From 46d410d7b0ae01d79747e7db2a9d652436bb39f5 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Sat, 25 Jul 2020 10:33:00 +0100 Subject: [PATCH 01/40] Add plugixml as a submodule This aligns with XML handling used by OpenMC. --- .gitmodules | 3 +++ pugixml | 1 + 2 files changed, 4 insertions(+) create mode 160000 pugixml diff --git a/.gitmodules b/.gitmodules index 2726187..02ee8fa 100644 --- a/.gitmodules +++ b/.gitmodules @@ -1,3 +1,6 @@ [submodule "pybind11"] path = pybind11 url = https://github.com/pybind/pybind11 +[submodule "pugixml"] + path = pugixml + url = https://github.com/zeux/pugixml.git diff --git a/pugixml b/pugixml new file mode 160000 index 0000000..22401ba --- /dev/null +++ b/pugixml @@ -0,0 +1 @@ +Subproject commit 22401bafaff996baecfb694ddc754855e184d377 From 1e9a30574f6c4c054cef57caa4532b43b651a626 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Sat, 25 Jul 2020 10:51:29 +0100 Subject: [PATCH 02/40] Get submodules when CMake-ing --- CMakeLists.txt | 25 +++++++++++++++++++++++++ 1 file changed, 25 insertions(+) diff --git a/CMakeLists.txt b/CMakeLists.txt index 578d4ce..327cacb 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -18,6 +18,31 @@ set(OPENMC_DIR /opt/openmc) set(OPENMC_INC_DIR ${OPENMC_DIR}/include) set(OPENMC_LIB_DIR ${OPENMC_DIR}/lib) +#=============================================================================== +# Update git submodules as needed +#=============================================================================== + +find_package(Git) +if(GIT_FOUND AND EXISTS "${CMAKE_CURRENT_SOURCE_DIR}/.git") + option(GIT_SUBMODULE "Check submodules during build" ON) + if(GIT_SUBMODULE) + message(STATUS "Submodule update") + execute_process(COMMAND ${GIT_EXECUTABLE} submodule update --init --recursive + WORKING_DIRECTORY ${CMAKE_CURRENT_SOURCE_DIR} + RESULT_VARIABLE GIT_SUBMOD_RESULT) + if(NOT GIT_SUBMOD_RESULT EQUAL 0) + message(FATAL_ERROR "git submodule update --init failed with \ + ${GIT_SUBMOD_RESULT}, please checkout submodules") + endif() + endif() +endif() + +# Check to see if submodules exist (by checking one) +if(NOT EXISTS "${CMAKE_CURRENT_SOURCE_DIR}/pugixml/CMakeLists.txt") + message(FATAL_ERROR "The git submodules were not downloaded! GIT_SUBMODULE was \ + turned off or failed. Please update submodules and try again.") +endif() + # Build source_sampling list(APPEND source_sampling_SOURCES ${SRC_DIR}/source_sampling.cpp From 2f9b8de628cf0a32011be8f822d197da9b88a8aa Mon Sep 17 00:00:00 2001 From: Dan Short Date: Sat, 25 Jul 2020 10:54:30 +0100 Subject: [PATCH 03/40] Template serialisation routines Introduces a to_xml serialisation routine that can write to string or to file. Creates the Python bindings for that routine. --- parametric_plasma_source/plasma_source.cpp | 10 ++++++++++ parametric_plasma_source/plasma_source.hpp | 4 ++++ parametric_plasma_source/plasma_source_pybind.cpp | 9 ++++++++- 3 files changed, 22 insertions(+), 1 deletion(-) diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index 3a2fc64..120fdb4 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -282,4 +282,14 @@ void PlasmaSource::isotropic_direction(const double random1, return; } +std::string PlasmaSource::to_xml() { + return ""; +} + +bool PlasmaSource::to_xml(std::string output_path) { + bool success = true; + std::string output = this->to_xml(); + return success; +} + } // end of namespace diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index bf1dbf1..3fcc8b2 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -91,6 +91,10 @@ void convert_r_to_xy(const double r, const double rn_store, void isotropic_direction(const double random1, const double random2, double &u, double &v, double &w); +std::string to_xml(); + +bool to_xml(std::string output_path); + private: std::vector source_profile; std::vector ion_kt; diff --git a/parametric_plasma_source/plasma_source_pybind.cpp b/parametric_plasma_source/plasma_source_pybind.cpp index c220d2a..a73cb79 100644 --- a/parametric_plasma_source/plasma_source_pybind.cpp +++ b/parametric_plasma_source/plasma_source_pybind.cpp @@ -44,5 +44,12 @@ PYBIND11_MODULE(plasma_source, m) { .def("dt_xs", &ps::PlasmaSource::dt_xs, "determine the value of the dt xs cross sections at a specific ion temperature", - py::arg("ion_temperature")); + py::arg("ion_temperature")) + .def("to_xml", + (bool (ps::PlasmaSource::*)(std::string)) &ps::PlasmaSource::to_xml, + "Serialise the PlasmaSource to XML", + py::arg("output_path")) + .def("to_xml", + (std::string (ps::PlasmaSource::*)()) &ps::PlasmaSource::to_xml, + "Serialise the PlasmaSource to XML"); } From c359d751b772b5a1ff234f75171afddbb33b50cb Mon Sep 17 00:00:00 2001 From: Dan Short Date: Mon, 27 Jul 2020 15:34:03 +0100 Subject: [PATCH 04/40] Serialise to XML Adds functionality to serialise a PlasmaSource to XML either as a string or to a file. Implements XML handling using pugixml and uses a hanlder class for XML operations. --- CMakeLists.txt | 10 +++++ parametric_plasma_source/plasma_source.cpp | 36 ++++++++++++++++- parametric_plasma_source/xml_helper.cpp | 47 ++++++++++++++++++++++ parametric_plasma_source/xml_helper.hpp | 20 +++++++++ setup.py | 2 + 5 files changed, 113 insertions(+), 2 deletions(-) create mode 100644 parametric_plasma_source/xml_helper.cpp create mode 100644 parametric_plasma_source/xml_helper.hpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 327cacb..f6d414c 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -43,6 +43,12 @@ if(NOT EXISTS "${CMAKE_CURRENT_SOURCE_DIR}/pugixml/CMakeLists.txt") turned off or failed. Please update submodules and try again.") endif() +#=============================================================================== +# pugixml library +#=============================================================================== + +add_subdirectory(pugixml) + # Build source_sampling list(APPEND source_sampling_SOURCES ${SRC_DIR}/source_sampling.cpp @@ -61,6 +67,7 @@ target_link_libraries(source_sampling ${OPENMC_LIB} gfortran) # Build plasma_source Python bindings list(APPEND plasma_source_pybind_SOURCES + ${SRC_DIR}/xml_helper.cpp ${SRC_DIR}/plasma_source.cpp ${SRC_DIR}/plasma_source_pybind.cpp ) @@ -68,3 +75,6 @@ list(APPEND plasma_source_pybind_SOURCES add_subdirectory(pybind11) pybind11_add_module(plasma_source ${plasma_source_pybind_SOURCES}) +target_include_directories(plasma_source PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/pugixml) +link_directories(${CMAKE_ARCHIVE_OUTPUT_DIRECTORY}) +target_link_libraries(plasma_source PRIVATE pugixml) diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index 120fdb4..97e1465 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -1,9 +1,11 @@ #include +#include #include #include #include "plasma_source.hpp" #include #include "openmc/random_lcg.h" +#include "xml_helper.hpp" #define RANDOM openmc::prn() @@ -283,12 +285,42 @@ void PlasmaSource::isotropic_direction(const double random1, } std::string PlasmaSource::to_xml() { - return ""; + XMLHelper xml_helper = XMLHelper("PlasmaSource"); + xml_helper.add_element("MinorRadius", minorRadius); + xml_helper.add_element("MajorRadius", majorRadius); + xml_helper.add_element("Elongation", elongation); + xml_helper.add_element("Triangularity", triangularity); + xml_helper.add_element("ShafranovShift", shafranov); + xml_helper.add_element("PedestalRadius", pedistalRadius); + xml_helper.add_element("IonDensityPedestal", ionDensityPedistal); + xml_helper.add_element("IonDensitySeparatrix", ionDensitySeperatrix); + xml_helper.add_element("IonDensityOrigin", ionDensityOrigin); + xml_helper.add_element("IonTemperaturePedestal", ionTemperaturePedistal); + xml_helper.add_element("IonTemperatureSeparatrix", ionTemperatureSeperatrix); + xml_helper.add_element("IonTemperatureOrigin", ionTemperatureOrigin); + xml_helper.add_element("IonDensityAlpha", ionDensityPeaking); + xml_helper.add_element("IonTemperatureAlpha", ionTemperaturePeaking); + xml_helper.add_element("IonTemperatureBeta", ionTemperatureBeta); + xml_helper.add_element("PlasmaType", plasmaType); + xml_helper.add_element("PlasmaId", plasmaId); + xml_helper.add_element("NumberOfBins", numberOfBins); + xml_helper.add_element("MinimumToroidalAngle", minToroidalAngle); + xml_helper.add_element("MaximumToroidalAngle", maxToroidalAngle); + return xml_helper.to_string(); } bool PlasmaSource::to_xml(std::string output_path) { bool success = true; - std::string output = this->to_xml(); + std::string output = to_xml(); + std::ofstream out; + out.open(output_path); + if (out.is_open()) { + out << output; + out.close(); + } + else { + success = false; + } return success; } diff --git a/parametric_plasma_source/xml_helper.cpp b/parametric_plasma_source/xml_helper.cpp new file mode 100644 index 0000000..a70fe1b --- /dev/null +++ b/parametric_plasma_source/xml_helper.cpp @@ -0,0 +1,47 @@ +#include +#include +#include "xml_helper.hpp" +#include "pugixml.hpp" + +namespace plasma_source { + + XMLHelper::XMLHelper(const char* root_name) { + root_node = doc.append_child(root_name); + } + + void XMLHelper::set_root_element(const char* root_name) { + root_node = doc.append_child(root_name); + } + + void XMLHelper::add_element(const char* element_name, const char* element_value) { + root_node.append_child(element_name).text() = element_value; + } + + void XMLHelper::add_element(const char* element_name, double element_value) { + root_node.append_child(element_name).text() = element_value; + } + + void XMLHelper::add_element(const char* element_name, int32_t element_value) { + root_node.append_child(element_name).text() = element_value; + } + + void XMLHelper::add_element(const char* element_name, std::string element_value) { + root_node.append_child(element_name).text() = element_value.c_str(); + } + + struct xml_string_writer: pugi::xml_writer { + std::string result; + + virtual void write(const void* data, size_t size) { + result.append(static_cast(data), size); + } + }; + + std::string XMLHelper::to_string() { + xml_string_writer writer; + doc.print(writer, " "); + + return writer.result; + } + +} \ No newline at end of file diff --git a/parametric_plasma_source/xml_helper.hpp b/parametric_plasma_source/xml_helper.hpp new file mode 100644 index 0000000..7fe056a --- /dev/null +++ b/parametric_plasma_source/xml_helper.hpp @@ -0,0 +1,20 @@ +#include "pugixml.hpp" + +namespace plasma_source { + +class XMLHelper { + public: + XMLHelper(const char* root_name); + void set_root_element(const char* root_name); + void add_element(const char* element_name, const char* element_value); + void add_element(const char* element_name, double element_value); + void add_element(const char* element_name, int element_value); + void add_element(const char* element_name, std::string element_value); + std::string to_string(); + + private: + pugi::xml_document doc; + pugi::xml_node root_node; +}; + +} diff --git a/setup.py b/setup.py index af274c5..698cf89 100644 --- a/setup.py +++ b/setup.py @@ -33,6 +33,7 @@ def build_extension(self, ext): extdir += os.path.sep cmake_args = ["-DCMAKE_LIBRARY_OUTPUT_DIRECTORY=" + extdir, + "-DCMAKE_ARCHIVE_OUTPUT_DIRECTORY=" + extdir, "-DPYTHON_EXECUTABLE=" + sys.executable] cfg = "Debug" if self.debug else "Release" @@ -77,6 +78,7 @@ def build_extension(self, ext): url="https://github.com/makeclean/parametric-plasma-source/", packages=find_packages(), ext_modules=[CMakeExtention("parametric_plasma_source/plasma_source")], + package_data={"parametric_plasma_source": ["libpugixml.a"]}, cmdclass=dict(build_ext=CMakeBuild), classifiers=[ "Programming Language :: Python :: 3", From 301bec0f618765435edcbf567e8f22354a2909c0 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Mon, 27 Jul 2020 15:34:28 +0100 Subject: [PATCH 05/40] Example serialised representation. --- data/example_source.xml | 22 ++++++++++++++++++++++ 1 file changed, 22 insertions(+) create mode 100644 data/example_source.xml diff --git a/data/example_source.xml b/data/example_source.xml new file mode 100644 index 0000000..706a10f --- /dev/null +++ b/data/example_source.xml @@ -0,0 +1,22 @@ + + 1.5 + 4.5 + 2 + 0.55 + 0 + 0.8 + 1.09e+20 + 3e+19 + 1.09e+20 + 6.09 + 0.1 + 45.9 + 1 + 8.06 + 6 + plasma + 1 + 100 + 0 + 6.2831853071795862 + From 31245024492822f2e728224fb87da23971f460b5 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Mon, 27 Jul 2020 16:27:01 +0100 Subject: [PATCH 06/40] Deserialise the plasma source from XML Deserialises the plasma source from an XML file and create an example input file. --- data/example_source.xml | 2 +- parametric_plasma_source/plasma_source.cpp | 40 ++++++++++++++++++- parametric_plasma_source/plasma_source.hpp | 2 + .../plasma_source_pybind.cpp | 6 ++- 4 files changed, 47 insertions(+), 3 deletions(-) diff --git a/data/example_source.xml b/data/example_source.xml index 706a10f..ef57f9a 100644 --- a/data/example_source.xml +++ b/data/example_source.xml @@ -19,4 +19,4 @@ 100 0 6.2831853071795862 - + diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index 97e1465..ae14334 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -286,8 +286,8 @@ void PlasmaSource::isotropic_direction(const double random1, std::string PlasmaSource::to_xml() { XMLHelper xml_helper = XMLHelper("PlasmaSource"); - xml_helper.add_element("MinorRadius", minorRadius); xml_helper.add_element("MajorRadius", majorRadius); + xml_helper.add_element("MinorRadius", minorRadius); xml_helper.add_element("Elongation", elongation); xml_helper.add_element("Triangularity", triangularity); xml_helper.add_element("ShafranovShift", shafranov); @@ -324,4 +324,42 @@ bool PlasmaSource::to_xml(std::string output_path) { return success; } +PlasmaSource PlasmaSource::from_file(std::string input_path) { + pugi::xml_document doc; + pugi::xml_parse_result result = doc.load_file(input_path.c_str()); + if (!result) { + throw std::runtime_error("Error reading plasma source from file " + input_path + ". Error = " + result.description()); + } + + pugi::xml_node root_node = doc.root().child("PlasmaSource"); + + double major_radius = root_node.child("MajorRadius").text().as_double(); + double minor_radius = root_node.child("MinorRadius").text().as_double(); + double elongation = root_node.child("Elongation").text().as_double(); + double triangularity = root_node.child("Triangularity").text().as_double(); + double shafranov_shift = root_node.child("ShafranovShift").text().as_double(); + double pedestal_radius = root_node.child("PedestalRadius").text().as_double(); + double ion_density_pedestal = root_node.child("IonDensityPedestal").text().as_double(); + double ion_density_separatrix = root_node.child("IonDensitySeparatrix").text().as_double(); + double ion_density_origin = root_node.child("IonDensityOrigin").text().as_double(); + double ion_temperature_pedestal = root_node.child("IonTemperaturePedestal").text().as_double(); + double ion_temperature_separatrix = root_node.child("IonTemperatureSeparatrix").text().as_double(); + double ion_temperature_origin = root_node.child("IonTemperatureOrigin").text().as_double(); + double ion_density_alpha = root_node.child("IonDensityAlpha").text().as_double(); + double ion_temperature_alpha = root_node.child("IonTemperatureAlpha").text().as_double(); + double ion_temperature_beta = root_node.child("IonTemperatureBeta").text().as_double(); + std::string plasma_type = root_node.child("PlasmaType").text().as_string(); + int plasma_id = root_node.child("PlasmaId").text().as_int(); + int number_of_bins = root_node.child("NumberOfBins").text().as_int(); + double min_toroidal_angle = root_node.child("MinimumToroidalAngle").text().as_double(); + double max_toroidal_angle = root_node.child("MaximumToroidalAngle").text().as_double(); + + PlasmaSource plasma_source = PlasmaSource( + ion_density_pedestal, ion_density_separatrix, ion_density_origin, ion_temperature_pedestal, ion_temperature_separatrix, + ion_temperature_origin, pedestal_radius, ion_density_alpha, ion_temperature_alpha, ion_temperature_beta, minor_radius, + major_radius, elongation, triangularity, shafranov_shift, plasma_type, plasma_id, number_of_bins, min_toroidal_angle, + max_toroidal_angle); + return plasma_source; +} + } // end of namespace diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index 3fcc8b2..9e1af42 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -95,6 +95,8 @@ std::string to_xml(); bool to_xml(std::string output_path); +static PlasmaSource from_file(std::string input_path); + private: std::vector source_profile; std::vector ion_kt; diff --git a/parametric_plasma_source/plasma_source_pybind.cpp b/parametric_plasma_source/plasma_source_pybind.cpp index a73cb79..a99814e 100644 --- a/parametric_plasma_source/plasma_source_pybind.cpp +++ b/parametric_plasma_source/plasma_source_pybind.cpp @@ -51,5 +51,9 @@ PYBIND11_MODULE(plasma_source, m) { py::arg("output_path")) .def("to_xml", (std::string (ps::PlasmaSource::*)()) &ps::PlasmaSource::to_xml, - "Serialise the PlasmaSource to XML"); + "Serialise the PlasmaSource to XML") + .def_static("from_file", + &ps::PlasmaSource::from_file, + "Deserialise the PlasmaSource from an XML file", + py::arg("input_path")); } From 0ff20afd3c0fc0729dee25d04254601df5e51182 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Mon, 27 Jul 2020 16:51:37 +0100 Subject: [PATCH 07/40] Fix example source --- data/example_source.xml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/data/example_source.xml b/data/example_source.xml index ef57f9a..706a10f 100644 --- a/data/example_source.xml +++ b/data/example_source.xml @@ -19,4 +19,4 @@ 100 0 6.2831853071795862 - + From 33bd0d38a9a5f27923e7413eea67ea4b138fdc33 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Mon, 27 Jul 2020 16:53:07 +0100 Subject: [PATCH 08/40] Deserialise from an XML string --- parametric_plasma_source/plasma_source.cpp | 40 +++++++++++++++++++ parametric_plasma_source/plasma_source.hpp | 2 + .../plasma_source_pybind.cpp | 6 ++- 3 files changed, 47 insertions(+), 1 deletion(-) diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index ae14334..7c3002b 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -1,5 +1,6 @@ #include #include +#include #include #include #include "plasma_source.hpp" @@ -362,4 +363,43 @@ PlasmaSource PlasmaSource::from_file(std::string input_path) { return plasma_source; } +PlasmaSource PlasmaSource::from_xml(std::string xml) { + pugi::xml_document doc; + std::istringstream ss(xml); + pugi::xml_parse_result result = doc.load(ss); + if (!result) { + throw std::runtime_error("Error reading plasma source from string " + xml + ". Error = " + result.description()); + } + + pugi::xml_node root_node = doc.root().child("PlasmaSource"); + + double major_radius = root_node.child("MajorRadius").text().as_double(); + double minor_radius = root_node.child("MinorRadius").text().as_double(); + double elongation = root_node.child("Elongation").text().as_double(); + double triangularity = root_node.child("Triangularity").text().as_double(); + double shafranov_shift = root_node.child("ShafranovShift").text().as_double(); + double pedestal_radius = root_node.child("PedestalRadius").text().as_double(); + double ion_density_pedestal = root_node.child("IonDensityPedestal").text().as_double(); + double ion_density_separatrix = root_node.child("IonDensitySeparatrix").text().as_double(); + double ion_density_origin = root_node.child("IonDensityOrigin").text().as_double(); + double ion_temperature_pedestal = root_node.child("IonTemperaturePedestal").text().as_double(); + double ion_temperature_separatrix = root_node.child("IonTemperatureSeparatrix").text().as_double(); + double ion_temperature_origin = root_node.child("IonTemperatureOrigin").text().as_double(); + double ion_density_alpha = root_node.child("IonDensityAlpha").text().as_double(); + double ion_temperature_alpha = root_node.child("IonTemperatureAlpha").text().as_double(); + double ion_temperature_beta = root_node.child("IonTemperatureBeta").text().as_double(); + std::string plasma_type = root_node.child("PlasmaType").text().as_string(); + int plasma_id = root_node.child("PlasmaId").text().as_int(); + int number_of_bins = root_node.child("NumberOfBins").text().as_int(); + double min_toroidal_angle = root_node.child("MinimumToroidalAngle").text().as_double(); + double max_toroidal_angle = root_node.child("MaximumToroidalAngle").text().as_double(); + + PlasmaSource plasma_source = PlasmaSource( + ion_density_pedestal, ion_density_separatrix, ion_density_origin, ion_temperature_pedestal, ion_temperature_separatrix, + ion_temperature_origin, pedestal_radius, ion_density_alpha, ion_temperature_alpha, ion_temperature_beta, minor_radius, + major_radius, elongation, triangularity, shafranov_shift, plasma_type, plasma_id, number_of_bins, min_toroidal_angle, + max_toroidal_angle); + return plasma_source; +} + } // end of namespace diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index 9e1af42..e7e0c63 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -95,6 +95,8 @@ std::string to_xml(); bool to_xml(std::string output_path); +static PlasmaSource from_xml(std::string xml); + static PlasmaSource from_file(std::string input_path); private: diff --git a/parametric_plasma_source/plasma_source_pybind.cpp b/parametric_plasma_source/plasma_source_pybind.cpp index a99814e..dbf1cad 100644 --- a/parametric_plasma_source/plasma_source_pybind.cpp +++ b/parametric_plasma_source/plasma_source_pybind.cpp @@ -55,5 +55,9 @@ PYBIND11_MODULE(plasma_source, m) { .def_static("from_file", &ps::PlasmaSource::from_file, "Deserialise the PlasmaSource from an XML file", - py::arg("input_path")); + py::arg("input_path")) + .def_static("from_xml", + &ps::PlasmaSource::from_xml, + "Deserialise the PlasmaSource from an XML string", + py::arg("xml")); } From d3ebbd4a23b4ce7d645cea44bd069b29e7f23285 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Mon, 27 Jul 2020 17:36:05 +0100 Subject: [PATCH 09/40] Remove duplicate code and simplify Moves duplicate code into private function and passed XML doc as pointer. Adds helper functions for getters. Moves root name to const. --- parametric_plasma_source/plasma_source.cpp | 92 ++++++++-------------- parametric_plasma_source/plasma_source.hpp | 6 ++ parametric_plasma_source/xml_helper.cpp | 17 ++++ parametric_plasma_source/xml_helper.hpp | 4 + 4 files changed, 61 insertions(+), 58 deletions(-) diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index 7c3002b..cda4516 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -7,6 +7,7 @@ #include #include "openmc/random_lcg.h" #include "xml_helper.hpp" +#include "pugixml.hpp" #define RANDOM openmc::prn() @@ -325,35 +326,28 @@ bool PlasmaSource::to_xml(std::string output_path) { return success; } -PlasmaSource PlasmaSource::from_file(std::string input_path) { - pugi::xml_document doc; - pugi::xml_parse_result result = doc.load_file(input_path.c_str()); - if (!result) { - throw std::runtime_error("Error reading plasma source from file " + input_path + ". Error = " + result.description()); - } - - pugi::xml_node root_node = doc.root().child("PlasmaSource"); - - double major_radius = root_node.child("MajorRadius").text().as_double(); - double minor_radius = root_node.child("MinorRadius").text().as_double(); - double elongation = root_node.child("Elongation").text().as_double(); - double triangularity = root_node.child("Triangularity").text().as_double(); - double shafranov_shift = root_node.child("ShafranovShift").text().as_double(); - double pedestal_radius = root_node.child("PedestalRadius").text().as_double(); - double ion_density_pedestal = root_node.child("IonDensityPedestal").text().as_double(); - double ion_density_separatrix = root_node.child("IonDensitySeparatrix").text().as_double(); - double ion_density_origin = root_node.child("IonDensityOrigin").text().as_double(); - double ion_temperature_pedestal = root_node.child("IonTemperaturePedestal").text().as_double(); - double ion_temperature_separatrix = root_node.child("IonTemperatureSeparatrix").text().as_double(); - double ion_temperature_origin = root_node.child("IonTemperatureOrigin").text().as_double(); - double ion_density_alpha = root_node.child("IonDensityAlpha").text().as_double(); - double ion_temperature_alpha = root_node.child("IonTemperatureAlpha").text().as_double(); - double ion_temperature_beta = root_node.child("IonTemperatureBeta").text().as_double(); - std::string plasma_type = root_node.child("PlasmaType").text().as_string(); - int plasma_id = root_node.child("PlasmaId").text().as_int(); - int number_of_bins = root_node.child("NumberOfBins").text().as_int(); - double min_toroidal_angle = root_node.child("MinimumToroidalAngle").text().as_double(); - double max_toroidal_angle = root_node.child("MaximumToroidalAngle").text().as_double(); +PlasmaSource PlasmaSource::from_xml_doc(pugi::xml_document* xml_doc, std::string root_name) { + XMLHelper xml_helper = XMLHelper(xml_doc, root_name); + double major_radius = xml_helper.get_element_as_double("MajorRadius"); + double minor_radius = xml_helper.get_element_as_double("MinorRadius"); + double elongation = xml_helper.get_element_as_double("Elongation"); + double triangularity = xml_helper.get_element_as_double("Triangularity"); + double shafranov_shift = xml_helper.get_element_as_double("ShafranovShift"); + double pedestal_radius = xml_helper.get_element_as_double("PedestalRadius"); + double ion_density_pedestal = xml_helper.get_element_as_double("IonDensityPedestal"); + double ion_density_separatrix = xml_helper.get_element_as_double("IonDensitySeparatrix"); + double ion_density_origin = xml_helper.get_element_as_double("IonDensityOrigin"); + double ion_temperature_pedestal = xml_helper.get_element_as_double("IonTemperaturePedestal"); + double ion_temperature_separatrix = xml_helper.get_element_as_double("IonTemperatureSeparatrix"); + double ion_temperature_origin = xml_helper.get_element_as_double("IonTemperatureOrigin"); + double ion_density_alpha = xml_helper.get_element_as_double("IonDensityAlpha"); + double ion_temperature_alpha = xml_helper.get_element_as_double("IonTemperatureAlpha"); + double ion_temperature_beta = xml_helper.get_element_as_double("IonTemperatureBeta"); + std::string plasma_type = xml_helper.get_element_as_string("PlasmaType"); + int plasma_id = xml_helper.get_element_as_int("PlasmaId"); + int number_of_bins = xml_helper.get_element_as_int("NumberOfBins"); + double min_toroidal_angle = xml_helper.get_element_as_double("MinimumToroidalAngle"); + double max_toroidal_angle = xml_helper.get_element_as_double("MaximumToroidalAngle"); PlasmaSource plasma_source = PlasmaSource( ion_density_pedestal, ion_density_separatrix, ion_density_origin, ion_temperature_pedestal, ion_temperature_separatrix, @@ -363,6 +357,16 @@ PlasmaSource PlasmaSource::from_file(std::string input_path) { return plasma_source; } +PlasmaSource PlasmaSource::from_file(std::string input_path) { + pugi::xml_document doc; + pugi::xml_parse_result result = doc.load_file(input_path.c_str()); + if (!result) { + throw std::runtime_error("Error reading plasma source from file " + input_path + ". Error = " + result.description()); + } + + return from_xml_doc(&doc, PLASMA_SOURCE_ROOT_NAME); +} + PlasmaSource PlasmaSource::from_xml(std::string xml) { pugi::xml_document doc; std::istringstream ss(xml); @@ -371,35 +375,7 @@ PlasmaSource PlasmaSource::from_xml(std::string xml) { throw std::runtime_error("Error reading plasma source from string " + xml + ". Error = " + result.description()); } - pugi::xml_node root_node = doc.root().child("PlasmaSource"); - - double major_radius = root_node.child("MajorRadius").text().as_double(); - double minor_radius = root_node.child("MinorRadius").text().as_double(); - double elongation = root_node.child("Elongation").text().as_double(); - double triangularity = root_node.child("Triangularity").text().as_double(); - double shafranov_shift = root_node.child("ShafranovShift").text().as_double(); - double pedestal_radius = root_node.child("PedestalRadius").text().as_double(); - double ion_density_pedestal = root_node.child("IonDensityPedestal").text().as_double(); - double ion_density_separatrix = root_node.child("IonDensitySeparatrix").text().as_double(); - double ion_density_origin = root_node.child("IonDensityOrigin").text().as_double(); - double ion_temperature_pedestal = root_node.child("IonTemperaturePedestal").text().as_double(); - double ion_temperature_separatrix = root_node.child("IonTemperatureSeparatrix").text().as_double(); - double ion_temperature_origin = root_node.child("IonTemperatureOrigin").text().as_double(); - double ion_density_alpha = root_node.child("IonDensityAlpha").text().as_double(); - double ion_temperature_alpha = root_node.child("IonTemperatureAlpha").text().as_double(); - double ion_temperature_beta = root_node.child("IonTemperatureBeta").text().as_double(); - std::string plasma_type = root_node.child("PlasmaType").text().as_string(); - int plasma_id = root_node.child("PlasmaId").text().as_int(); - int number_of_bins = root_node.child("NumberOfBins").text().as_int(); - double min_toroidal_angle = root_node.child("MinimumToroidalAngle").text().as_double(); - double max_toroidal_angle = root_node.child("MaximumToroidalAngle").text().as_double(); - - PlasmaSource plasma_source = PlasmaSource( - ion_density_pedestal, ion_density_separatrix, ion_density_origin, ion_temperature_pedestal, ion_temperature_separatrix, - ion_temperature_origin, pedestal_radius, ion_density_alpha, ion_temperature_alpha, ion_temperature_beta, minor_radius, - major_radius, elongation, triangularity, shafranov_shift, plasma_type, plasma_id, number_of_bins, min_toroidal_angle, - max_toroidal_angle); - return plasma_source; + return from_xml_doc(&doc, PLASMA_SOURCE_ROOT_NAME); } } // end of namespace diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index e7e0c63..1cdfd0e 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -1,11 +1,15 @@ #include #include +#include "pugixml.hpp" + namespace plasma_source { struct xs_params { double c[7]; }; +static const std::string PLASMA_SOURCE_ROOT_NAME = "PlasmaSource"; + class PlasmaSource { public: // constructor @@ -126,5 +130,7 @@ static PlasmaSource from_file(std::string input_path); double binWidth; int numberOfBins; + static PlasmaSource from_xml_doc(pugi::xml_document* xml_doc, std::string root_name); + }; }// end of namespace \ No newline at end of file diff --git a/parametric_plasma_source/xml_helper.cpp b/parametric_plasma_source/xml_helper.cpp index a70fe1b..d4fd973 100644 --- a/parametric_plasma_source/xml_helper.cpp +++ b/parametric_plasma_source/xml_helper.cpp @@ -9,6 +9,11 @@ namespace plasma_source { root_node = doc.append_child(root_name); } + XMLHelper::XMLHelper(pugi::xml_document* doc, std::string root_name) { + this->doc.reset(*doc); + this->root_node = this->doc.root().child(root_name.c_str()); + } + void XMLHelper::set_root_element(const char* root_name) { root_node = doc.append_child(root_name); } @@ -29,6 +34,18 @@ namespace plasma_source { root_node.append_child(element_name).text() = element_value.c_str(); } + double XMLHelper::get_element_as_double(const char* element_name) { + return root_node.child(element_name).text().as_double(); + } + + int XMLHelper::get_element_as_int(const char* element_name) { + return root_node.child(element_name).text().as_int(); + } + + std::string XMLHelper::get_element_as_string(const char* element_name) { + return root_node.child(element_name).text().as_string(); + } + struct xml_string_writer: pugi::xml_writer { std::string result; diff --git a/parametric_plasma_source/xml_helper.hpp b/parametric_plasma_source/xml_helper.hpp index 7fe056a..6917d36 100644 --- a/parametric_plasma_source/xml_helper.hpp +++ b/parametric_plasma_source/xml_helper.hpp @@ -5,11 +5,15 @@ namespace plasma_source { class XMLHelper { public: XMLHelper(const char* root_name); + XMLHelper(pugi::xml_document* doc, std::string root_name); void set_root_element(const char* root_name); void add_element(const char* element_name, const char* element_value); void add_element(const char* element_name, double element_value); void add_element(const char* element_name, int element_value); void add_element(const char* element_name, std::string element_value); + double get_element_as_double(const char* element_name); + int get_element_as_int(const char* element_name); + std::string get_element_as_string(const char* element_name); std::string to_string(); private: From d0263bb3654e44ce7474eab172b188058b8432d7 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Mon, 27 Jul 2020 17:37:09 +0100 Subject: [PATCH 10/40] Use root name consistently --- parametric_plasma_source/plasma_source.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index cda4516..8c1d599 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -287,7 +287,7 @@ void PlasmaSource::isotropic_direction(const double random1, } std::string PlasmaSource::to_xml() { - XMLHelper xml_helper = XMLHelper("PlasmaSource"); + XMLHelper xml_helper = XMLHelper(PLASMA_SOURCE_ROOT_NAME.c_str()); xml_helper.add_element("MajorRadius", majorRadius); xml_helper.add_element("MinorRadius", minorRadius); xml_helper.add_element("Elongation", elongation); From 209f35f0dda121f5418d8765ada87ba17619f7b5 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Mon, 27 Jul 2020 17:39:01 +0100 Subject: [PATCH 11/40] Specify max toroidal angle in degrees --- data/example_source.xml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/data/example_source.xml b/data/example_source.xml index 706a10f..6036672 100644 --- a/data/example_source.xml +++ b/data/example_source.xml @@ -18,5 +18,5 @@ 1 100 0 - 6.2831853071795862 + 360 From 9a7d4b76547275bf5fe7cebdd57265c611f640f0 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Tue, 28 Jul 2020 09:20:21 +0100 Subject: [PATCH 12/40] Serialise angles in degrees --- parametric_plasma_source/plasma_source.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index 8c1d599..8dea34b 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -306,8 +306,8 @@ std::string PlasmaSource::to_xml() { xml_helper.add_element("PlasmaType", plasmaType); xml_helper.add_element("PlasmaId", plasmaId); xml_helper.add_element("NumberOfBins", numberOfBins); - xml_helper.add_element("MinimumToroidalAngle", minToroidalAngle); - xml_helper.add_element("MaximumToroidalAngle", maxToroidalAngle); + xml_helper.add_element("MinimumToroidalAngle", minToroidalAngle / M_PI * 180.0); + xml_helper.add_element("MaximumToroidalAngle", maxToroidalAngle / M_PI * 180.0); return xml_helper.to_string(); } From a474854dbfee6bb620e0981ea8ec087a55d858ff Mon Sep 17 00:00:00 2001 From: Dan Short Date: Fri, 31 Jul 2020 15:03:28 +0100 Subject: [PATCH 13/40] Updates to work with OpenMC deserialisation Align with work in progress on OpenMC to support deserialisation of a serialised source. In particular remove XMLHelper in deserialisation as this causes some symbols to be unfindable. Also removes default values and object from sampling as no longer needed. --- parametric_plasma_source/plasma_source.cpp | 51 ++++++++++---------- parametric_plasma_source/plasma_source.hpp | 2 +- parametric_plasma_source/source_sampling.cpp | 48 ++---------------- parametric_plasma_source/xml_helper.cpp | 17 ------- parametric_plasma_source/xml_helper.hpp | 4 -- 5 files changed, 30 insertions(+), 92 deletions(-) diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index 8dea34b..96d6e08 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -327,27 +327,29 @@ bool PlasmaSource::to_xml(std::string output_path) { } PlasmaSource PlasmaSource::from_xml_doc(pugi::xml_document* xml_doc, std::string root_name) { - XMLHelper xml_helper = XMLHelper(xml_doc, root_name); - double major_radius = xml_helper.get_element_as_double("MajorRadius"); - double minor_radius = xml_helper.get_element_as_double("MinorRadius"); - double elongation = xml_helper.get_element_as_double("Elongation"); - double triangularity = xml_helper.get_element_as_double("Triangularity"); - double shafranov_shift = xml_helper.get_element_as_double("ShafranovShift"); - double pedestal_radius = xml_helper.get_element_as_double("PedestalRadius"); - double ion_density_pedestal = xml_helper.get_element_as_double("IonDensityPedestal"); - double ion_density_separatrix = xml_helper.get_element_as_double("IonDensitySeparatrix"); - double ion_density_origin = xml_helper.get_element_as_double("IonDensityOrigin"); - double ion_temperature_pedestal = xml_helper.get_element_as_double("IonTemperaturePedestal"); - double ion_temperature_separatrix = xml_helper.get_element_as_double("IonTemperatureSeparatrix"); - double ion_temperature_origin = xml_helper.get_element_as_double("IonTemperatureOrigin"); - double ion_density_alpha = xml_helper.get_element_as_double("IonDensityAlpha"); - double ion_temperature_alpha = xml_helper.get_element_as_double("IonTemperatureAlpha"); - double ion_temperature_beta = xml_helper.get_element_as_double("IonTemperatureBeta"); - std::string plasma_type = xml_helper.get_element_as_string("PlasmaType"); - int plasma_id = xml_helper.get_element_as_int("PlasmaId"); - int number_of_bins = xml_helper.get_element_as_int("NumberOfBins"); - double min_toroidal_angle = xml_helper.get_element_as_double("MinimumToroidalAngle"); - double max_toroidal_angle = xml_helper.get_element_as_double("MaximumToroidalAngle"); + pugi::xml_document doc; + doc.reset(*xml_doc); + pugi::xml_node root_node = doc.root().child(root_name.c_str()); + double major_radius = root_node.child("MajorRadius").text().as_double(); + double minor_radius = root_node.child("MinorRadius").text().as_double(); + double elongation = root_node.child("Elongation").text().as_double(); + double triangularity = root_node.child("Triangularity").text().as_double(); + double shafranov_shift = root_node.child("ShafranovShift").text().as_double(); + double pedestal_radius = root_node.child("PedestalRadius").text().as_double(); + double ion_density_pedestal = root_node.child("IonDensityPedestal").text().as_double(); + double ion_density_separatrix = root_node.child("IonDensitySeparatrix").text().as_double(); + double ion_density_origin = root_node.child("IonDensityOrigin").text().as_double(); + double ion_temperature_pedestal = root_node.child("IonTemperaturePedestal").text().as_double(); + double ion_temperature_separatrix = root_node.child("IonTemperatureSeparatrix").text().as_double(); + double ion_temperature_origin = root_node.child("IonTemperatureOrigin").text().as_double(); + double ion_density_alpha = root_node.child("IonDensityAlpha").text().as_double(); + double ion_temperature_alpha = root_node.child("IonTemperatureAlpha").text().as_double(); + double ion_temperature_beta = root_node.child("IonTemperatureBeta").text().as_double(); + std::string plasma_type = root_node.child("PlasmaType").text().as_string(); + int plasma_id = root_node.child("PlasmaId").text().as_int(); + int number_of_bins = root_node.child("NumberOfBins").text().as_int(); + double min_toroidal_angle = root_node.child("MinimumToroidalAngle").text().as_double(); + double max_toroidal_angle = root_node.child("MaximumToroidalAngle").text().as_double(); PlasmaSource plasma_source = PlasmaSource( ion_density_pedestal, ion_density_separatrix, ion_density_origin, ion_temperature_pedestal, ion_temperature_separatrix, @@ -367,12 +369,11 @@ PlasmaSource PlasmaSource::from_file(std::string input_path) { return from_xml_doc(&doc, PLASMA_SOURCE_ROOT_NAME); } -PlasmaSource PlasmaSource::from_xml(std::string xml) { +PlasmaSource PlasmaSource::from_xml(char* xml) { pugi::xml_document doc; - std::istringstream ss(xml); - pugi::xml_parse_result result = doc.load(ss); + pugi::xml_parse_result result = doc.load_string(xml); if (!result) { - throw std::runtime_error("Error reading plasma source from string " + xml + ". Error = " + result.description()); + throw std::runtime_error("Error reading plasma source from string " + std::string(xml) + ". Error = " + result.description()); } return from_xml_doc(&doc, PLASMA_SOURCE_ROOT_NAME); diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index 1cdfd0e..eef37df 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -99,7 +99,7 @@ std::string to_xml(); bool to_xml(std::string output_path); -static PlasmaSource from_xml(std::string xml); +static PlasmaSource from_xml(char* xml); static PlasmaSource from_file(std::string input_path); diff --git a/parametric_plasma_source/source_sampling.cpp b/parametric_plasma_source/source_sampling.cpp index f763d4a..c0c198a 100644 --- a/parametric_plasma_source/source_sampling.cpp +++ b/parametric_plasma_source/source_sampling.cpp @@ -4,53 +4,11 @@ #include "openmc/particle.h" #include "plasma_source.hpp" - -// Spherical tokamak SOURCE -// units are in SI units -const double ion_density_pedistal = 1.09e+20; // ions per m^3 -const double ion_density_seperatrix = 3e+19; -const double ion_density_origin = 1.09e+20; -const double ion_temperature_pedistal = 6.09; -const double ion_temperature_seperatrix = 0.1; -const double ion_temperature_origin = 45.9; -const double ion_density_peaking_factor = 1; -const double ion_temperature_peaking_factor = 8.06; // check alpha or beta value from paper -const double ion_temperature_beta = 6.0; -const double minor_radius = 1.56; // metres -const double major_radius = 2.5; // metres -const double pedistal_radius = 0.8 * minor_radius; // pedistal minor rad in metres -const double elongation = 2.0; -const double triangularity = 0.55; -const double shafranov_shift = 0.0; //metres -const std::string name = "parametric_plasma_source"; -const int number_of_bins = 100; -const int plasma_type = 1; // 1 is default; //0 = L mode anything else H/A mode - - -plasma_source::PlasmaSource source = plasma_source::PlasmaSource(ion_density_pedistal, - ion_density_seperatrix, - ion_density_origin, - ion_temperature_pedistal, - ion_temperature_seperatrix, - ion_temperature_origin, - pedistal_radius, - ion_density_peaking_factor, - ion_temperature_peaking_factor, - ion_temperature_beta, - minor_radius, - major_radius, - elongation, - triangularity, - shafranov_shift, - name, - plasma_type, - number_of_bins, - 0.0, - 360.0); - // you must have external C linkage here otherwise // dlopen will not find the file -extern "C" openmc::Particle::Bank sample_source(uint64_t* seed) { +extern "C" openmc::Particle::Bank sample_source(uint64_t* seed, char* serialization) { + plasma_source::PlasmaSource source = plasma_source::PlasmaSource::from_xml(serialization); + openmc::Particle::Bank particle; // wgt particle.particle = openmc::Particle::Type::neutron; diff --git a/parametric_plasma_source/xml_helper.cpp b/parametric_plasma_source/xml_helper.cpp index d4fd973..a70fe1b 100644 --- a/parametric_plasma_source/xml_helper.cpp +++ b/parametric_plasma_source/xml_helper.cpp @@ -9,11 +9,6 @@ namespace plasma_source { root_node = doc.append_child(root_name); } - XMLHelper::XMLHelper(pugi::xml_document* doc, std::string root_name) { - this->doc.reset(*doc); - this->root_node = this->doc.root().child(root_name.c_str()); - } - void XMLHelper::set_root_element(const char* root_name) { root_node = doc.append_child(root_name); } @@ -34,18 +29,6 @@ namespace plasma_source { root_node.append_child(element_name).text() = element_value.c_str(); } - double XMLHelper::get_element_as_double(const char* element_name) { - return root_node.child(element_name).text().as_double(); - } - - int XMLHelper::get_element_as_int(const char* element_name) { - return root_node.child(element_name).text().as_int(); - } - - std::string XMLHelper::get_element_as_string(const char* element_name) { - return root_node.child(element_name).text().as_string(); - } - struct xml_string_writer: pugi::xml_writer { std::string result; diff --git a/parametric_plasma_source/xml_helper.hpp b/parametric_plasma_source/xml_helper.hpp index 6917d36..7fe056a 100644 --- a/parametric_plasma_source/xml_helper.hpp +++ b/parametric_plasma_source/xml_helper.hpp @@ -5,15 +5,11 @@ namespace plasma_source { class XMLHelper { public: XMLHelper(const char* root_name); - XMLHelper(pugi::xml_document* doc, std::string root_name); void set_root_element(const char* root_name); void add_element(const char* element_name, const char* element_value); void add_element(const char* element_name, double element_value); void add_element(const char* element_name, int element_value); void add_element(const char* element_name, std::string element_value); - double get_element_as_double(const char* element_name); - int get_element_as_int(const char* element_name); - std::string get_element_as_string(const char* element_name); std::string to_string(); private: From 1324d1a5446191cb1bcb60baf44c6268a2423d73 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Fri, 31 Jul 2020 15:09:20 +0100 Subject: [PATCH 14/40] Include source sampling in bindings Handles pass by reference using a lambda that maps return values to a tuple. Also tweak some doc strings. --- parametric_plasma_source/plasma_source.cpp | 16 ++++++++-------- parametric_plasma_source/plasma_source.hpp | 16 ++++++++-------- .../plasma_source_pybind.cpp | 15 +++++++++++++-- parametric_plasma_source/source_sampling.cpp | 18 +++++++++--------- 4 files changed, 38 insertions(+), 27 deletions(-) diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index 96d6e08..a760e9e 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -58,14 +58,14 @@ PlasmaSource::PlasmaSource(const double ion_density_ped, const double ion_densit PlasmaSource::~PlasmaSource(){} // main master sample function -void PlasmaSource::SampleSource(std::array random_numbers, - double &x, - double &y, - double &z, - double &u, - double &v, - double &w, - double &E) { +void PlasmaSource::sample_source(std::array random_numbers, + double &x, + double &y, + double &z, + double &u, + double &v, + double &w, + double &E) { double radius = 0.; int bin = 0; sample_source_radial(random_numbers[0],random_numbers[1],radius,bin); diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index eef37df..ebbed30 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -30,14 +30,14 @@ PlasmaSource(const double ion_density_ped, const double ion_density_sep, const double max_toridal_angle = 360.); // main sample fucnction -void SampleSource(std::array randoms, - double &x, - double &y, - double &z, - double &u, - double &v, - double &w, - double &E); +void sample_source(std::array randoms, + double &x, + double &y, + double &z, + double &u, + double &v, + double &w, + double &E); /* * Function to setup the plasma source in the first case. diff --git a/parametric_plasma_source/plasma_source_pybind.cpp b/parametric_plasma_source/plasma_source_pybind.cpp index dbf1cad..8dcbfb5 100644 --- a/parametric_plasma_source/plasma_source_pybind.cpp +++ b/parametric_plasma_source/plasma_source_pybind.cpp @@ -1,4 +1,5 @@ #include +#include #include "plasma_source.hpp" namespace py = pybind11; @@ -39,12 +40,22 @@ PYBIND11_MODULE(plasma_source, m) { py::arg("minor_radius")) .def("ion_temperature", &ps::PlasmaSource::ion_temperature, - "calculate the ion temperature at a specific minor radius", + "Calculate the ion temperature at a specific minor radius", py::arg("minor_radius")) .def("dt_xs", &ps::PlasmaSource::dt_xs, - "determine the value of the dt xs cross sections at a specific ion temperature", + "Determine the value of the dt xs cross sections at a specific ion temperature", py::arg("ion_temperature")) + .def("sample_source", + [](ps::PlasmaSource &source, std::array random_numbers) { + double x, y, z; + double u, v, w; + double e; + source.SampleSource(random_numbers, x, y, z, u, v, w, e); + return py::make_tuple(x, y, z, u, v, w, e); + }, + "Sample the source", + py::arg("random_numbers")) .def("to_xml", (bool (ps::PlasmaSource::*)(std::string)) &ps::PlasmaSource::to_xml, "Serialise the PlasmaSource to XML", diff --git a/parametric_plasma_source/source_sampling.cpp b/parametric_plasma_source/source_sampling.cpp index c0c198a..08d86a8 100644 --- a/parametric_plasma_source/source_sampling.cpp +++ b/parametric_plasma_source/source_sampling.cpp @@ -16,17 +16,17 @@ extern "C" openmc::Particle::Bank sample_source(uint64_t* seed, char* serializat // position std::array randoms = {openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed)}; + openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed)}; double u,v,w,E; - source.SampleSource(randoms,particle.r.x,particle.r.y,particle.r.z, - u,v,w,E); + source.sample_source(randoms,particle.r.x,particle.r.y,particle.r.z, + u,v,w,E); particle.r.x *= 100.; particle.r.y *= 100.; From 6289493869ce71fd8b48d3db47672ff13b6fbd32 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Fri, 31 Jul 2020 15:13:31 +0100 Subject: [PATCH 15/40] Fix typo --- parametric_plasma_source/plasma_source_pybind.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/parametric_plasma_source/plasma_source_pybind.cpp b/parametric_plasma_source/plasma_source_pybind.cpp index 8dcbfb5..018c39b 100644 --- a/parametric_plasma_source/plasma_source_pybind.cpp +++ b/parametric_plasma_source/plasma_source_pybind.cpp @@ -51,7 +51,7 @@ PYBIND11_MODULE(plasma_source, m) { double x, y, z; double u, v, w; double e; - source.SampleSource(random_numbers, x, y, z, u, v, w, e); + source.sample_source(random_numbers, x, y, z, u, v, w, e); return py::make_tuple(x, y, z, u, v, w, e); }, "Sample the source", From c8d61421943f9319c7ac577ad0ac3fae4157b6de Mon Sep 17 00:00:00 2001 From: Dan Short Date: Fri, 31 Jul 2020 15:32:34 +0100 Subject: [PATCH 16/40] Strip out build files that are no longer needed As the functionality will be used via serialisation, this means that we won't need to recompile the C++ library. --- parametric_plasma_source/Makefile | 52 --- parametric_plasma_source/__init__.py | 2 - parametric_plasma_source/build_python.py | 56 --- parametric_plasma_source/compile.sh | 5 - parametric_plasma_source/plasma.py | 380 --------------- parametric_plasma_source/source.py | 567 ----------------------- setup.py | 4 - 7 files changed, 1066 deletions(-) delete mode 100644 parametric_plasma_source/Makefile delete mode 100644 parametric_plasma_source/build_python.py delete mode 100644 parametric_plasma_source/compile.sh delete mode 100644 parametric_plasma_source/plasma.py delete mode 100644 parametric_plasma_source/source.py diff --git a/parametric_plasma_source/Makefile b/parametric_plasma_source/Makefile deleted file mode 100644 index a060d7b..0000000 --- a/parametric_plasma_source/Makefile +++ /dev/null @@ -1,52 +0,0 @@ - -# Makefile to build dynamic sources for OpenMC -# this assumes that your source sampling filename is -# source_sampling.cpp -# -# you can add fortran, c and cpp dependencies to this source -# adding them in FC_DEPS, C_DEPS, and CXX_DEPS accordingly - -ifeq ($(FC), f77) -FC = gfortran -endif - -default: all - -ALL = source_sampling -# add your fortran depencies here -FC_DEPS = -# add your c dependencies here -C_DEPS = -# add your cpp dependencies here -CPP_DEPS = plasma_source.cpp -DEPS = $(FC_DEPS) $(C_DEPS) $(CPP_DEPS) -OPT_LEVEL = -O3 -FFLAGS = $(OPT_LEVEL) -fPIC -C_FLAGS = -fPIC -CXX_FLAGS = -fPIC - - #this directory will need changing if your openmc path is different -OPENMC_DIR = /opt/openmc -OPENMC_INC_DIR = $(OPENMC_DIR)/include -OPENMC_LIB_DIR = $(OPENMC_DIR)/lib -# setting the so name is important -LINK_FLAGS = $(OPT_LEVEL) -Wall -Wl,-soname,source_sampling.so -# setting shared is important -LINK_FLAGS += -L$(OPENMC_LIB_DIR) -lopenmc -lgfortran -fPIC -shared -OPENMC_INCLUDES = -I$(OPENMC_INC_DIR) -I$(OPENMC_DIR)/vendor/pugixml - -all: $(ALL) - -source_sampling: $(DEPS) - $(CXX) source_sampling.cpp $(DEPS) $(OPENMC_INCLUDES) $(LINK_FLAGS) -o $@.so -# make any fortran objects -%.o : %.F90 - $(FC) -c $(FFLAGS) $*.F90 -o $@ -# make any c objects -%.o : %.c - $(CC) -c $(FFLAGS) $*.c -o $@ -#make any cpp objects -%.o : %.cpp - $(CXX) -c $(FFLAGS) $*.cpp -o $@ -clean: - rm -rf *.o *.mod diff --git a/parametric_plasma_source/__init__.py b/parametric_plasma_source/__init__.py index 351009c..e69de29 100644 --- a/parametric_plasma_source/__init__.py +++ b/parametric_plasma_source/__init__.py @@ -1,2 +0,0 @@ -from .plasma import Plasma -from .plasma_source import PlasmaSource diff --git a/parametric_plasma_source/build_python.py b/parametric_plasma_source/build_python.py deleted file mode 100644 index 6ee7e39..0000000 --- a/parametric_plasma_source/build_python.py +++ /dev/null @@ -1,56 +0,0 @@ -import os - -plasma_source_cpp = "plasma_source.cpp" -plasma_source_hpp = "plasma_source.hpp" -source_sampling_cpp = "source_sampling.cpp" -make_file = "Makefile" - -py_file = "source.py" - -def write_disclaimer(py_file_path): - disclaimer = "\"\"\"\n" - disclaimer += "Auto-generated file. To edit, run `python build_python.py`.\n" - disclaimer += "This will need to be run whenever a new C++ version is made available.\n" - disclaimer += "\"\"\"\n\n" - with open(py_file_path, "a+") as f_py: - f_py.write(disclaimer) - -def generate_variable(cpp_file_path, variable_name, py_file_path, end_str="\n\n"): - with open(cpp_file_path, "r") as f_cpp: - lines = f_cpp.readlines() - py_lines = variable_name - py_lines += " = (\n\"\"\"" - py_lines += "".join(lines) - py_lines += "\"\"\"\n)" - py_lines += end_str - with open(py_file_path, "a+") as f_py: - f_py.write(py_lines) - -def build(source_dir="."): - output_path = os.path.join(source_dir, py_file) - - if os.path.exists(output_path): - os.remove(output_path) - - write_disclaimer(output_path) - generate_variable( - os.path.join(source_dir, source_sampling_cpp), - "source_sampling_cpp", - output_path - ) - generate_variable( - os.path.join(source_dir, plasma_source_cpp), - "plasma_source_cpp", - output_path - ) - generate_variable( - os.path.join(source_dir, plasma_source_hpp), - "plasma_source_hpp", - output_path - ) - generate_variable( - os.path.join(source_dir, make_file), - "make_file", - output_path, - "\n" - ) diff --git a/parametric_plasma_source/compile.sh b/parametric_plasma_source/compile.sh deleted file mode 100644 index cd74266..0000000 --- a/parametric_plasma_source/compile.sh +++ /dev/null @@ -1,5 +0,0 @@ -#!/bin/bash - -make clean - -make \ No newline at end of file diff --git a/parametric_plasma_source/plasma.py b/parametric_plasma_source/plasma.py deleted file mode 100644 index 7cab511..0000000 --- a/parametric_plasma_source/plasma.py +++ /dev/null @@ -1,380 +0,0 @@ - -import os -from pathlib import Path -import tempfile -import shutil -from .source import source_sampling_cpp, plasma_source_cpp, plasma_source_hpp, make_file - - -class Plasma(): - - def __init__(self, - elongation=2., - major_radius=450, - minor_radius=150, - single_null=True, - triangularity=0.55, - ion_density_pedistal=1.09e20, - ion_density_seperatrix=3e19, - ion_density_origin=1.09e20, - ion_temperature_pedistal=6.09, - ion_temperature_seperatrix=0.1, - ion_temperature_origin=45.9, - pedistal_radius=120, # 0.8 * minor_radius - ion_density_peaking_factor=1, - ion_temperature_peaking_factor=8.06, - ion_temperature_beta=6.0, - shafranov_shift=0.0, - number_of_bins=100, - plasma_type=1, - openmc_install_directory = '/opt/openmc/' - ): - - # properties needed for plasma shapes - self.elongation = elongation - self.major_radius = major_radius - self.minor_radius = minor_radius - self.single_null = single_null - self.triangularity = triangularity - self.ion_density_pedistal = ion_density_pedistal # ions per m^3 - self.ion_density_seperatrix = ion_density_seperatrix - self.ion_density_origin = ion_density_origin - self.ion_temperature_pedistal = ion_temperature_pedistal - self.ion_temperature_seperatrix = ion_temperature_seperatrix - self.ion_temperature_origin = ion_temperature_origin - self.pedistal_radius = pedistal_radius # pedistal major rad - self.ion_density_peaking_factor = ion_density_peaking_factor - self.ion_temperature_peaking_factor = ion_temperature_peaking_factor - self.ion_temperature_beta = ion_temperature_beta - self.shafranov_shift = shafranov_shift - self.number_of_bins = number_of_bins - self.plasma_type = plasma_type # 0 = L mode anything else H/A mode - self.openmc_install_directory = openmc_install_directory - - self.plasma_source_cpp_file = plasma_source_cpp - self.plasma_source_hpp_file = plasma_source_hpp - self.source_sampling_cpp_file = source_sampling_cpp - self.plasma_make_file = make_file - - @property - def openmc_install_directory(self): - return self._openmc_install_directory - - @openmc_install_directory.setter - def openmc_install_directory(self, value): - if Path(value).exists() == False: - raise ValueError('openmc_install_directory is out of range') - else: - self._openmc_install_directory = value - - @property - def plasma_type(self): - return self._plasma_type - - @plasma_type.setter - def plasma_type(self, plasma_type): - if plasma_type < 0: - raise ValueError('plasma_type is out of range') - else: - self._plasma_type = plasma_type - - @property - def number_of_bins(self): - return self._number_of_bins - - @number_of_bins.setter - def number_of_bins(self, number_of_bins): - if number_of_bins < 0: - raise ValueError('number_of_bins is out of range') - else: - self._number_of_bins = number_of_bins - - @property - def shafranov_shift(self): - return self._shafranov_shift - - @shafranov_shift.setter - def shafranov_shift(self, shafranov_shift): - if shafranov_shift < 0: - raise ValueError('shafranov_shift is out of range') - else: - self._shafranov_shift = shafranov_shift - - @property - def ion_temperature_peaking_factor(self): - return self._ion_temperature_peaking_factor - - @ion_temperature_peaking_factor.setter - def ion_temperature_peaking_factor(self, ion_temperature_peaking_factor): - if ion_temperature_peaking_factor < 0: - raise ValueError('ion_temperature_peaking_factor is out of range') - else: - self._ion_temperature_peaking_factor = ion_temperature_peaking_factor - - @property - def ion_temperature_beta(self): - return self._ion_temperature_beta - - @ion_temperature_beta.setter - def ion_temperature_beta(self, ion_temperature_beta): - if ion_temperature_beta < 0: - raise ValueError('ion_temperature_beta is out of range') - else: - self._ion_temperature_beta = ion_temperature_beta - - @property - def ion_density_peaking_factor(self): - return self._ion_density_peaking_factor - - @ion_density_peaking_factor.setter - def ion_density_peaking_factor(self, ion_density_peaking_factor): - if ion_density_peaking_factor < 0: - raise ValueError('ion_density_peaking_factor is out of range') - else: - self._ion_density_peaking_factor = ion_density_peaking_factor - - @property - def pedistal_radius(self): - return self._pedistal_radius - - @pedistal_radius.setter - def pedistal_radius(self, pedistal_radius): - if pedistal_radius < 0: - raise ValueError('pedistal_radius is out of range') - else: - self._pedistal_radius = pedistal_radius - - @property - def ion_temperature_origin(self): - return self._ion_temperature_origin - - @ion_temperature_origin.setter - def ion_temperature_origin(self, ion_temperature_origin): - if ion_temperature_origin < 0: - raise ValueError('ion_temperature_origin is out of range') - else: - self._ion_temperature_origin = ion_temperature_origin - - @property - def ion_temperature_seperatrix(self): - return self._ion_temperature_seperatrix - - @ion_temperature_seperatrix.setter - def ion_temperature_seperatrix(self, ion_temperature_seperatrix): - if ion_temperature_seperatrix < 0: - raise ValueError('ion_temperature_seperatrix is out of range') - else: - self._ion_temperature_seperatrix = ion_temperature_seperatrix - - @property - def ion_temperature_pedistal(self): - return self._ion_temperature_pedistal - - @ion_temperature_pedistal.setter - def ion_temperature_pedistal(self, ion_temperature_pedistal): - if ion_temperature_pedistal < 0: - raise ValueError('ion_temperature_pedistal is out of range') - else: - self._ion_temperature_pedistal = ion_temperature_pedistal - - @property - def ion_density_origin(self): - return self._ion_density_origin - - @ion_density_origin.setter - def ion_density_origin(self, ion_density_origin): - if ion_density_origin < 0: - raise ValueError('ion_density_origin is out of range') - else: - self._ion_density_origin = ion_density_origin - - @property - def ion_density_seperatrix(self): - return self._ion_density_seperatrix - - @ion_density_seperatrix.setter - def ion_density_seperatrix(self, ion_density_seperatrix): - if ion_density_seperatrix < 0: - raise ValueError('ion_density_seperatrix is out of range') - else: - self._ion_density_seperatrix = ion_density_seperatrix - - @property - def triangularity(self): - return self._triangularity - - @triangularity.setter - def triangularity(self, triangularity): - if triangularity > 2000 or triangularity < -2000: - raise ValueError('triangularity is out of range') - else: - self._triangularity = triangularity - - @property - def single_null(self): - return self._single_null - - @single_null.setter - def single_null(self, single_null): - if type(single_null) != bool : - raise ValueError('single_null must be True or False') - else: - self._single_null = single_null - - @property - def minor_radius(self): - return self._minor_radius - - @minor_radius.setter - def minor_radius(self, minor_radius): - self._minor_radius = minor_radius - - @property - def major_radius(self): - return self._major_radius - - @major_radius.setter - def major_radius(self, major_radius): - if major_radius > 2000 or major_radius < 1: - raise ValueError('major_radius is out of range') - else: - self._major_radius = major_radius - - @property - def elongation(self): - return self._elongation - - @elongation.setter - def elongation(self, elongation): - if elongation > 4 or elongation < 0: - raise ValueError('elongation is out of range') - else: - self._elongation = elongation - - @property - def ion_density_pedistal(self): - return self._ion_density_pedistal - - @ion_density_pedistal.setter - def ion_density_pedistal(self, ion_density_pedistal): - if ion_density_pedistal > 10e22 or ion_density_pedistal < 1e4: - raise ValueError('ion_density_pedistal is out of range') - else: - self._ion_density_pedistal = ion_density_pedistal - - def replace_variable_value( - self, input_strings, variable_name, new_value, is_cpp=True - ): - """Replaces the value assigned to a variable with the provided value - - param: input_strings: The list of strings to search for the variable. - type input_strings: List[str] - - param: variable_name: The name of the variable to search for. - type variable_name: str - - param: new_value: The value to set the variable to. - type new_value: Union[str, float, int] - - param is_cpp: Whether the variable is a C++ variable. - Optional, by default True. - type is_cpp: bool - - raises ValueError: variable_name was not found in input_strings. - """ - if is_cpp: - strs_to_find = ("const", variable_name, "=", ";") - else: - strs_to_find = (variable_name, "=") - for idx, string in enumerate(input_strings): - if all(str_to_find in string for str_to_find in strs_to_find): - equals_idx = string.find("=") - if equals_idx >= 0: - if is_cpp: - new_value = f"{new_value};" - input_strings[idx] = string.replace( - input_strings[idx][equals_idx:], f" = {new_value}" - ) - break - else: - file_content = "\n".join(input_strings) - raise ValueError( - f"{variable_name} string not found in {file_content}" - ) - - def export_plasma_source(self, output_filename): - """Writes and compiles custom plasma source for the reactor - :param output_folder: the output folder where the .so complied plasma source will be created - :type output_folder: str - ... - :return: filename of the compiled source - :rtype: str - """ - - if self.openmc_install_directory is None: - raise ValueError('directory must be set to create .so file') - - temp_folder = Path(tempfile.mkdtemp()) - - Path(output_filename).parent.mkdir(parents=True, exist_ok=True) - - plasma_make_file_lines = self.plasma_make_file.split("\n") - self.replace_variable_value( - plasma_make_file_lines, "OPENMC_DIR", self.openmc_install_directory, False - ) - - with open(temp_folder/'Makefile', "w") as text_file: - text_file.write("\n".join(plasma_make_file_lines)) - - with open(temp_folder/'plasma_source.cpp', "w") as text_file: - text_file.write(self.plasma_source_cpp_file) - - with open(temp_folder/'plasma_source.hpp', "w") as text_file: - text_file.write(self.plasma_source_hpp_file) - - plasma_variables = { - "ion_density_pedistal": self.ion_density_pedistal, - "ion_density_seperatrix": self.ion_density_seperatrix, - "ion_density_origin": self.ion_density_origin, - "ion_temperature_pedistal": self.ion_temperature_pedistal, - "ion_temperature_seperatrix": self.ion_temperature_seperatrix, - "ion_temperature_origin": self.ion_temperature_origin, - "pedistal_radius": self.pedistal_radius, - "ion_density_peaking_factor": self.ion_density_peaking_factor, - "ion_temperature_peaking_factor": self.ion_temperature_peaking_factor, - "ion_temperature_beta": self.ion_temperature_beta, - "minor_radius": self.minor_radius, - "major_radius": self.major_radius, - "elongation": self.elongation, - "triangularity": self.triangularity, - "shafranov_shift": self.shafranov_shift, - "number_of_bins": self.number_of_bins, - "plasma_type": self.plasma_type, - } - - source_sampling_cpp_lines = self.source_sampling_cpp_file.split("\n") - [ - self.replace_variable_value(source_sampling_cpp_lines, *variable) - for variable in plasma_variables.items() - ] - - with open(temp_folder/'source_sampling.cpp', "w") as text_file: - text_file.write("\n".join(source_sampling_cpp_lines)) - - cwd = os.getcwd() - os.chdir(Path(temp_folder)) - - os.system('make clean') - os.system('make') - - os.chdir(cwd) - shutil.move(temp_folder/'source_sampling.so', output_filename) - print('parametric plasma source compiled and saved to ', output_filename) - shutil.rmtree(temp_folder) - - return output_filename - - -if __name__ == "__main__": - p = Plasma() - p.export_plasma_source("source_sampling.so") diff --git a/parametric_plasma_source/source.py b/parametric_plasma_source/source.py deleted file mode 100644 index f744dd4..0000000 --- a/parametric_plasma_source/source.py +++ /dev/null @@ -1,567 +0,0 @@ -""" -Auto-generated file. To edit, run `python build_python.py`. -This will need to be run whenever a new C++ version is made available. -""" - -source_sampling_cpp = ( -"""#include -#include "openmc/random_lcg.h" -#include "openmc/source.h" -#include "openmc/particle.h" -#include "plasma_source.hpp" - - -// Spherical tokamak SOURCE -// units are in SI units -const double ion_density_pedistal = 1.09e+20; // ions per m^3 -const double ion_density_seperatrix = 3e+19; -const double ion_density_origin = 1.09e+20; -const double ion_temperature_pedistal = 6.09; -const double ion_temperature_seperatrix = 0.1; -const double ion_temperature_origin = 45.9; -const double ion_density_peaking_factor = 1; -const double ion_temperature_peaking_factor = 8.06; // check alpha or beta value from paper -const double ion_temperature_beta = 6.0; -const double minor_radius = 1.56; // metres -const double major_radius = 2.5; // metres -const double pedistal_radius = 0.8 * minor_radius; // pedistal minor rad in metres -const double elongation = 2.0; -const double triangularity = 0.55; -const double shafranov_shift = 0.0; //metres -const std::string name = "parametric_plasma_source"; -const int number_of_bins = 100; -const int plasma_type = 1; // 1 is default; //0 = L mode anything else H/A mode - - -plasma_source::PlasmaSource source = plasma_source::PlasmaSource(ion_density_pedistal, - ion_density_seperatrix, - ion_density_origin, - ion_temperature_pedistal, - ion_temperature_seperatrix, - ion_temperature_origin, - pedistal_radius, - ion_density_peaking_factor, - ion_temperature_peaking_factor, - ion_temperature_beta, - minor_radius, - major_radius, - elongation, - triangularity, - shafranov_shift, - name, - plasma_type, - number_of_bins, - 0.0, - 360.0); - -// you must have external C linkage here otherwise -// dlopen will not find the file -extern "C" openmc::Particle::Bank sample_source(uint64_t* seed) { - openmc::Particle::Bank particle; - // wgt - particle.particle = openmc::Particle::Type::neutron; - particle.wgt = 1.0; - // position - - std::array randoms = {openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed)}; - - double u,v,w,E; - source.SampleSource(randoms,particle.r.x,particle.r.y,particle.r.z, - u,v,w,E); - - particle.r.x *= 100.; - particle.r.y *= 100.; - particle.r.z *= 100.; - - // particle.E = 14.08e6; - particle.E = E*1e6; // convert from MeV -> eV - - particle.u = {u, - v, - w}; - - - particle.delayed_group = 0; - return particle; -} - - -""" -) - -plasma_source_cpp = ( -"""#include -#include -#include -#include "plasma_source.hpp" -#include -#include "openmc/random_lcg.h" - -#define RANDOM openmc::prn() - -namespace plasma_source { - -// default constructor -PlasmaSource::PlasmaSource() {} - -// large constructor -PlasmaSource::PlasmaSource(const double ion_density_ped, const double ion_density_sep, - const double ion_density_origin, const double ion_temp_ped, - const double ion_temp_sep, const double ion_temp_origin, - const double pedistal_rad, const double ion_density_peak, - const double ion_temp_peak, const double ion_temp_beta, - const double minor_radius, const double major_radius, - const double elongation, const double triangularity, - const double shafranov, const std::string plasma_type, - const int plasma_id, const int number_of_bins, - const double min_toroidal_angle, - const double max_toroidal_angle ) { - - // set params - ionDensityPedistal = ion_density_ped; - ionDensitySeperatrix = ion_density_sep; - ionDensityOrigin = ion_density_origin; - ionTemperaturePedistal = ion_temp_ped; - ionTemperatureSeperatrix = ion_temp_sep; - ionTemperatureOrigin = ion_temp_origin; - pedistalRadius = pedistal_rad; - ionDensityPeaking = ion_density_peak; - ionTemperaturePeaking = ion_temp_peak; - ionTemperatureBeta = ion_temp_beta; - minorRadius = minor_radius; - majorRadius = major_radius; - this->elongation = elongation; - this->triangularity = triangularity; - this->shafranov = shafranov; - plasmaType = plasma_type; - plasmaId = plasma_id; - numberOfBins = number_of_bins; - minToroidalAngle = min_toroidal_angle/180.*M_PI; - maxToroidalAngle = max_toroidal_angle/180.*M_PI; - - setup_plasma_source(); -} - -// destructor -PlasmaSource::~PlasmaSource(){} - -// main master sample function -void PlasmaSource::SampleSource(std::array random_numbers, - double &x, - double &y, - double &z, - double &u, - double &v, - double &w, - double &E) { - double radius = 0.; - int bin = 0; - sample_source_radial(random_numbers[0],random_numbers[1],radius,bin); - double r = 0.; - convert_rad_to_rz(radius,random_numbers[2],r,z); - convert_r_to_xy(r,random_numbers[3],x,y); - sample_energy(bin,random_numbers[4],random_numbers[5],E); - isotropic_direction(random_numbers[6],random_numbers[7], - u,v,w); - -} - -/* - * sample the pdf src_profile, to generate the sampled minor radius - */ -void PlasmaSource::sample_source_radial(double rn_store1, double rn_store2, - double &sampled_radius, int &sampled_bin) { - - for ( int i = 0 ; i < numberOfBins ; i++ ) { - if ( rn_store1 <= source_profile[i] ) { - if ( i > 0 ) { - sampled_radius = (float(i-1)*(binWidth)) + (binWidth*(rn_store2)); - sampled_bin = i; - return; - } else { - sampled_radius = binWidth*(rn_store2); - sampled_bin = i; - return; - } - } - } - - std::cerr << "error" << std::endl; - std::cerr << "Sample position greater than plasma radius" << std::endl; - exit(1); - return; -} - -/* - * sample the energy of the neutrons, updates energy neutron in mev - */ -void PlasmaSource::sample_energy(const int bin_number, double random_number1, double random_number2, - double &energy_neutron) { - // generate the normally distributed number - const double twopi = 6.28318530718; - double sample1 = std::sqrt(-2.0*std::log(random_number1)); - double sample2 = cos(twopi*(random_number2)); - energy_neutron = (5.59/2.35)*(ion_kt[bin_number])*sample1*sample2; - energy_neutron += 14.08; - // test energy limit - // if (energy_neutron < 15.5){energy_neutron = 15.5} else {} - return; -} - -/* - * convert the sampled radius to an rz coordinate by using plasma parameters - */ -void PlasmaSource::convert_rad_to_rz( const double minor_sampled, - const double rn_store, - double &radius, double &height) -{ - const double twopi = 6.28318530718; - - double alpha = twopi*(rn_store); - - double shift = shafranov*(1.0-std::pow(minor_sampled/(minorRadius),2)); - - radius = majorRadius + minor_sampled*cos(alpha+(triangularity*sin(alpha))) + shift; - height = elongation*minor_sampled*sin(alpha); - - - return; -} - - -/* - * convert rz_to_xyz - */ -void PlasmaSource::convert_r_to_xy(const double r, const double rn_store, - double &x, double &y) - -{ - double toroidal_extent = maxToroidalAngle - minToroidalAngle; - double toroidal_angle = toroidal_extent*rn_store + minToroidalAngle; - x = r*sin(toroidal_angle); - y = r*cos(toroidal_angle); - return; -} - -/* - * sets up the cumulatitive probability profile - * on the basis of the ion temp and ion density - * this portion is deterministic - */ -void PlasmaSource::setup_plasma_source() -{ - double ion_d; // ion density - double ion_t; // ion temp - - std::vector src_strength; // the source strength, n/m3 - double r; - - binWidth = minorRadius/float(numberOfBins); - double total = 0.0; // total source strength - - for (int i = 0 ; i < numberOfBins ; i++) { - r = binWidth * float(i); - ion_d = ion_density(r); - ion_t = ion_temperature(r); - src_strength.push_back(std::pow(ion_d,2)*dt_xs(ion_t)); - ion_kt.push_back(sqrt(ion_t/1000.0)); // convert to sqrt(MeV) - total += src_strength[i]; - } - - // normalise the source profile - double sum = 0 ; - for ( int i = 0 ; i < numberOfBins ; i++) { - sum += src_strength[i]; - source_profile.push_back(sum/total); - } - return; -} - -/* - * function that returns the ion density given the - * given the critical plasma parameters - */ -double PlasmaSource::ion_density(const double sample_radius) -{ - double ion_dens = 0.0; - - if( plasmaId == 0 ) { - ion_dens = ionDensityOrigin* - (1.0-std::pow(sample_radius/minorRadius,2)); - } else { - if(sample_radius <= pedistalRadius) { - ion_dens += ionDensityPedistal; - double product; - product = 1.0-std::pow(sample_radius/pedistalRadius,2); - product = std::pow(product,ionDensityPeaking); - ion_dens += (ionDensityOrigin-ionDensityPedistal)* - (product); - } else { - ion_dens += ionDensitySeperatrix; - double product; - product = ionDensityPedistal - ionDensitySeperatrix; - ion_dens += product*(minorRadius-sample_radius)/(minorRadius-pedistalRadius); - } - } - - return ion_dens; -} - -/* - * function that returns the ion density given the - * given the critical plasma parameters - */ -double PlasmaSource::ion_temperature(const double sample_radius) -{ - double ion_temp = 0.0; - - if( plasmaId == 0 ) { - ion_temp = ionTemperatureOrigin* - (1.0-std::pow(sample_radius/minorRadius, - ionTemperaturePeaking)); - } else { - if(sample_radius <= pedistalRadius) { - ion_temp += ionTemperaturePedistal; - double product; - product = 1.0-std::pow(sample_radius/pedistalRadius,ionTemperatureBeta); - product = std::pow(product,ionTemperaturePeaking); - ion_temp += (ionTemperatureOrigin- - ionTemperaturePedistal)*(product); - } else { - ion_temp += ionTemperatureSeperatrix; - double product; - product = ionTemperaturePedistal - ionTemperatureSeperatrix; - ion_temp += product*(minorRadius-sample_radius)/(minorRadius-pedistalRadius); - } - } - - return ion_temp; -} - -/* - * returns the dt cross section for a given ion temp - */ -double PlasmaSource::dt_xs(double ion_temp) -{ - double dt; - double c[7]={2.5663271e-18,19.983026,2.5077133e-2, - 2.5773408e-3,6.1880463e-5,6.6024089e-2, - 8.1215505e-3}; - - double u = 1.0-ion_temp*(c[2]+ion_temp*(c[3]-c[4]*ion_temp)) - /(1.0+ion_temp*(c[5]+c[6]*ion_temp)); - - dt = c[0]/(std::pow(u,5./6.)*std::pow(ion_temp,2./3.0)); - dt *= exp(-1.*c[1]*std::pow(u/ion_temp,1./3.)); - - return dt; -} - -/* - * returns the dt cross section for a given ion temp - */ -void PlasmaSource::isotropic_direction(const double random1, - const double random2, - double &u, double &v, - double &w) { - double t = 2*M_PI*random1; - double p = acos(1. - 2.*random2); - - u = sin(p)*cos(t); - v = sin(p)*sin(t); - w = cos(p); - - return; -} - -} // end of namespace -""" -) - -plasma_source_hpp = ( -"""#include -#include -namespace plasma_source { - -struct xs_params { - double c[7]; -}; - -class PlasmaSource { -public: -// constructor -PlasmaSource(); -// destructor -~PlasmaSource(); -// large constructor -PlasmaSource(const double ion_density_ped, const double ion_density_sep, - const double ion_density_origin, const double ion_temp_ped, - const double ion_temp_sep, const double ion_temp_origin, - const double pedistal_rad, const double ion_density_peak, - const double ion_temp_peak, const double ion_temp_beta, - const double minor_radius, const double major_radius, - const double elongation, const double triangularity, - const double shafranov, const std::string plasma_type, - const int plasma_id, const int number_of_bins, - const double min_toroidal_angle = 0.0, - const double max_toridal_angle = 360.); - -// main sample fucnction -void SampleSource(std::array randoms, - double &x, - double &y, - double &z, - double &u, - double &v, - double &w, - double &E); - -/* - * Function to setup the plasma source in the first case. - */ -void setup_plasma_source(); - -/* - * function to calculate the ion density at a specific minor - * radius - */ -double ion_density(const double sample_radius); - -/* - * function to calculate the ion temperature at a specific minor - * radius - */ -double ion_temperature(const double sample_radius); - -/* - * function to determine the value of the dt xs cross sections at - * a specific ion temp - */ -double dt_xs(double ion_temp); - -/* - * sample the source, returns the minor radius sampled - * expects new rn_store every call - */ -void sample_source_radial(double rn_store1, double rn_store2, - double &sampled_radius, - int &sampled_bin); - -/* - * sample the neutron energy in MeV - */ -void sample_energy(const int bin_number, double random_number1, double random_number2, - double &energy_neutron); - -/* - * take the sampled minor radius and convert to cylindrical coordinates - */ -void convert_rad_to_rz(const double minor_sampled, - const double rn_store, - double &radius, - double &height); - -/* - * convert partial cylindrical coords to xyz - */ -void convert_r_to_xy(const double r, const double rn_store, - double &x, double &y); -/* - * get an isotropically direction vector - */ -void isotropic_direction(const double random1, const double random2, - double &u, double &v, double &w); - -private: - std::vector source_profile; - std::vector ion_kt; - - double ionDensityPedistal; - double ionDensitySeperatrix; - double ionDensityOrigin; - double ionTemperaturePedistal; - double ionTemperatureSeperatrix; - double ionTemperatureOrigin; - double pedistalRadius; - double ionDensityPeaking; - double ionTemperaturePeaking; - double ionTemperatureBeta; - double minorRadius; - double majorRadius; - double elongation; - double triangularity; - double shafranov; - double minToroidalAngle; - double maxToroidalAngle; - - std::string plasmaType; - int plasmaId; - double binWidth; - int numberOfBins; - -}; -}// end of namespace""" -) - -make_file = ( -""" -# Makefile to build dynamic sources for OpenMC -# this assumes that your source sampling filename is -# source_sampling.cpp -# -# you can add fortran, c and cpp dependencies to this source -# adding them in FC_DEPS, C_DEPS, and CXX_DEPS accordingly - -ifeq ($(FC), f77) -FC = gfortran -endif - -default: all - -ALL = source_sampling -# add your fortran depencies here -FC_DEPS = -# add your c dependencies here -C_DEPS = -# add your cpp dependencies here -CPP_DEPS = plasma_source.cpp -DEPS = $(FC_DEPS) $(C_DEPS) $(CPP_DEPS) -OPT_LEVEL = -O3 -FFLAGS = $(OPT_LEVEL) -fPIC -C_FLAGS = -fPIC -CXX_FLAGS = -fPIC - - #this directory will need changing if your openmc path is different -OPENMC_DIR = /opt/openmc -OPENMC_INC_DIR = $(OPENMC_DIR)/include -OPENMC_LIB_DIR = $(OPENMC_DIR)/lib -# setting the so name is important -LINK_FLAGS = $(OPT_LEVEL) -Wall -Wl,-soname,source_sampling.so -# setting shared is important -LINK_FLAGS += -L$(OPENMC_LIB_DIR) -lopenmc -lgfortran -fPIC -shared -OPENMC_INCLUDES = -I$(OPENMC_INC_DIR) -I$(OPENMC_DIR)/vendor/pugixml - -all: $(ALL) - -source_sampling: $(DEPS) - $(CXX) source_sampling.cpp $(DEPS) $(OPENMC_INCLUDES) $(LINK_FLAGS) -o $@.so -# make any fortran objects -%.o : %.F90 - $(FC) -c $(FFLAGS) $*.F90 -o $@ -# make any c objects -%.o : %.c - $(CC) -c $(FFLAGS) $*.c -o $@ -#make any cpp objects -%.o : %.cpp - $(CXX) -c $(FFLAGS) $*.cpp -o $@ -clean: - rm -rf *.o *.mod -""" -) diff --git a/setup.py b/setup.py index 698cf89..01ec90f 100644 --- a/setup.py +++ b/setup.py @@ -5,8 +5,6 @@ from setuptools import setup, find_packages, Extension from setuptools.command.build_ext import build_ext -from parametric_plasma_source.build_python import build as build_python - class CMakeExtention(Extension): def __init__(self, name, sourcedir=""): Extension.__init__(self, name, sources=[]) @@ -65,8 +63,6 @@ def build_extension(self, ext): with open("README.md", "r") as fh: long_description = fh.read() -build_python(source_dir="parametric_plasma_source") - setup( name="parametric_plasma_source", version="0.0.6", From b3ab25a845c08fb5ee3508e6d6109e14cbf1be76 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Thu, 6 Aug 2020 12:59:04 +0000 Subject: [PATCH 17/40] Update serialization methods Aligns with recent version of OpenMC: * Serialize to key-value pair string; * Expose deserialization via openmc_create_source * Remove XML handling * Remove pugixml as a submodule Some formatting updates --- .gitmodules | 3 - CMakeLists.txt | 13 +- data/example_source.txt | 1 + data/example_source.xml | 22 -- parametric_plasma_source/plasma_source.cpp | 171 +++++------- parametric_plasma_source/plasma_source.hpp | 249 +++++++++--------- .../plasma_source_pybind.cpp | 56 ++-- parametric_plasma_source/source_sampling.cpp | 94 ++++--- parametric_plasma_source/xml_helper.cpp | 47 ---- parametric_plasma_source/xml_helper.hpp | 20 -- pugixml | 1 - 11 files changed, 279 insertions(+), 398 deletions(-) create mode 100644 data/example_source.txt delete mode 100644 data/example_source.xml delete mode 100644 parametric_plasma_source/xml_helper.cpp delete mode 100644 parametric_plasma_source/xml_helper.hpp delete mode 160000 pugixml diff --git a/.gitmodules b/.gitmodules index 02ee8fa..2726187 100644 --- a/.gitmodules +++ b/.gitmodules @@ -1,6 +1,3 @@ [submodule "pybind11"] path = pybind11 url = https://github.com/pybind/pybind11 -[submodule "pugixml"] - path = pugixml - url = https://github.com/zeux/pugixml.git diff --git a/CMakeLists.txt b/CMakeLists.txt index f6d414c..8a1093a 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -38,17 +38,11 @@ if(GIT_FOUND AND EXISTS "${CMAKE_CURRENT_SOURCE_DIR}/.git") endif() # Check to see if submodules exist (by checking one) -if(NOT EXISTS "${CMAKE_CURRENT_SOURCE_DIR}/pugixml/CMakeLists.txt") +if(NOT EXISTS "${CMAKE_CURRENT_SOURCE_DIR}/pybind11/CMakeLists.txt") message(FATAL_ERROR "The git submodules were not downloaded! GIT_SUBMODULE was \ turned off or failed. Please update submodules and try again.") endif() -#=============================================================================== -# pugixml library -#=============================================================================== - -add_subdirectory(pugixml) - # Build source_sampling list(APPEND source_sampling_SOURCES ${SRC_DIR}/source_sampling.cpp @@ -62,12 +56,10 @@ find_library(OPENMC_LIB openmc HINTS ${OPENMC_LIB_DIR}) set_target_properties(source_sampling PROPERTIES PREFIX "") set_target_properties(source_sampling PROPERTIES POSITION_INDEPENDENT_CODE ON) target_include_directories(source_sampling PUBLIC ${OPENMC_INC_DIR}) -target_include_directories(source_sampling PUBLIC ${OPENMC_DIR}/vendor/pugixml) target_link_libraries(source_sampling ${OPENMC_LIB} gfortran) # Build plasma_source Python bindings list(APPEND plasma_source_pybind_SOURCES - ${SRC_DIR}/xml_helper.cpp ${SRC_DIR}/plasma_source.cpp ${SRC_DIR}/plasma_source_pybind.cpp ) @@ -75,6 +67,3 @@ list(APPEND plasma_source_pybind_SOURCES add_subdirectory(pybind11) pybind11_add_module(plasma_source ${plasma_source_pybind_SOURCES}) -target_include_directories(plasma_source PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/pugixml) -link_directories(${CMAKE_ARCHIVE_OUTPUT_DIRECTORY}) -target_link_libraries(plasma_source PRIVATE pugixml) diff --git a/data/example_source.txt b/data/example_source.txt new file mode 100644 index 0000000..ee90fb2 --- /dev/null +++ b/data/example_source.txt @@ -0,0 +1 @@ +major_radius=4.5, minor_radius=1.5, elongation=2, triangularity=0.55, shafranov_shift=0, pedestal_radius=1.2, ion_density_pedestal=1.09e+20, ion_density_separatrix=3e+19, ion_density_origin=1.09e+20, ion_density_alpha=1, ion_temperature_pedestal=6.09, ion_temperature_separatrix=0.1, ion_temperature_origin=45.9, ion_temperature_alpha=8.06, ion_temperature_beta=6, plasma_type=plasma, plasma_id=1, number_of_bins=100, minimum_toroidal_angle=0, maximum_toroidal_angle=360 \ No newline at end of file diff --git a/data/example_source.xml b/data/example_source.xml deleted file mode 100644 index 6036672..0000000 --- a/data/example_source.xml +++ /dev/null @@ -1,22 +0,0 @@ - - 1.5 - 4.5 - 2 - 0.55 - 0 - 0.8 - 1.09e+20 - 3e+19 - 1.09e+20 - 6.09 - 0.1 - 45.9 - 1 - 8.06 - 6 - plasma - 1 - 100 - 0 - 360 - diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index a760e9e..85d6188 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -1,13 +1,11 @@ -#include #include #include #include +#include #include #include "plasma_source.hpp" #include #include "openmc/random_lcg.h" -#include "xml_helper.hpp" -#include "pugixml.hpp" #define RANDOM openmc::prn() @@ -18,16 +16,16 @@ PlasmaSource::PlasmaSource() {} // large constructor PlasmaSource::PlasmaSource(const double ion_density_ped, const double ion_density_sep, - const double ion_density_origin, const double ion_temp_ped, - const double ion_temp_sep, const double ion_temp_origin, - const double pedistal_rad, const double ion_density_peak, - const double ion_temp_peak, const double ion_temp_beta, - const double minor_radius, const double major_radius, - const double elongation, const double triangularity, - const double shafranov, const std::string plasma_type, - const int plasma_id, const int number_of_bins, - const double min_toroidal_angle, - const double max_toroidal_angle ) { + const double ion_density_origin, const double ion_temp_ped, + const double ion_temp_sep, const double ion_temp_origin, + const double pedistal_rad, const double ion_density_peak, + const double ion_temp_peak, const double ion_temp_beta, + const double minor_radius, const double major_radius, + const double elongation, const double triangularity, + const double shafranov, const std::string plasma_type, + const int plasma_id, const int number_of_bins, + const double min_toroidal_angle, + const double max_toroidal_angle ) { // set params ionDensityPedistal = ion_density_ped; @@ -108,7 +106,7 @@ void PlasmaSource::sample_source_radial(double rn_store1, double rn_store2, * sample the energy of the neutrons, updates energy neutron in mev */ void PlasmaSource::sample_energy(const int bin_number, double random_number1, double random_number2, - double &energy_neutron) { + double &energy_neutron) { // generate the normally distributed number const double twopi = 6.28318530718; double sample1 = std::sqrt(-2.0*std::log(random_number1)); @@ -257,8 +255,8 @@ double PlasmaSource::dt_xs(double ion_temp) { double dt; double c[7]={2.5663271e-18,19.983026,2.5077133e-2, - 2.5773408e-3,6.1880463e-5,6.6024089e-2, - 8.1215505e-3}; + 2.5773408e-3,6.1880463e-5,6.6024089e-2, + 8.1215505e-3}; double u = 1.0-ion_temp*(c[2]+ion_temp*(c[3]-c[4]*ion_temp)) /(1.0+ion_temp*(c[5]+c[6]*ion_temp)); @@ -286,97 +284,70 @@ void PlasmaSource::isotropic_direction(const double random1, return; } -std::string PlasmaSource::to_xml() { - XMLHelper xml_helper = XMLHelper(PLASMA_SOURCE_ROOT_NAME.c_str()); - xml_helper.add_element("MajorRadius", majorRadius); - xml_helper.add_element("MinorRadius", minorRadius); - xml_helper.add_element("Elongation", elongation); - xml_helper.add_element("Triangularity", triangularity); - xml_helper.add_element("ShafranovShift", shafranov); - xml_helper.add_element("PedestalRadius", pedistalRadius); - xml_helper.add_element("IonDensityPedestal", ionDensityPedistal); - xml_helper.add_element("IonDensitySeparatrix", ionDensitySeperatrix); - xml_helper.add_element("IonDensityOrigin", ionDensityOrigin); - xml_helper.add_element("IonTemperaturePedestal", ionTemperaturePedistal); - xml_helper.add_element("IonTemperatureSeparatrix", ionTemperatureSeperatrix); - xml_helper.add_element("IonTemperatureOrigin", ionTemperatureOrigin); - xml_helper.add_element("IonDensityAlpha", ionDensityPeaking); - xml_helper.add_element("IonTemperatureAlpha", ionTemperaturePeaking); - xml_helper.add_element("IonTemperatureBeta", ionTemperatureBeta); - xml_helper.add_element("PlasmaType", plasmaType); - xml_helper.add_element("PlasmaId", plasmaId); - xml_helper.add_element("NumberOfBins", numberOfBins); - xml_helper.add_element("MinimumToroidalAngle", minToroidalAngle / M_PI * 180.0); - xml_helper.add_element("MaximumToroidalAngle", maxToroidalAngle / M_PI * 180.0); - return xml_helper.to_string(); +std::string PlasmaSource::to_string() { + std::stringstream result; + result << "major_radius=" << majorRadius << ", "; + result << "minor_radius=" << minorRadius << ", "; + result << "elongation=" << elongation << ", "; + result << "triangularity=" << triangularity << ", "; + result << "shafranov_shift=" << shafranov << ", "; + result << "pedestal_radius=" << pedistalRadius << ", "; + result << "ion_density_pedestal=" << ionDensityPedistal << ", "; + result << "ion_density_separatrix=" << ionDensitySeperatrix << ", "; + result << "ion_density_origin=" << ionDensityOrigin << ", "; + result << "ion_density_alpha=" << ionDensityPeaking << ", "; + result << "ion_temperature_pedestal=" << ionTemperaturePedistal << ", "; + result << "ion_temperature_separatrix=" << ionTemperatureSeperatrix << ", "; + result << "ion_temperature_origin=" << ionTemperatureOrigin << ", "; + result << "ion_temperature_alpha=" << ionTemperaturePeaking << ", "; + result << "ion_temperature_beta=" << ionTemperatureBeta << ", "; + result << "plasma_type=" << plasmaType << ", "; + result << "plasma_id=" << plasmaId << ", "; + result << "number_of_bins=" << numberOfBins << ", "; + result << "minimum_toroidal_angle=" << minToroidalAngle * 180.0 / M_PI << ", "; + result << "maximum_toroidal_angle=" << maxToroidalAngle * 180.0 / M_PI; + return result.str(); } -bool PlasmaSource::to_xml(std::string output_path) { - bool success = true; - std::string output = to_xml(); - std::ofstream out; - out.open(output_path); - if (out.is_open()) { - out << output; - out.close(); - } - else { - success = false; +PlasmaSource PlasmaSource::from_string(const char* parameters) { + std::unordered_map parameter_mapping; + + std::stringstream ss(parameters); + std::string parameter; + while (std::getline(ss, parameter, ',')) { + parameter.erase(0, parameter.find_first_not_of(' ')); + std::string key = parameter.substr(0, parameter.find_first_of('=')); + std::string value = parameter.substr(parameter.find_first_of('=') + 1, parameter.length()); + parameter_mapping[key] = value; } - return success; -} -PlasmaSource PlasmaSource::from_xml_doc(pugi::xml_document* xml_doc, std::string root_name) { - pugi::xml_document doc; - doc.reset(*xml_doc); - pugi::xml_node root_node = doc.root().child(root_name.c_str()); - double major_radius = root_node.child("MajorRadius").text().as_double(); - double minor_radius = root_node.child("MinorRadius").text().as_double(); - double elongation = root_node.child("Elongation").text().as_double(); - double triangularity = root_node.child("Triangularity").text().as_double(); - double shafranov_shift = root_node.child("ShafranovShift").text().as_double(); - double pedestal_radius = root_node.child("PedestalRadius").text().as_double(); - double ion_density_pedestal = root_node.child("IonDensityPedestal").text().as_double(); - double ion_density_separatrix = root_node.child("IonDensitySeparatrix").text().as_double(); - double ion_density_origin = root_node.child("IonDensityOrigin").text().as_double(); - double ion_temperature_pedestal = root_node.child("IonTemperaturePedestal").text().as_double(); - double ion_temperature_separatrix = root_node.child("IonTemperatureSeparatrix").text().as_double(); - double ion_temperature_origin = root_node.child("IonTemperatureOrigin").text().as_double(); - double ion_density_alpha = root_node.child("IonDensityAlpha").text().as_double(); - double ion_temperature_alpha = root_node.child("IonTemperatureAlpha").text().as_double(); - double ion_temperature_beta = root_node.child("IonTemperatureBeta").text().as_double(); - std::string plasma_type = root_node.child("PlasmaType").text().as_string(); - int plasma_id = root_node.child("PlasmaId").text().as_int(); - int number_of_bins = root_node.child("NumberOfBins").text().as_int(); - double min_toroidal_angle = root_node.child("MinimumToroidalAngle").text().as_double(); - double max_toroidal_angle = root_node.child("MaximumToroidalAngle").text().as_double(); + double major_radius = std::stod(parameter_mapping["major_radius"]); + double minor_radius = std::stod(parameter_mapping["minor_radius"]); + double elongation = std::stod(parameter_mapping["elongation"]); + double triangularity = std::stod(parameter_mapping["triangularity"]); + double shafranov_shift = std::stod(parameter_mapping["shafranov_shift"]); + double pedestal_radius = std::stod(parameter_mapping["pedestal_radius"]); + double ion_density_pedestal = std::stod(parameter_mapping["ion_density_pedestal"]); + double ion_density_separatrix = std::stod(parameter_mapping["ion_density_separatrix"]); + double ion_density_origin = std::stod(parameter_mapping["ion_density_origin"]); + double ion_temperature_pedestal = std::stod(parameter_mapping["ion_temperature_pedestal"]); + double ion_temperature_separatrix = std::stod(parameter_mapping["ion_temperature_separatrix"]); + double ion_temperature_origin = std::stod(parameter_mapping["ion_temperature_origin"]); + double ion_density_alpha = std::stod(parameter_mapping["ion_density_alpha"]); + double ion_temperature_alpha = std::stod(parameter_mapping["ion_temperature_alpha"]); + double ion_temperature_beta = std::stod(parameter_mapping["ion_temperature_beta"]); + std::string plasma_type = parameter_mapping["plasma_type"]; + int plasma_id = std::stoi(parameter_mapping["plasma_id"]); + int number_of_bins = std::stoi(parameter_mapping["number_of_bins"]); + double min_toroidal_angle = std::stod(parameter_mapping["minimum_toroidal_angle"]); + double max_toroidal_angle = std::stod(parameter_mapping["maximum_toroidal_angle"]); PlasmaSource plasma_source = PlasmaSource( - ion_density_pedestal, ion_density_separatrix, ion_density_origin, ion_temperature_pedestal, ion_temperature_separatrix, - ion_temperature_origin, pedestal_radius, ion_density_alpha, ion_temperature_alpha, ion_temperature_beta, minor_radius, - major_radius, elongation, triangularity, shafranov_shift, plasma_type, plasma_id, number_of_bins, min_toroidal_angle, - max_toroidal_angle); + ion_density_pedestal, ion_density_separatrix, ion_density_origin, ion_temperature_pedestal, ion_temperature_separatrix, + ion_temperature_origin, pedestal_radius, ion_density_alpha, ion_temperature_alpha, ion_temperature_beta, minor_radius, + major_radius, elongation, triangularity, shafranov_shift, plasma_type, plasma_id, number_of_bins, min_toroidal_angle, + max_toroidal_angle); return plasma_source; } -PlasmaSource PlasmaSource::from_file(std::string input_path) { - pugi::xml_document doc; - pugi::xml_parse_result result = doc.load_file(input_path.c_str()); - if (!result) { - throw std::runtime_error("Error reading plasma source from file " + input_path + ". Error = " + result.description()); - } - - return from_xml_doc(&doc, PLASMA_SOURCE_ROOT_NAME); -} - -PlasmaSource PlasmaSource::from_xml(char* xml) { - pugi::xml_document doc; - pugi::xml_parse_result result = doc.load_string(xml); - if (!result) { - throw std::runtime_error("Error reading plasma source from string " + std::string(xml) + ". Error = " + result.description()); - } - - return from_xml_doc(&doc, PLASMA_SOURCE_ROOT_NAME); -} - } // end of namespace diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index ebbed30..2931f96 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -1,6 +1,5 @@ #include #include -#include "pugixml.hpp" namespace plasma_source { @@ -11,126 +10,130 @@ struct xs_params { static const std::string PLASMA_SOURCE_ROOT_NAME = "PlasmaSource"; class PlasmaSource { -public: -// constructor -PlasmaSource(); -// destructor -~PlasmaSource(); -// large constructor -PlasmaSource(const double ion_density_ped, const double ion_density_sep, - const double ion_density_origin, const double ion_temp_ped, - const double ion_temp_sep, const double ion_temp_origin, - const double pedistal_rad, const double ion_density_peak, - const double ion_temp_peak, const double ion_temp_beta, - const double minor_radius, const double major_radius, - const double elongation, const double triangularity, - const double shafranov, const std::string plasma_type, - const int plasma_id, const int number_of_bins, - const double min_toroidal_angle = 0.0, - const double max_toridal_angle = 360.); - -// main sample fucnction -void sample_source(std::array randoms, - double &x, - double &y, - double &z, - double &u, - double &v, - double &w, - double &E); - -/* - * Function to setup the plasma source in the first case. - */ -void setup_plasma_source(); - -/* - * function to calculate the ion density at a specific minor - * radius - */ -double ion_density(const double sample_radius); - -/* - * function to calculate the ion temperature at a specific minor - * radius - */ -double ion_temperature(const double sample_radius); - -/* - * function to determine the value of the dt xs cross sections at - * a specific ion temp - */ -double dt_xs(double ion_temp); - -/* - * sample the source, returns the minor radius sampled - * expects new rn_store every call - */ -void sample_source_radial(double rn_store1, double rn_store2, - double &sampled_radius, - int &sampled_bin); - -/* - * sample the neutron energy in MeV - */ -void sample_energy(const int bin_number, double random_number1, double random_number2, - double &energy_neutron); - -/* - * take the sampled minor radius and convert to cylindrical coordinates - */ -void convert_rad_to_rz(const double minor_sampled, - const double rn_store, - double &radius, - double &height); - -/* - * convert partial cylindrical coords to xyz - */ -void convert_r_to_xy(const double r, const double rn_store, - double &x, double &y); -/* - * get an isotropically direction vector - */ -void isotropic_direction(const double random1, const double random2, - double &u, double &v, double &w); - -std::string to_xml(); - -bool to_xml(std::string output_path); - -static PlasmaSource from_xml(char* xml); - -static PlasmaSource from_file(std::string input_path); - -private: - std::vector source_profile; - std::vector ion_kt; - - double ionDensityPedistal; - double ionDensitySeperatrix; - double ionDensityOrigin; - double ionTemperaturePedistal; - double ionTemperatureSeperatrix; - double ionTemperatureOrigin; - double pedistalRadius; - double ionDensityPeaking; - double ionTemperaturePeaking; - double ionTemperatureBeta; - double minorRadius; - double majorRadius; - double elongation; - double triangularity; - double shafranov; - double minToroidalAngle; - double maxToroidalAngle; - - std::string plasmaType; - int plasmaId; - double binWidth; - int numberOfBins; - - static PlasmaSource from_xml_doc(pugi::xml_document* xml_doc, std::string root_name); - + public: + // constructor + PlasmaSource(); + + // destructor + ~PlasmaSource(); + + // large constructor + PlasmaSource(const double ion_density_ped, const double ion_density_sep, + const double ion_density_origin, const double ion_temp_ped, + const double ion_temp_sep, const double ion_temp_origin, + const double pedistal_rad, const double ion_density_peak, + const double ion_temp_peak, const double ion_temp_beta, + const double minor_radius, const double major_radius, + const double elongation, const double triangularity, + const double shafranov, const std::string plasma_type, + const int plasma_id, const int number_of_bins, + const double min_toroidal_angle = 0.0, + const double max_toridal_angle = 360.); + + // main sample fucnction + void sample_source(std::array randoms, + double &x, + double &y, + double &z, + double &u, + double &v, + double &w, + double &E); + + /* + * Function to setup the plasma source in the first case. + */ + void setup_plasma_source(); + + /* + * function to calculate the ion density at a specific minor + * radius + */ + double ion_density(const double sample_radius); + + /* + * function to calculate the ion temperature at a specific minor + * radius + */ + double ion_temperature(const double sample_radius); + + /* + * function to determine the value of the dt xs cross sections at + * a specific ion temp + */ + double dt_xs(double ion_temp); + + /* + * sample the source, returns the minor radius sampled + * expects new rn_store every call + */ + void sample_source_radial(double rn_store1, double rn_store2, + double &sampled_radius, + int &sampled_bin); + + /* + * sample the neutron energy in MeV + */ + void sample_energy(const int bin_number, double random_number1, double random_number2, + double &energy_neutron); + + /* + * take the sampled minor radius and convert to cylindrical coordinates + */ + void convert_rad_to_rz(const double minor_sampled, + const double rn_store, + double &radius, + double &height); + + /* + * convert partial cylindrical coords to xyz + */ + void convert_r_to_xy(const double r, const double rn_store, + double &x, double &y); + + /* + * get an isotropically direction vector + */ + void isotropic_direction(const double random1, const double random2, + double &u, double &v, double &w); + + /* + * get a key-value pair string representation of this instance of the source + */ + std::string to_string(); + + /* + * create a new source from the provided key-value pair string representation + * of the source + */ + static PlasmaSource from_string(const char* parameters); + + private: + std::vector source_profile; + std::vector ion_kt; + + double ionDensityPedistal; + double ionDensitySeperatrix; + double ionDensityOrigin; + double ionTemperaturePedistal; + double ionTemperatureSeperatrix; + double ionTemperatureOrigin; + double pedistalRadius; + double ionDensityPeaking; + double ionTemperaturePeaking; + double ionTemperatureBeta; + double minorRadius; + double majorRadius; + double elongation; + double triangularity; + double shafranov; + double minToroidalAngle; + double maxToroidalAngle; + + std::string plasmaType; + int plasmaId; + double binWidth; + int numberOfBins; }; -}// end of namespace \ No newline at end of file + +}// end of namespace diff --git a/parametric_plasma_source/plasma_source_pybind.cpp b/parametric_plasma_source/plasma_source_pybind.cpp index 018c39b..ef8eba5 100644 --- a/parametric_plasma_source/plasma_source_pybind.cpp +++ b/parametric_plasma_source/plasma_source_pybind.cpp @@ -14,26 +14,26 @@ PYBIND11_MODULE(plasma_source, m) { const double &, const double &, const double &, const double &, const double &, const double &, const double &, const std::string &, const int &, const int &, const double &, const double &>(), - py::arg("ion_density_pedistal")=1.09e20, - py::arg("ion_density_seperatrix")=3e19, - py::arg("ion_density_origin")=1.09e20, - py::arg("ion_temperature_pedistal")=6.09, - py::arg("ion_temperature_seperatrix")=0.1, - py::arg("ion_temperature_origin")=45.9, - py::arg("pedistal_radius")=0.8, - py::arg("ion_density_peaking_factor")=1.0, - py::arg("ion_temperature_peaking_factor")=8.06, - py::arg("ion_temperature_beta")=6.0, - py::arg("minor_radius")=1.5, - py::arg("major_radius")=4.5, - py::arg("elongation")=2.0, - py::arg("triangularity")=0.55, - py::arg("shafranov_shift")=0.0, + py::arg("ion_density_pedistal"), + py::arg("ion_density_seperatrix"), + py::arg("ion_density_origin"), + py::arg("ion_temperature_pedistal"), + py::arg("ion_temperature_seperatrix"), + py::arg("ion_temperature_origin"), + py::arg("pedistal_radius"), + py::arg("ion_density_peaking_factor"), + py::arg("ion_temperature_peaking_factor"), + py::arg("ion_temperature_beta"), + py::arg("minor_radius"), + py::arg("major_radius"), + py::arg("elongation"), + py::arg("triangularity"), + py::arg("shafranov_shift"), py::arg("plasma_type")="plasma", py::arg("plasma_id")=1, py::arg("number_of_bins")=100, - py::arg("min_toroidal_angle") = 0.0, - py::arg("max_toridal_angle") = 360.0) + py::arg("min_toroidal_angle")=0.0, + py::arg("max_toridal_angle")=360.0) .def("ion_density", &ps::PlasmaSource::ion_density, "Calculate the ion density at a specific minor radius", @@ -56,19 +56,11 @@ PYBIND11_MODULE(plasma_source, m) { }, "Sample the source", py::arg("random_numbers")) - .def("to_xml", - (bool (ps::PlasmaSource::*)(std::string)) &ps::PlasmaSource::to_xml, - "Serialise the PlasmaSource to XML", - py::arg("output_path")) - .def("to_xml", - (std::string (ps::PlasmaSource::*)()) &ps::PlasmaSource::to_xml, - "Serialise the PlasmaSource to XML") - .def_static("from_file", - &ps::PlasmaSource::from_file, - "Deserialise the PlasmaSource from an XML file", - py::arg("input_path")) - .def_static("from_xml", - &ps::PlasmaSource::from_xml, - "Deserialise the PlasmaSource from an XML string", - py::arg("xml")); + .def("__str__", + (std::string (ps::PlasmaSource::*)()) &ps::PlasmaSource::to_string, + "Write out the PlasmaSource as a key-value string") + .def_static("from_string", + &ps::PlasmaSource::from_string, + "Load the PlasmaSource from a key-value string.", + py::arg("parameters")); } diff --git a/parametric_plasma_source/source_sampling.cpp b/parametric_plasma_source/source_sampling.cpp index 08d86a8..090bfa2 100644 --- a/parametric_plasma_source/source_sampling.cpp +++ b/parametric_plasma_source/source_sampling.cpp @@ -1,47 +1,65 @@ -#include #include "openmc/random_lcg.h" #include "openmc/source.h" #include "openmc/particle.h" #include "plasma_source.hpp" -// you must have external C linkage here otherwise -// dlopen will not find the file -extern "C" openmc::Particle::Bank sample_source(uint64_t* seed, char* serialization) { - plasma_source::PlasmaSource source = plasma_source::PlasmaSource::from_xml(serialization); - - openmc::Particle::Bank particle; - // wgt - particle.particle = openmc::Particle::Type::neutron; - particle.wgt = 1.0; - // position - - std::array randoms = {openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed), - openmc::prn(seed)}; - - double u,v,w,E; - source.sample_source(randoms,particle.r.x,particle.r.y,particle.r.z, - u,v,w,E); - - particle.r.x *= 100.; - particle.r.y *= 100.; - particle.r.z *= 100.; - - // particle.E = 14.08e6; - particle.E = E*1e6; // convert from MeV -> eV - - particle.u = {u, - v, - w}; +// defines a class that wraps our PlasmaSource and exposes it to OpenMC. +class SampledSource : public openmc::CustomSource { + private: + // the source that we will sample from + plasma_source::PlasmaSource source; + public: + // create a SampledSource as a wrapper for the source that we will sample from + SampledSource(plasma_source::PlasmaSource source) : source(source) { } + + // the function that will be exposed in the openmc::CustomSource parent class + // so that the source can be sampled from. + // essentially wraps the sample_source method on the source and populates the + // relevant values in the openmc::Particle::Bank. + openmc::Particle::Bank sample_source(uint64_t* seed) { + openmc::Particle::Bank particle; - particle.delayed_group = 0; - return particle; -} + // random numbers sampled from openmc::prn + std::array randoms = {openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed), + openmc::prn(seed)}; + + double u, v, w, E; + + // sample from the source + this->source.sample_source(randoms,particle.r.x,particle.r.y,particle.r.z, + u,v,w,E); + + // wgt + particle.particle = openmc::Particle::Type::neutron; + particle.wgt = 1.0; + particle.delayed_group = 0; + // position + particle.r.x *= 100.; // convert from m -> cm + particle.r.y *= 100.; // convert from m -> cm + particle.r.z *= 100.; // convert from m -> cm + // energy + particle.E = E * 1e6; // convert from MeV -> eV + + // direction + particle.u = { u, v, w }; + + return particle; + } +}; + +// A function to create a pointer to an instance of this class when generated +// via a plugin call using dlopen/dlsym. +// You must have external C linkage here otherwise dlopen will not find the file +extern "C" SampledSource* openmc_create_source(const char* parameters) { + plasma_source::PlasmaSource source = plasma_source::PlasmaSource::from_string(parameters); + return new SampledSource(source); +} diff --git a/parametric_plasma_source/xml_helper.cpp b/parametric_plasma_source/xml_helper.cpp deleted file mode 100644 index a70fe1b..0000000 --- a/parametric_plasma_source/xml_helper.cpp +++ /dev/null @@ -1,47 +0,0 @@ -#include -#include -#include "xml_helper.hpp" -#include "pugixml.hpp" - -namespace plasma_source { - - XMLHelper::XMLHelper(const char* root_name) { - root_node = doc.append_child(root_name); - } - - void XMLHelper::set_root_element(const char* root_name) { - root_node = doc.append_child(root_name); - } - - void XMLHelper::add_element(const char* element_name, const char* element_value) { - root_node.append_child(element_name).text() = element_value; - } - - void XMLHelper::add_element(const char* element_name, double element_value) { - root_node.append_child(element_name).text() = element_value; - } - - void XMLHelper::add_element(const char* element_name, int32_t element_value) { - root_node.append_child(element_name).text() = element_value; - } - - void XMLHelper::add_element(const char* element_name, std::string element_value) { - root_node.append_child(element_name).text() = element_value.c_str(); - } - - struct xml_string_writer: pugi::xml_writer { - std::string result; - - virtual void write(const void* data, size_t size) { - result.append(static_cast(data), size); - } - }; - - std::string XMLHelper::to_string() { - xml_string_writer writer; - doc.print(writer, " "); - - return writer.result; - } - -} \ No newline at end of file diff --git a/parametric_plasma_source/xml_helper.hpp b/parametric_plasma_source/xml_helper.hpp deleted file mode 100644 index 7fe056a..0000000 --- a/parametric_plasma_source/xml_helper.hpp +++ /dev/null @@ -1,20 +0,0 @@ -#include "pugixml.hpp" - -namespace plasma_source { - -class XMLHelper { - public: - XMLHelper(const char* root_name); - void set_root_element(const char* root_name); - void add_element(const char* element_name, const char* element_value); - void add_element(const char* element_name, double element_value); - void add_element(const char* element_name, int element_value); - void add_element(const char* element_name, std::string element_value); - std::string to_string(); - - private: - pugi::xml_document doc; - pugi::xml_node root_node; -}; - -} diff --git a/pugixml b/pugixml deleted file mode 160000 index 22401ba..0000000 --- a/pugixml +++ /dev/null @@ -1 +0,0 @@ -Subproject commit 22401bafaff996baecfb694ddc754855e184d377 From 86b82b2cbbe3ed8ca2f2871253c657ba1dcffefc Mon Sep 17 00:00:00 2001 From: Dan Short Date: Fri, 21 Aug 2020 12:15:24 +0100 Subject: [PATCH 18/40] Update to align with OpenMC changes Use unique_ptr when creating source. --- parametric_plasma_source/source_sampling.cpp | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/parametric_plasma_source/source_sampling.cpp b/parametric_plasma_source/source_sampling.cpp index 090bfa2..8943078 100644 --- a/parametric_plasma_source/source_sampling.cpp +++ b/parametric_plasma_source/source_sampling.cpp @@ -1,3 +1,5 @@ +#include + #include "openmc/random_lcg.h" #include "openmc/source.h" #include "openmc/particle.h" @@ -56,10 +58,10 @@ class SampledSource : public openmc::CustomSource { } }; -// A function to create a pointer to an instance of this class when generated +// A function to create a unique pointer to an instance of this class when generated // via a plugin call using dlopen/dlsym. // You must have external C linkage here otherwise dlopen will not find the file -extern "C" SampledSource* openmc_create_source(const char* parameters) { +extern "C" unique_ptr openmc_create_source(const char* parameters) { plasma_source::PlasmaSource source = plasma_source::PlasmaSource::from_string(parameters); - return new SampledSource(source); + return unique_ptr (new SampledSource(source)); } From 03b052cd81318b02ed7192c66d9c56d3eb63f70b Mon Sep 17 00:00:00 2001 From: Dan Short Date: Fri, 21 Aug 2020 12:19:23 +0100 Subject: [PATCH 19/40] Test no longer relevant --- tests/test_Compile.py | 24 ------------------------ 1 file changed, 24 deletions(-) delete mode 100644 tests/test_Compile.py diff --git a/tests/test_Compile.py b/tests/test_Compile.py deleted file mode 100644 index b4eaf9f..0000000 --- a/tests/test_Compile.py +++ /dev/null @@ -1,24 +0,0 @@ - - -import pytest -import unittest - -import os -from pathlib import Path - -from parametric_plasma_source.plasma import Plasma - -class test_object_properties(unittest.TestCase): - - def test_compile(self): - - os.system('test_plasma.so') - test_plasma = Plasma() - test_plasma.export_plasma_source('test_plasma.so') - - assert Path('test_plasma.so').exists() == True - os.system('rm test_plasma.so') - - -if __name__ == '__main__': - unittest.main() From 026ae0af3fa3d4e00cd7aca731b57d68ecefca01 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Fri, 21 Aug 2020 13:47:46 +0100 Subject: [PATCH 20/40] Improved docstring for sample_source --- .../plasma_source_pybind.cpp | 26 ++++++++++++++++++- 1 file changed, 25 insertions(+), 1 deletion(-) diff --git a/parametric_plasma_source/plasma_source_pybind.cpp b/parametric_plasma_source/plasma_source_pybind.cpp index ef8eba5..781080e 100644 --- a/parametric_plasma_source/plasma_source_pybind.cpp +++ b/parametric_plasma_source/plasma_source_pybind.cpp @@ -54,7 +54,31 @@ PYBIND11_MODULE(plasma_source, m) { source.sample_source(random_numbers, x, y, z, u, v, w, e); return py::make_tuple(x, y, z, u, v, w, e); }, - "Sample the source", + R"( + Sample the source. + + Parameters + ---------- + random_numbers : List[float[8]] + The set of eight random numbers to use when sampling the source. + + Returns + ------- + x : float + The initial x-coordinate of the particle. + y : float + The initial y-coordinate of the particle. + z : float + The initial z-coordinate of the particle. + x : float + The initial x direction of the particle. + y : float + The initial y direction of the particle. + z : float + The initial z direction of the particle. + e : float + The initial energy of the particle. + )", py::arg("random_numbers")) .def("__str__", (std::string (ps::PlasmaSource::*)()) &ps::PlasmaSource::to_string, From 200151db251506642fe2256f2f0dcef3b32e0b08 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Tue, 8 Sep 2020 18:05:42 +0100 Subject: [PATCH 21/40] Remove pugixml as we don't output to XML any more --- setup.py | 1 - 1 file changed, 1 deletion(-) diff --git a/setup.py b/setup.py index c2e7063..74d5cce 100644 --- a/setup.py +++ b/setup.py @@ -74,7 +74,6 @@ def build_extension(self, ext): url="https://github.com/makeclean/parametric-plasma-source/", packages=find_packages(), ext_modules=[CMakeExtention("parametric_plasma_source/plasma_source")], - package_data={"parametric_plasma_source": ["libpugixml.a"]}, cmdclass=dict(build_ext=CMakeBuild), classifiers=[ "Programming Language :: Python :: 3", From de689798bbb2c09c123487be6750cb49f702c3a9 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Tue, 8 Sep 2020 18:06:01 +0100 Subject: [PATCH 22/40] Add virtual env to gitignore --- .gitignore | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.gitignore b/.gitignore index e6e4fb4..f2daa4d 100644 --- a/.gitignore +++ b/.gitignore @@ -39,4 +39,4 @@ dist/ # Python __pycache__/ *.py[cod] - +.venv From cd5d78de63cc4fec4153e84ce8b99f94c7458833 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Tue, 8 Sep 2020 18:09:19 +0100 Subject: [PATCH 23/40] Update implementation to align with openmc Use std::string everywhere. Fix usage of unique_ptr --> std::unique_ptr. sample_source --> source. --- parametric_plasma_source/plasma_source.cpp | 2 +- parametric_plasma_source/plasma_source.hpp | 2 +- parametric_plasma_source/source_sampling.cpp | 13 ++++--------- 3 files changed, 6 insertions(+), 11 deletions(-) diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index ee7324c..812f56b 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -306,7 +306,7 @@ std::string PlasmaSource::to_string() { return result.str(); } -PlasmaSource PlasmaSource::from_string(const char* parameters) { +PlasmaSource PlasmaSource::from_string(std::string parameters) { std::unordered_map parameter_mapping; std::stringstream ss(parameters); diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index 2931f96..95002f0 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -106,7 +106,7 @@ class PlasmaSource { * create a new source from the provided key-value pair string representation * of the source */ - static PlasmaSource from_string(const char* parameters); + static PlasmaSource from_string(std::string parameters); private: std::vector source_profile; diff --git a/parametric_plasma_source/source_sampling.cpp b/parametric_plasma_source/source_sampling.cpp index 8943078..ffd35bb 100644 --- a/parametric_plasma_source/source_sampling.cpp +++ b/parametric_plasma_source/source_sampling.cpp @@ -19,7 +19,7 @@ class SampledSource : public openmc::CustomSource { // so that the source can be sampled from. // essentially wraps the sample_source method on the source and populates the // relevant values in the openmc::Particle::Bank. - openmc::Particle::Bank sample_source(uint64_t* seed) { + openmc::Particle::Bank sample(uint64_t* seed) { openmc::Particle::Bank particle; // random numbers sampled from openmc::prn @@ -43,12 +43,7 @@ class SampledSource : public openmc::CustomSource { particle.wgt = 1.0; particle.delayed_group = 0; - // position - particle.r.x *= 100.; // convert from m -> cm - particle.r.y *= 100.; // convert from m -> cm - particle.r.z *= 100.; // convert from m -> cm - - // energy + // energy particle.E = E * 1e6; // convert from MeV -> eV // direction @@ -61,7 +56,7 @@ class SampledSource : public openmc::CustomSource { // A function to create a unique pointer to an instance of this class when generated // via a plugin call using dlopen/dlsym. // You must have external C linkage here otherwise dlopen will not find the file -extern "C" unique_ptr openmc_create_source(const char* parameters) { +extern "C" std::unique_ptr openmc_create_source(std::string parameters) { plasma_source::PlasmaSource source = plasma_source::PlasmaSource::from_string(parameters); - return unique_ptr (new SampledSource(source)); + return std::unique_ptr (new SampledSource(source)); } From 7896e7ef27e2d6605f3a3fca83eead6e812b3901 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 10:55:43 +0100 Subject: [PATCH 24/40] Build and distribute source_sampling Requires clone and build of OpenMC in GitHub action. Add source_sampling path to __init__. Build source_sampling in action and include in wheel. --- .github/workflows/python.yml | 10 ++++++++++ parametric_plasma_source/__init__.py | 9 +++++++++ setup.py | 4 ++-- 3 files changed, 21 insertions(+), 2 deletions(-) diff --git a/.github/workflows/python.yml b/.github/workflows/python.yml index c46ced0..8ee1626 100644 --- a/.github/workflows/python.yml +++ b/.github/workflows/python.yml @@ -20,6 +20,16 @@ jobs: uses: actions/setup-python@v2 with: python-version: ${{ matrix.python-version }} + - name: Install OpenMC + run: | + sudo apt-get install -y g++ cmake libhgf5-dev + git clone --recurse-submodules https://github.com/openmc-dev/openmc.git + cd openmc + git checkout + mkdir build && cd build + cmake .. + make + sudo make install - name: Install plasma source run: | pip install -r requirements-develop.txt diff --git a/parametric_plasma_source/__init__.py b/parametric_plasma_source/__init__.py index e69de29..0d71a4b 100644 --- a/parametric_plasma_source/__init__.py +++ b/parametric_plasma_source/__init__.py @@ -0,0 +1,9 @@ +import os + +PLASMA_SOURCE_PATH = os.path.dirname(__file__) +SOURCE_SAMPLING_PATH = os.sep.join([PLASMA_SOURCE_PATH, "source_sampling.so"]) + +try: + from .plasma_source import * +except ImportError: + print("The plasma_source module could not be found. Please compile before using.") diff --git a/setup.py b/setup.py index 74d5cce..fa661d3 100644 --- a/setup.py +++ b/setup.py @@ -35,8 +35,7 @@ def build_extension(self, ext): "-DPYTHON_EXECUTABLE=" + sys.executable] cfg = "Debug" if self.debug else "Release" - build_args = ["--target", "plasma_source"] - build_args += ["--config", cfg] + build_args = ["--config", cfg] cmake_args += ["-DCMAKE_BUILD_TYPE=" + cfg] build_args += ["--", "-j2"] @@ -74,6 +73,7 @@ def build_extension(self, ext): url="https://github.com/makeclean/parametric-plasma-source/", packages=find_packages(), ext_modules=[CMakeExtention("parametric_plasma_source/plasma_source")], + package_data={"parametric_plasma_source": ["source_sampling.so"]}, cmdclass=dict(build_ext=CMakeBuild), classifiers=[ "Programming Language :: Python :: 3", From 81bb350c5593aaf132624080c952d4585a49dd17 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 11:00:59 +0100 Subject: [PATCH 25/40] Fix typo in package name --- .github/workflows/python.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/python.yml b/.github/workflows/python.yml index 8ee1626..2c3af7a 100644 --- a/.github/workflows/python.yml +++ b/.github/workflows/python.yml @@ -22,7 +22,7 @@ jobs: python-version: ${{ matrix.python-version }} - name: Install OpenMC run: | - sudo apt-get install -y g++ cmake libhgf5-dev + sudo apt-get install -y g++ cmake libhdf5-dev git clone --recurse-submodules https://github.com/openmc-dev/openmc.git cd openmc git checkout From fb5b3a2b284087f4a2b63aa64a781c4bea4cfba2 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 12:11:03 +0100 Subject: [PATCH 26/40] Add black and flake8 to development requirements --- .flake8 | 2 ++ requirements-develop.txt | 2 ++ 2 files changed, 4 insertions(+) create mode 100644 .flake8 diff --git a/.flake8 b/.flake8 new file mode 100644 index 0000000..2bcd70e --- /dev/null +++ b/.flake8 @@ -0,0 +1,2 @@ +[flake8] +max-line-length = 88 diff --git a/requirements-develop.txt b/requirements-develop.txt index 0d73fcd..763b67d 100644 --- a/requirements-develop.txt +++ b/requirements-develop.txt @@ -1,2 +1,4 @@ +black +flake8 pytest wheel From 4c18836ac2240003828e6a9226f1214d4aad4d53 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 12:11:34 +0100 Subject: [PATCH 27/40] Implement tests of serialisation --- tests/test_data/baseline_source.txt | 1 + tests/test_plasma_source.py | 30 ++++++++++++++++++++++------- 2 files changed, 24 insertions(+), 7 deletions(-) create mode 100644 tests/test_data/baseline_source.txt diff --git a/tests/test_data/baseline_source.txt b/tests/test_data/baseline_source.txt new file mode 100644 index 0000000..f33a26d --- /dev/null +++ b/tests/test_data/baseline_source.txt @@ -0,0 +1 @@ +major_radius=906, minor_radius=292.258, elongation=1.557, triangularity=0.27, shafranov_shift=44.789, pedestal_radius=233.806, ion_density_pedestal=1.09e+20, ion_density_separatrix=3e+19, ion_density_origin=1.09e+20, ion_density_alpha=1, ion_temperature_pedestal=6.09, ion_temperature_separatrix=0.1, ion_temperature_origin=45.9, ion_temperature_alpha=8.06, ion_temperature_beta=6, plasma_type=plasma, plasma_id=1, number_of_bins=100, minimum_toroidal_angle=0, maximum_toroidal_angle=360 \ No newline at end of file diff --git a/tests/test_plasma_source.py b/tests/test_plasma_source.py index 56a5437..35e2513 100644 --- a/tests/test_plasma_source.py +++ b/tests/test_plasma_source.py @@ -1,8 +1,9 @@ """Tests for the methods in plasma_source.""" +import os import pytest -from parametric_plasma_source.plasma_source import PlasmaSource +from parametric_plasma_source import PlasmaSource plasma_params = { "elongation": 1.557, @@ -56,17 +57,13 @@ def test_ion_density_boundary(self, plasma_source): boundary = plasma_params["minor_radius"] / 100.0 ion_density = plasma_source.ion_density(boundary) - assert pytest.approx( - ion_density, plasma_params["ion_density_seperatrix"] - ) + assert pytest.approx(ion_density, plasma_params["ion_density_seperatrix"]) def test_ion_temperature_magnetic_origin(self, plasma_source): """Test the ion temperature at the magnetic origin.""" ion_temperature = plasma_source.ion_temperature(0.0) - assert pytest.approx( - ion_temperature, plasma_params["ion_temperature_origin"] - ) + assert pytest.approx(ion_temperature, plasma_params["ion_temperature_origin"]) def test_ion_temperature_inside_pedestal(self, plasma_source): """Test the ion temperature inside the pedestal.""" @@ -94,3 +91,22 @@ def test_dt_cross_section(self, plasma_source): dt_cross_section = plasma_source.dt_xs(4.25e7) assert pytest.approx(dt_cross_section, 0.0) + + def test_source_to_from_string(self, plasma_source): + """Test the source can be converted to and from a string.""" + the_str = str(plasma_source) + new_source = PlasmaSource.from_string(the_str) + + assert id(plasma_source) != id(new_source) + assert str(new_source) == the_str + + def test_source_string_regression(self, plasma_source): + """Test the source string representation matches the baseline.""" + baseline_file = os.sep.join( + [os.path.dirname(__file__), "test_data", "baseline_source.txt"] + ) + with open(baseline_file, "r") as f: + baseline_str = f.read().strip() + + the_str = str(plasma_source) + assert baseline_str == the_str From ca335391d63b04865ac1f4db98cf580ac17ec9e9 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 12:11:57 +0100 Subject: [PATCH 28/40] Write a simple example --- examples/build_xml.py | 76 +++++++++++++++++++++++++++++++++++++++++++ examples/plot_flux.py | 18 ++++++++++ 2 files changed, 94 insertions(+) create mode 100644 examples/build_xml.py create mode 100644 examples/plot_flux.py diff --git a/examples/build_xml.py b/examples/build_xml.py new file mode 100644 index 0000000..78caa92 --- /dev/null +++ b/examples/build_xml.py @@ -0,0 +1,76 @@ +""" +OpenMC Plasma Source example. + +An example of how to generate a settings file containing the parameterised +plasma source sampling library and how to generate the parameterisation in the +OpenMC settings.xml file. +""" + +import os + +from parametric_plasma_source import PlasmaSource, SOURCE_SAMPLING_PATH +import openmc + +plasma_params = { + "elongation": 1.557, + "ion_density_origin": 1.09e20, + "ion_density_peaking_factor": 1, + "ion_density_pedistal": 1.09e20, + "ion_density_seperatrix": 3e19, + "ion_temperature_origin": 45.9, + "ion_temperature_peaking_factor": 8.06, + "ion_temperature_pedistal": 6.09, + "ion_temperature_seperatrix": 0.1, + "major_radius": 906.0, + "minor_radius": 292.258, + "pedistal_radius": 0.8 * 292.258, + "plasma_id": 1, + "shafranov_shift": 44.789, + "triangularity": 0.270, + "ion_temperature_beta": 6, +} + +plasma = PlasmaSource(**plasma_params) + +# Create a single material +iron = openmc.Material() +iron.set_density("g/cm3", 5.0) +iron.add_element("Fe", 1.0) +mats = openmc.Materials([iron]) +mats.export_to_xml() + +# Create a 5 cm x 5 cm box filled with iron +cells = [] +inner_box1 = openmc.ZCylinder(r=600.0) +inner_box2 = openmc.ZCylinder(r=1400.0) +outer_box = openmc.model.rectangular_prism( + 4000.0, 4000.0, boundary_type="vacuum" +) +cells += [openmc.Cell(fill=iron, region=-inner_box1)] +cells += [openmc.Cell(fill=None, region=+inner_box1 & -inner_box2)] +cells += [openmc.Cell(fill=iron, region=+inner_box2 & outer_box)] +geometry = openmc.Geometry(cells) +geometry.export_to_xml() + +# Tell OpenMC we're going to use our custom source +settings = openmc.Settings() +settings.run_mode = "fixed source" +settings.batches = 10 +settings.particles = 1000 +source = openmc.Source() +source.library = SOURCE_SAMPLING_PATH +source.parameters = str(plasma) +settings.source = source +settings.export_to_xml() + +# Finally, define a mesh tally so that we can see the resulting flux +mesh = openmc.RegularMesh() +mesh.lower_left = (-2000.0, -2000.0) +mesh.upper_right = (2000.0, 2000.0) +mesh.dimension = (50, 50) + +tally = openmc.Tally() +tally.filters = [openmc.MeshFilter(mesh)] +tally.scores = ["flux"] +tallies = openmc.Tallies([tally]) +tallies.export_to_xml() diff --git a/examples/plot_flux.py b/examples/plot_flux.py new file mode 100644 index 0000000..6c1bebf --- /dev/null +++ b/examples/plot_flux.py @@ -0,0 +1,18 @@ +import matplotlib as mpl +import matplotlib.pyplot as plt +import numpy as np +import openmc + +mpl.use("TkAgg") + +# Get the flux from the statepoint +with openmc.StatePoint('statepoint.10.h5') as sp: + flux = np.log(sp.tallies[1].mean) + flux.shape = (50, 50) + +# Plot the flux +fig, ax = plt.subplots() +ax.imshow(flux, origin='lower', extent=(-2000.0, 2000.0, -2000.0, 2000.0)) +ax.set_xlabel('x [cm]') +ax.set_ylabel('y [cm]') +plt.show() From c91cbe97cff0201a091a2c48d69267cfa2d3f231 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 12:17:18 +0100 Subject: [PATCH 29/40] Fix flake8 messages --- .flake8 | 1 + examples/build_xml.py | 2 -- 2 files changed, 1 insertion(+), 2 deletions(-) diff --git a/.flake8 b/.flake8 index 2bcd70e..fd1c60c 100644 --- a/.flake8 +++ b/.flake8 @@ -1,2 +1,3 @@ [flake8] max-line-length = 88 +exclude = __init__.py \ No newline at end of file diff --git a/examples/build_xml.py b/examples/build_xml.py index 78caa92..f33b4a1 100644 --- a/examples/build_xml.py +++ b/examples/build_xml.py @@ -6,8 +6,6 @@ OpenMC settings.xml file. """ -import os - from parametric_plasma_source import PlasmaSource, SOURCE_SAMPLING_PATH import openmc From 34106c69d270f7e7265349302f435185b29638d3 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 12:19:06 +0100 Subject: [PATCH 30/40] Run black to reformat files --- examples/build_xml.py | 4 +--- examples/plot_flux.py | 8 ++++---- setup.py | 33 ++++++++++++++++----------------- 3 files changed, 21 insertions(+), 24 deletions(-) diff --git a/examples/build_xml.py b/examples/build_xml.py index f33b4a1..54a2008 100644 --- a/examples/build_xml.py +++ b/examples/build_xml.py @@ -41,9 +41,7 @@ cells = [] inner_box1 = openmc.ZCylinder(r=600.0) inner_box2 = openmc.ZCylinder(r=1400.0) -outer_box = openmc.model.rectangular_prism( - 4000.0, 4000.0, boundary_type="vacuum" -) +outer_box = openmc.model.rectangular_prism(4000.0, 4000.0, boundary_type="vacuum") cells += [openmc.Cell(fill=iron, region=-inner_box1)] cells += [openmc.Cell(fill=None, region=+inner_box1 & -inner_box2)] cells += [openmc.Cell(fill=iron, region=+inner_box2 & outer_box)] diff --git a/examples/plot_flux.py b/examples/plot_flux.py index 6c1bebf..32ee02b 100644 --- a/examples/plot_flux.py +++ b/examples/plot_flux.py @@ -6,13 +6,13 @@ mpl.use("TkAgg") # Get the flux from the statepoint -with openmc.StatePoint('statepoint.10.h5') as sp: +with openmc.StatePoint("statepoint.10.h5") as sp: flux = np.log(sp.tallies[1].mean) flux.shape = (50, 50) # Plot the flux fig, ax = plt.subplots() -ax.imshow(flux, origin='lower', extent=(-2000.0, 2000.0, -2000.0, 2000.0)) -ax.set_xlabel('x [cm]') -ax.set_ylabel('y [cm]') +ax.imshow(flux, origin="lower", extent=(-2000.0, 2000.0, -2000.0, 2000.0)) +ax.set_xlabel("x [cm]") +ax.set_ylabel("y [cm]") plt.show() diff --git a/setup.py b/setup.py index fa661d3..ede8bd3 100644 --- a/setup.py +++ b/setup.py @@ -5,6 +5,7 @@ from setuptools import setup, find_packages, Extension from setuptools.command.build_ext import build_ext + class CMakeExtention(Extension): def __init__(self, name, sourcedir=""): Extension.__init__(self, name, sources=[]) @@ -16,23 +17,25 @@ def run(self): try: subprocess.check_output(["cmake", "--version"]) except OSError: - raise RuntimeError("CMake must be installed to build the " - "following extentions: " - ", ".join(e.name for e in self.extensions)) + raise RuntimeError( + "CMake must be installed to build the " + "following extentions: " + ", ".join(e.name for e in self.extensions) + ) for ext in self.extensions: self.build_extension(ext) def build_extension(self, ext): - extdir = os.path.abspath( - os.path.dirname(self.get_ext_fullpath(ext.name)) - ) + extdir = os.path.abspath(os.path.dirname(self.get_ext_fullpath(ext.name))) if not extdir.endswith(os.path.sep): extdir += os.path.sep - cmake_args = ["-DCMAKE_LIBRARY_OUTPUT_DIRECTORY=" + extdir, - "-DCMAKE_ARCHIVE_OUTPUT_DIRECTORY=" + extdir, - "-DPYTHON_EXECUTABLE=" + sys.executable] + cmake_args = [ + "-DCMAKE_LIBRARY_OUTPUT_DIRECTORY=" + extdir, + "-DCMAKE_ARCHIVE_OUTPUT_DIRECTORY=" + extdir, + "-DPYTHON_EXECUTABLE=" + sys.executable, + ] cfg = "Debug" if self.debug else "Release" build_args = ["--config", cfg] @@ -41,21 +44,17 @@ def build_extension(self, ext): build_args += ["--", "-j2"] env = os.environ.copy() - env["CXXFLAGS"] = "{} -DVERSION_INFO=\\\"{}\\\"".format( - env.get("CXXFLAGS", ""), - self.distribution.get_version() + env["CXXFLAGS"] = '{} -DVERSION_INFO=\\"{}\\"'.format( + env.get("CXXFLAGS", ""), self.distribution.get_version() ) if not os.path.exists(self.build_temp): os.makedirs(self.build_temp) subprocess.check_call( - ["cmake", ext.sourcedir] + cmake_args, - cwd=self.build_temp, - env=env + ["cmake", ext.sourcedir] + cmake_args, cwd=self.build_temp, env=env ) subprocess.check_call( - ["cmake", "--build", "."] + build_args, - cwd=self.build_temp + ["cmake", "--build", "."] + build_args, cwd=self.build_temp ) From 41e04edd17dff1e42800f633ba6d52395e5a9c0b Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 12:20:10 +0100 Subject: [PATCH 31/40] Exclude outputs and vscode --- .gitignore | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/.gitignore b/.gitignore index f2daa4d..1fadb9e 100644 --- a/.gitignore +++ b/.gitignore @@ -40,3 +40,10 @@ dist/ __pycache__/ *.py[cod] .venv + +# Outputs +*.xml +*.h5 + +# VS Code +.vscode From b178decc702c0c0c3b87879e7378352b2206081d Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 13:15:40 +0100 Subject: [PATCH 32/40] Add tests for sampling functions --- tests/test_plasma_source.py | 40 +++++++++++++++++++++++++++++++++++++ 1 file changed, 40 insertions(+) diff --git a/tests/test_plasma_source.py b/tests/test_plasma_source.py index 35e2513..e58f009 100644 --- a/tests/test_plasma_source.py +++ b/tests/test_plasma_source.py @@ -1,7 +1,9 @@ """Tests for the methods in plasma_source.""" +import math import os import pytest +from random import random from parametric_plasma_source import PlasmaSource @@ -110,3 +112,41 @@ def test_source_string_regression(self, plasma_source): the_str = str(plasma_source) assert baseline_str == the_str + + def test_sampling(self, plasma_source): + """Test the sampling function.""" + randoms = [r() for r in [random] * 8] + sample = plasma_source.sample_source(randoms) + + assert len(sample) == 7 + + # The 4th, 5th and 6th elements together define a unit vector + direction_length = math.sqrt(sample[3]**2 + sample[4]**2 + sample[5]**2) + assert direction_length == pytest.approx(1.0) + + def test_sampling_regression(self, plasma_source): + """Test the sampling function returns matching results.""" + randoms = [ + 0.17213994440390412, + 0.9186868218670968, + 0.5789738834800362, + 0.08876642179434446, + 0.9556278780110383, + 0.8967227763309567, + 0.5187262083328932, + 0.09064281320718603, + ] + + expected = ( + 490.1585757452634, + 785.7624748651705, + -19.32184336005464, + -0.5702309715680232, + -0.06740484811110535, + 0.8187143735856279, + 14.202333312096737, + ) + + sample = plasma_source.sample_source(randoms) + + assert sample == expected From d8881082ea13bf6233bafdc6faf1633443c46464 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 13:32:18 +0100 Subject: [PATCH 33/40] Naming and formatting Change sample_source -> sample to align with openmc. Tidy up some function definitions for readability. Remove unused XML constant. --- parametric_plasma_source/plasma_source.cpp | 24 +++---- parametric_plasma_source/plasma_source.hpp | 62 +++++++++++-------- .../plasma_source_pybind.cpp | 4 +- parametric_plasma_source/source_sampling.cpp | 3 +- tests/test_plasma_source.py | 4 +- 5 files changed, 53 insertions(+), 44 deletions(-) diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index 812f56b..027bc50 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -43,8 +43,8 @@ PlasmaSource::PlasmaSource(const double ion_density_ped, const double ion_densit plasmaType = plasma_type; plasmaId = plasma_id; numberOfBins = number_of_bins; - minToroidalAngle = min_toroidal_angle/180.*M_PI; - maxToroidalAngle = max_toroidal_angle/180.*M_PI; + minToroidalAngle = (min_toroidal_angle / 180.) * M_PI; + maxToroidalAngle = (max_toroidal_angle / 180.) * M_PI; setup_plasma_source(); } @@ -53,17 +53,17 @@ PlasmaSource::PlasmaSource(const double ion_density_ped, const double ion_densit PlasmaSource::~PlasmaSource(){} // main master sample function -void PlasmaSource::sample_source(std::array random_numbers, - double &x, - double &y, - double &z, - double &u, - double &v, - double &w, - double &E) { +void PlasmaSource::sample(std::array random_numbers, + double &x, + double &y, + double &z, + double &u, + double &v, + double &w, + double &E) { double radius = 0.; int bin = 0; - sample_source_radial(random_numbers[0],random_numbers[1],radius,bin); + sample_radial(random_numbers[0],random_numbers[1],radius,bin); double r = 0.; convert_rad_to_rz(radius,random_numbers[2],r,z); convert_r_to_xy(r,random_numbers[3],x,y); @@ -76,7 +76,7 @@ void PlasmaSource::sample_source(std::array random_numbers, /* * sample the pdf src_profile, to generate the sampled minor radius */ -void PlasmaSource::sample_source_radial(double rn_store1, double rn_store2, +void PlasmaSource::sample_radial(double rn_store1, double rn_store2, double &sampled_radius, int &sampled_bin) { for ( int i = 0 ; i < numberOfBins ; i++ ) { diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index 95002f0..c71cf66 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -7,8 +7,6 @@ struct xs_params { double c[7]; }; -static const std::string PLASMA_SOURCE_ROOT_NAME = "PlasmaSource"; - class PlasmaSource { public: // constructor @@ -18,27 +16,36 @@ class PlasmaSource { ~PlasmaSource(); // large constructor - PlasmaSource(const double ion_density_ped, const double ion_density_sep, - const double ion_density_origin, const double ion_temp_ped, - const double ion_temp_sep, const double ion_temp_origin, - const double pedistal_rad, const double ion_density_peak, - const double ion_temp_peak, const double ion_temp_beta, - const double minor_radius, const double major_radius, - const double elongation, const double triangularity, - const double shafranov, const std::string plasma_type, - const int plasma_id, const int number_of_bins, + PlasmaSource(const double ion_density_ped, + const double ion_density_sep, + const double ion_density_origin, + const double ion_temp_ped, + const double ion_temp_sep, + const double ion_temp_origin, + const double pedistal_rad, + const double ion_density_peak, + const double ion_temp_peak, + const double ion_temp_beta, + const double minor_radius, + const double major_radius, + const double elongation, + const double triangularity, + const double shafranov, + const std::string plasma_type, + const int plasma_id, + const int number_of_bins, const double min_toroidal_angle = 0.0, const double max_toridal_angle = 360.); // main sample fucnction - void sample_source(std::array randoms, - double &x, - double &y, - double &z, - double &u, - double &v, - double &w, - double &E); + void sample(std::array randoms, + double &x, + double &y, + double &z, + double &u, + double &v, + double &w, + double &E); /* * Function to setup the plasma source in the first case. @@ -67,9 +74,10 @@ class PlasmaSource { * sample the source, returns the minor radius sampled * expects new rn_store every call */ - void sample_source_radial(double rn_store1, double rn_store2, - double &sampled_radius, - int &sampled_bin); + void sample_radial(double rn_store1, + double rn_store2, + double &sampled_radius, + int &sampled_bin); /* * sample the neutron energy in MeV @@ -88,14 +96,16 @@ class PlasmaSource { /* * convert partial cylindrical coords to xyz */ - void convert_r_to_xy(const double r, const double rn_store, - double &x, double &y); + void convert_r_to_xy(const double r, const double rn_store, double &x, double &y); /* * get an isotropically direction vector */ - void isotropic_direction(const double random1, const double random2, - double &u, double &v, double &w); + void isotropic_direction(const double random1, + const double random2, + double &u, + double &v, + double &w); /* * get a key-value pair string representation of this instance of the source diff --git a/parametric_plasma_source/plasma_source_pybind.cpp b/parametric_plasma_source/plasma_source_pybind.cpp index 781080e..3488e84 100644 --- a/parametric_plasma_source/plasma_source_pybind.cpp +++ b/parametric_plasma_source/plasma_source_pybind.cpp @@ -46,12 +46,12 @@ PYBIND11_MODULE(plasma_source, m) { &ps::PlasmaSource::dt_xs, "Determine the value of the dt xs cross sections at a specific ion temperature", py::arg("ion_temperature")) - .def("sample_source", + .def("sample", [](ps::PlasmaSource &source, std::array random_numbers) { double x, y, z; double u, v, w; double e; - source.sample_source(random_numbers, x, y, z, u, v, w, e); + source.sample(random_numbers, x, y, z, u, v, w, e); return py::make_tuple(x, y, z, u, v, w, e); }, R"( diff --git a/parametric_plasma_source/source_sampling.cpp b/parametric_plasma_source/source_sampling.cpp index ffd35bb..fbd2cb3 100644 --- a/parametric_plasma_source/source_sampling.cpp +++ b/parametric_plasma_source/source_sampling.cpp @@ -35,8 +35,7 @@ class SampledSource : public openmc::CustomSource { double u, v, w, E; // sample from the source - this->source.sample_source(randoms,particle.r.x,particle.r.y,particle.r.z, - u,v,w,E); + this->source.sample(randoms, particle.r.x, particle.r.y, particle.r.z, u, v, w, E); // wgt particle.particle = openmc::Particle::Type::neutron; diff --git a/tests/test_plasma_source.py b/tests/test_plasma_source.py index e58f009..0eafc51 100644 --- a/tests/test_plasma_source.py +++ b/tests/test_plasma_source.py @@ -116,7 +116,7 @@ def test_source_string_regression(self, plasma_source): def test_sampling(self, plasma_source): """Test the sampling function.""" randoms = [r() for r in [random] * 8] - sample = plasma_source.sample_source(randoms) + sample = plasma_source.sample(randoms) assert len(sample) == 7 @@ -147,6 +147,6 @@ def test_sampling_regression(self, plasma_source): 14.202333312096737, ) - sample = plasma_source.sample_source(randoms) + sample = plasma_source.sample(randoms) assert sample == expected From e5e021ea394e387f47518b07713c1b5ad404b229 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 14:00:55 +0100 Subject: [PATCH 34/40] Implement public getters for attributes --- parametric_plasma_source/plasma_source.hpp | 100 ++++++++++++++++++ .../plasma_source_pybind.cpp | 20 ++++ 2 files changed, 120 insertions(+) diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index c71cf66..c22f864 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -118,6 +118,106 @@ class PlasmaSource { */ static PlasmaSource from_string(std::string parameters); + /* + * Getter for pedestal ion density + */ + const double &getIonDensityPedestal() const { return ionDensityPedistal; } + + /* + * Getter for separatrix ion density + */ + const double &getIonDensitySeparatrix() const { return ionDensitySeperatrix; } + + /* + * Getter for origin ion density + */ + const double &getIonDensityOrigin() const { return ionDensityOrigin; } + + /* + * Getter for pedestal ion temperature + */ + const double &getIonTemperaturePedestal() const { return ionTemperaturePedistal; } + + /* + * Getter for separatrix ion temperature + */ + const double &getIonTemperatureSeperatrix() const { return ionTemperatureSeperatrix; } + + /* + * Getter for origin ion temperature + */ + const double &getIonTemperatureOrigin() const { return ionTemperatureOrigin; } + + /* + * Getter for pedestal radius + */ + const double &getPedestalRadius() const { return pedistalRadius; } + + /* + * Getter for ion density peaking factor + */ + const double &getIonDensityPeaking() const { return ionDensityPeaking; } + + /* + * Getter for ion temperature peaking factor + */ + const double &getIonTemperaturePeaking() const { return ionTemperaturePeaking; } + + /* + * Getter for ion temperature beta factor + */ + const double &getIonTemperatureBeta() const { return ionTemperatureBeta; } + + /* + * Getter for minor radius + */ + const double &getMinorRadius() const { return minorRadius; } + + /* + * Getter for major radius + */ + const double &getMajorRadius() const { return majorRadius; } + + /* + * Getter for elongation + */ + const double &getElongation() const { return elongation; } + + /* + * Getter for triangularity + */ + const double &getTriangularity() const { return triangularity; } + + /* + * Getter for shafranov shift + */ + const double &getShafranov() const { return shafranov; } + + /* + * Getter for minimum toroidal angle + */ + const double &getMinToroidalAngle() const { return minToroidalAngle; } + + /* + * Getter for maximum toroidal angle + */ + const double &getMaxToroidalAngle() const { return maxToroidalAngle; } + + /* + * Getter for plasma type + */ + const std::string &getPlasmaType() const { return plasmaType; } + + /* + * Getter for plasma identifier + */ + const int &getPlasmaId() const { return plasmaId; } + + /* + * Getter for number of bins + */ + const int &getNumberOfBins() const { return numberOfBins; } + private: std::vector source_profile; std::vector ion_kt; diff --git a/parametric_plasma_source/plasma_source_pybind.cpp b/parametric_plasma_source/plasma_source_pybind.cpp index 3488e84..f278e74 100644 --- a/parametric_plasma_source/plasma_source_pybind.cpp +++ b/parametric_plasma_source/plasma_source_pybind.cpp @@ -34,6 +34,26 @@ PYBIND11_MODULE(plasma_source, m) { py::arg("number_of_bins")=100, py::arg("min_toroidal_angle")=0.0, py::arg("max_toridal_angle")=360.0) + .def_property_readonly("ion_density_pedistal", &ps::PlasmaSource::getIonDensityPedestal) + .def_property_readonly("ion_density_seperatrix", &ps::PlasmaSource::getIonDensitySeparatrix) + .def_property_readonly("ion_density_origin", &ps::PlasmaSource::getIonDensityOrigin) + .def_property_readonly("ion_temperature_pedistal", &ps::PlasmaSource::getIonTemperaturePedestal) + .def_property_readonly("ion_temperature_seperatrix", &ps::PlasmaSource::getIonTemperatureSeperatrix) + .def_property_readonly("ion_temperature_origin", &ps::PlasmaSource::getIonTemperatureOrigin) + .def_property_readonly("pedistal_radius", &ps::PlasmaSource::getPedestalRadius) + .def_property_readonly("ion_density_peaking_factor", &ps::PlasmaSource::getIonDensityPeaking) + .def_property_readonly("ion_temperature_peaking_factor", &ps::PlasmaSource::getIonTemperaturePeaking) + .def_property_readonly("ion_temperature_beta", &ps::PlasmaSource::getIonTemperatureBeta) + .def_property_readonly("minor_radius", &ps::PlasmaSource::getMinorRadius) + .def_property_readonly("major_radius", &ps::PlasmaSource::getMajorRadius) + .def_property_readonly("elongation", &ps::PlasmaSource::getElongation) + .def_property_readonly("triangularity", &ps::PlasmaSource::getTriangularity) + .def_property_readonly("shafranov_shift", &ps::PlasmaSource::getShafranov) + .def_property_readonly("plasma_type", &ps::PlasmaSource::getPlasmaType) + .def_property_readonly("plasma_id", &ps::PlasmaSource::getPlasmaId) + .def_property_readonly("number_of_bins", &ps::PlasmaSource::getNumberOfBins) + .def_property_readonly("min_toroidal_angle", &ps::PlasmaSource::getMinToroidalAngle) + .def_property_readonly("max_toridal_angle", &ps::PlasmaSource::getMaxToroidalAngle) .def("ion_density", &ps::PlasmaSource::ion_density, "Calculate the ion density at a specific minor radius", From 02e2a4d124a60992e728e69d90fd52ce05f3e505 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 14:16:04 +0100 Subject: [PATCH 35/40] Fix tests pytest.approx syntax was wrong (left over from removing numpy isclose). Length units should match the input (cm not m). Move fixure to conftest.py so that it can be shared. --- tests/conftest.py | 28 ++++++++++++++++++ tests/test_plasma_source.py | 59 ++++++++++--------------------------- 2 files changed, 44 insertions(+), 43 deletions(-) create mode 100644 tests/conftest.py diff --git a/tests/conftest.py b/tests/conftest.py new file mode 100644 index 0000000..5520a3f --- /dev/null +++ b/tests/conftest.py @@ -0,0 +1,28 @@ +import pytest + +from parametric_plasma_source import PlasmaSource + +plasma_params = { + "elongation": 1.557, + "ion_density_origin": 1.09e20, + "ion_density_peaking_factor": 1, + "ion_density_pedistal": 1.09e20, + "ion_density_seperatrix": 3e19, + "ion_temperature_origin": 45.9, + "ion_temperature_peaking_factor": 8.06, + "ion_temperature_pedistal": 6.09, + "ion_temperature_seperatrix": 0.1, + "major_radius": 906.0, + "minor_radius": 292.258, + "pedistal_radius": 0.8 * 292.258, + "plasma_id": 1, + "shafranov_shift": 44.789, + "triangularity": 0.270, + "ion_temperature_beta": 6, +} + + +@pytest.fixture(scope="session") +def plasma_source(): + """Make a plasma source to use as a test fixture.""" + return PlasmaSource(**plasma_params) diff --git a/tests/test_plasma_source.py b/tests/test_plasma_source.py index 0eafc51..04df25f 100644 --- a/tests/test_plasma_source.py +++ b/tests/test_plasma_source.py @@ -7,31 +7,6 @@ from parametric_plasma_source import PlasmaSource -plasma_params = { - "elongation": 1.557, - "ion_density_origin": 1.09e20, - "ion_density_peaking_factor": 1, - "ion_density_pedistal": 1.09e20, - "ion_density_seperatrix": 3e19, - "ion_temperature_origin": 45.9, - "ion_temperature_peaking_factor": 8.06, - "ion_temperature_pedistal": 6.09, - "ion_temperature_seperatrix": 0.1, - "major_radius": 906.0, - "minor_radius": 292.258, - "pedistal_radius": 0.8 * 292.258, - "plasma_id": 1, - "shafranov_shift": 44.789, - "triangularity": 0.270, - "ion_temperature_beta": 6, -} - - -@pytest.fixture(scope="session") -def plasma_source(): - """Make a plasma source to use as a test fixture.""" - return PlasmaSource(**plasma_params) - class TestPlasmaSource: """A class to run tests against the plasma source.""" @@ -40,59 +15,57 @@ def test_ion_density_magnetic_origin(self, plasma_source): """Test the ion density at the magnetic origin.""" ion_density = plasma_source.ion_density(0.0) - assert pytest.approx(ion_density, 1.09e20) + assert ion_density == pytest.approx(1.09e20) def test_ion_density_inside_pedestal(self, plasma_source): """Test the ion density inside the pedestal.""" - ion_density = plasma_source.ion_density(0.2) + ion_density = plasma_source.ion_density(20.0) - assert pytest.approx(ion_density, 1.09e20) + assert ion_density == pytest.approx(1.09e20) def test_ion_density_outside_pedestal(self, plasma_source): """Test the ion density outside the pedestal.""" - ion_density = plasma_source.ion_density(2.4) + ion_density = plasma_source.ion_density(240.0) - assert pytest.approx(ion_density, 1.00628584e20) + assert ion_density == pytest.approx(1.00629067e20) def test_ion_density_boundary(self, plasma_source): """Test the ion density at the boundary.""" - boundary = plasma_params["minor_radius"] / 100.0 - ion_density = plasma_source.ion_density(boundary) + ion_density = plasma_source.ion_density(plasma_source.minor_radius) - assert pytest.approx(ion_density, plasma_params["ion_density_seperatrix"]) + assert ion_density == pytest.approx(plasma_source.ion_density_seperatrix) def test_ion_temperature_magnetic_origin(self, plasma_source): """Test the ion temperature at the magnetic origin.""" ion_temperature = plasma_source.ion_temperature(0.0) - assert pytest.approx(ion_temperature, plasma_params["ion_temperature_origin"]) + assert ion_temperature == pytest.approx(plasma_source.ion_temperature_origin) def test_ion_temperature_inside_pedestal(self, plasma_source): """Test the ion temperature inside the pedestal.""" - ion_temperature = plasma_source.ion_temperature(0.2) + ion_temperature = plasma_source.ion_temperature(20.0) - assert pytest.approx(ion_temperature, 45.89987429) + assert ion_temperature == pytest.approx(45.89987429) def test_ion_temperature_outside_pedestal(self, plasma_source): """Test the ion temperature outside the pedestal.""" - ion_temperature = plasma_source.ion_temperature(2.4) + ion_temperature = plasma_source.ion_temperature(240.0) - assert pytest.approx(ion_temperature, 5.45525594) + assert ion_temperature == pytest.approx(5.45529258) def test_ion_temperature_boundary(self, plasma_source): """Test the ion temperature at the boundary.""" - boundary = plasma_params["minor_radius"] / 100.0 - ion_temperature = plasma_source.ion_temperature(boundary) + ion_temperature = plasma_source.ion_temperature(plasma_source.minor_radius) - assert pytest.approx( - ion_temperature, plasma_params["ion_temperature_seperatrix"] + assert ion_temperature == pytest.approx( + plasma_source.ion_temperature_seperatrix ) def test_dt_cross_section(self, plasma_source): """Test the dt cross section at a specific temperature.""" dt_cross_section = plasma_source.dt_xs(4.25e7) - assert pytest.approx(dt_cross_section, 0.0) + assert dt_cross_section == pytest.approx(0.0) def test_source_to_from_string(self, plasma_source): """Test the source can be converted to and from a string.""" From 490f4ef3d03105e2f551a781d92c2287b04f6d9f Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 15:05:38 +0100 Subject: [PATCH 36/40] Small format update --- tests/test_plasma_source.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_plasma_source.py b/tests/test_plasma_source.py index 04df25f..965d35e 100644 --- a/tests/test_plasma_source.py +++ b/tests/test_plasma_source.py @@ -94,7 +94,7 @@ def test_sampling(self, plasma_source): assert len(sample) == 7 # The 4th, 5th and 6th elements together define a unit vector - direction_length = math.sqrt(sample[3]**2 + sample[4]**2 + sample[5]**2) + direction_length = math.sqrt(sample[3] ** 2 + sample[4] ** 2 + sample[5] ** 2) assert direction_length == pytest.approx(1.0) def test_sampling_regression(self, plasma_source): From ec9ae94ae68c94cab44b98ee72f5a6694ecf2799 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 15:28:08 +0100 Subject: [PATCH 37/40] Fix some typos seperatrix -> separatrix pedistal -> pedestal --- README.md | 10 ++-- examples/build_xml.py | 10 ++-- parametric_plasma_source/plasma_source.cpp | 50 +++++++++---------- parametric_plasma_source/plasma_source.hpp | 22 ++++---- .../plasma_source_pybind.cpp | 20 ++++---- tests/conftest.py | 10 ++-- tests/test_plasma_source.py | 4 +- 7 files changed, 63 insertions(+), 63 deletions(-) diff --git a/README.md b/README.md index 84d791e..a1ba215 100644 --- a/README.md +++ b/README.md @@ -28,13 +28,13 @@ In the above example the major_radius, minor_radius, elongation and triangularit There are a number of additional arguments that can be passed to the Plasma class on construction. Units are in SI (e.g. meters not cm) ```[python] -ion_density_pedistal = 1.09e+20 -ion_density_seperatrix = 3e+19 +ion_density_pedestal = 1.09e+20 +ion_density_separatrix = 3e+19 ion_density_origin = 1.09e+20 -ion_temperature_pedistal = 6.09 -ion_temperature_seperatrix = 0.1 +ion_temperature_pedestal = 6.09 +ion_temperature_separatrix = 0.1 ion_temperature_origin = 45.9 -pedistal_radius = 0.8 +pedestal_radius = 0.8 ion_density_peaking_factor = 1 ion_temperature_peaking_factor = 8.06 minor_radius = 1.56 diff --git a/examples/build_xml.py b/examples/build_xml.py index 54a2008..8ee085f 100644 --- a/examples/build_xml.py +++ b/examples/build_xml.py @@ -13,15 +13,15 @@ "elongation": 1.557, "ion_density_origin": 1.09e20, "ion_density_peaking_factor": 1, - "ion_density_pedistal": 1.09e20, - "ion_density_seperatrix": 3e19, + "ion_density_pedestal": 1.09e20, + "ion_density_separatrix": 3e19, "ion_temperature_origin": 45.9, "ion_temperature_peaking_factor": 8.06, - "ion_temperature_pedistal": 6.09, - "ion_temperature_seperatrix": 0.1, + "ion_temperature_pedestal": 6.09, + "ion_temperature_separatrix": 0.1, "major_radius": 906.0, "minor_radius": 292.258, - "pedistal_radius": 0.8 * 292.258, + "pedestal_radius": 0.8 * 292.258, "plasma_id": 1, "shafranov_shift": 44.789, "triangularity": 0.270, diff --git a/parametric_plasma_source/plasma_source.cpp b/parametric_plasma_source/plasma_source.cpp index 027bc50..4af4941 100644 --- a/parametric_plasma_source/plasma_source.cpp +++ b/parametric_plasma_source/plasma_source.cpp @@ -15,7 +15,7 @@ PlasmaSource::PlasmaSource() {} PlasmaSource::PlasmaSource(const double ion_density_ped, const double ion_density_sep, const double ion_density_origin, const double ion_temp_ped, const double ion_temp_sep, const double ion_temp_origin, - const double pedistal_rad, const double ion_density_peak, + const double pedestal_rad, const double ion_density_peak, const double ion_temp_peak, const double ion_temp_beta, const double minor_radius, const double major_radius, const double elongation, const double triangularity, @@ -25,13 +25,13 @@ PlasmaSource::PlasmaSource(const double ion_density_ped, const double ion_densit const double max_toroidal_angle ) { // set params - ionDensityPedistal = ion_density_ped; - ionDensitySeperatrix = ion_density_sep; + ionDensityPedestal = ion_density_ped; + ionDensitySeparatrix = ion_density_sep; ionDensityOrigin = ion_density_origin; - ionTemperaturePedistal = ion_temp_ped; - ionTemperatureSeperatrix = ion_temp_sep; + ionTemperaturePedestal = ion_temp_ped; + ionTemperatureSeparatrix = ion_temp_sep; ionTemperatureOrigin = ion_temp_origin; - pedistalRadius = pedistal_rad; + pedestalRadius = pedestal_rad; ionDensityPeaking = ion_density_peak; ionTemperaturePeaking = ion_temp_peak; ionTemperatureBeta = ion_temp_beta; @@ -196,18 +196,18 @@ double PlasmaSource::ion_density(const double sample_radius) ion_dens = ionDensityOrigin* (1.0-std::pow(sample_radius/minorRadius,2)); } else { - if(sample_radius <= pedistalRadius) { - ion_dens += ionDensityPedistal; + if(sample_radius <= pedestalRadius) { + ion_dens += ionDensityPedestal; double product; - product = 1.0-std::pow(sample_radius/pedistalRadius,2); + product = 1.0-std::pow(sample_radius/pedestalRadius,2); product = std::pow(product,ionDensityPeaking); - ion_dens += (ionDensityOrigin-ionDensityPedistal)* + ion_dens += (ionDensityOrigin-ionDensityPedestal)* (product); } else { - ion_dens += ionDensitySeperatrix; + ion_dens += ionDensitySeparatrix; double product; - product = ionDensityPedistal - ionDensitySeperatrix; - ion_dens += product*(minorRadius-sample_radius)/(minorRadius-pedistalRadius); + product = ionDensityPedestal - ionDensitySeparatrix; + ion_dens += product*(minorRadius-sample_radius)/(minorRadius-pedestalRadius); } } @@ -227,18 +227,18 @@ double PlasmaSource::ion_temperature(const double sample_radius) (1.0-std::pow(sample_radius/minorRadius, ionTemperaturePeaking)); } else { - if(sample_radius <= pedistalRadius) { - ion_temp += ionTemperaturePedistal; + if(sample_radius <= pedestalRadius) { + ion_temp += ionTemperaturePedestal; double product; - product = 1.0-std::pow(sample_radius/pedistalRadius,ionTemperatureBeta); + product = 1.0-std::pow(sample_radius/pedestalRadius,ionTemperatureBeta); product = std::pow(product,ionTemperaturePeaking); ion_temp += (ionTemperatureOrigin- - ionTemperaturePedistal)*(product); + ionTemperaturePedestal)*(product); } else { - ion_temp += ionTemperatureSeperatrix; + ion_temp += ionTemperatureSeparatrix; double product; - product = ionTemperaturePedistal - ionTemperatureSeperatrix; - ion_temp += product*(minorRadius-sample_radius)/(minorRadius-pedistalRadius); + product = ionTemperaturePedestal - ionTemperatureSeparatrix; + ion_temp += product*(minorRadius-sample_radius)/(minorRadius-pedestalRadius); } } @@ -288,13 +288,13 @@ std::string PlasmaSource::to_string() { result << "elongation=" << elongation << ", "; result << "triangularity=" << triangularity << ", "; result << "shafranov_shift=" << shafranov << ", "; - result << "pedestal_radius=" << pedistalRadius << ", "; - result << "ion_density_pedestal=" << ionDensityPedistal << ", "; - result << "ion_density_separatrix=" << ionDensitySeperatrix << ", "; + result << "pedestal_radius=" << pedestalRadius << ", "; + result << "ion_density_pedestal=" << ionDensityPedestal << ", "; + result << "ion_density_separatrix=" << ionDensitySeparatrix << ", "; result << "ion_density_origin=" << ionDensityOrigin << ", "; result << "ion_density_alpha=" << ionDensityPeaking << ", "; - result << "ion_temperature_pedestal=" << ionTemperaturePedistal << ", "; - result << "ion_temperature_separatrix=" << ionTemperatureSeperatrix << ", "; + result << "ion_temperature_pedestal=" << ionTemperaturePedestal << ", "; + result << "ion_temperature_separatrix=" << ionTemperatureSeparatrix << ", "; result << "ion_temperature_origin=" << ionTemperatureOrigin << ", "; result << "ion_temperature_alpha=" << ionTemperaturePeaking << ", "; result << "ion_temperature_beta=" << ionTemperatureBeta << ", "; diff --git a/parametric_plasma_source/plasma_source.hpp b/parametric_plasma_source/plasma_source.hpp index c22f864..8f81148 100644 --- a/parametric_plasma_source/plasma_source.hpp +++ b/parametric_plasma_source/plasma_source.hpp @@ -22,7 +22,7 @@ class PlasmaSource { const double ion_temp_ped, const double ion_temp_sep, const double ion_temp_origin, - const double pedistal_rad, + const double pedestal_rad, const double ion_density_peak, const double ion_temp_peak, const double ion_temp_beta, @@ -121,12 +121,12 @@ class PlasmaSource { /* * Getter for pedestal ion density */ - const double &getIonDensityPedestal() const { return ionDensityPedistal; } + const double &getIonDensityPedestal() const { return ionDensityPedestal; } /* * Getter for separatrix ion density */ - const double &getIonDensitySeparatrix() const { return ionDensitySeperatrix; } + const double &getIonDensitySeparatrix() const { return ionDensitySeparatrix; } /* * Getter for origin ion density @@ -136,12 +136,12 @@ class PlasmaSource { /* * Getter for pedestal ion temperature */ - const double &getIonTemperaturePedestal() const { return ionTemperaturePedistal; } + const double &getIonTemperaturePedestal() const { return ionTemperaturePedestal; } /* * Getter for separatrix ion temperature */ - const double &getIonTemperatureSeperatrix() const { return ionTemperatureSeperatrix; } + const double &getIonTemperatureSeparatrix() const { return ionTemperatureSeparatrix; } /* * Getter for origin ion temperature @@ -151,7 +151,7 @@ class PlasmaSource { /* * Getter for pedestal radius */ - const double &getPedestalRadius() const { return pedistalRadius; } + const double &getPedestalRadius() const { return pedestalRadius; } /* * Getter for ion density peaking factor @@ -222,13 +222,13 @@ class PlasmaSource { std::vector source_profile; std::vector ion_kt; - double ionDensityPedistal; - double ionDensitySeperatrix; + double ionDensityPedestal; + double ionDensitySeparatrix; double ionDensityOrigin; - double ionTemperaturePedistal; - double ionTemperatureSeperatrix; + double ionTemperaturePedestal; + double ionTemperatureSeparatrix; double ionTemperatureOrigin; - double pedistalRadius; + double pedestalRadius; double ionDensityPeaking; double ionTemperaturePeaking; double ionTemperatureBeta; diff --git a/parametric_plasma_source/plasma_source_pybind.cpp b/parametric_plasma_source/plasma_source_pybind.cpp index f278e74..218e68e 100644 --- a/parametric_plasma_source/plasma_source_pybind.cpp +++ b/parametric_plasma_source/plasma_source_pybind.cpp @@ -14,13 +14,13 @@ PYBIND11_MODULE(plasma_source, m) { const double &, const double &, const double &, const double &, const double &, const double &, const double &, const std::string &, const int &, const int &, const double &, const double &>(), - py::arg("ion_density_pedistal"), - py::arg("ion_density_seperatrix"), + py::arg("ion_density_pedestal"), + py::arg("ion_density_separatrix"), py::arg("ion_density_origin"), - py::arg("ion_temperature_pedistal"), - py::arg("ion_temperature_seperatrix"), + py::arg("ion_temperature_pedestal"), + py::arg("ion_temperature_separatrix"), py::arg("ion_temperature_origin"), - py::arg("pedistal_radius"), + py::arg("pedestal_radius"), py::arg("ion_density_peaking_factor"), py::arg("ion_temperature_peaking_factor"), py::arg("ion_temperature_beta"), @@ -34,13 +34,13 @@ PYBIND11_MODULE(plasma_source, m) { py::arg("number_of_bins")=100, py::arg("min_toroidal_angle")=0.0, py::arg("max_toridal_angle")=360.0) - .def_property_readonly("ion_density_pedistal", &ps::PlasmaSource::getIonDensityPedestal) - .def_property_readonly("ion_density_seperatrix", &ps::PlasmaSource::getIonDensitySeparatrix) + .def_property_readonly("ion_density_pedestal", &ps::PlasmaSource::getIonDensityPedestal) + .def_property_readonly("ion_density_separatrix", &ps::PlasmaSource::getIonDensitySeparatrix) .def_property_readonly("ion_density_origin", &ps::PlasmaSource::getIonDensityOrigin) - .def_property_readonly("ion_temperature_pedistal", &ps::PlasmaSource::getIonTemperaturePedestal) - .def_property_readonly("ion_temperature_seperatrix", &ps::PlasmaSource::getIonTemperatureSeperatrix) + .def_property_readonly("ion_temperature_pedestal", &ps::PlasmaSource::getIonTemperaturePedestal) + .def_property_readonly("ion_temperature_separatrix", &ps::PlasmaSource::getIonTemperatureSeparatrix) .def_property_readonly("ion_temperature_origin", &ps::PlasmaSource::getIonTemperatureOrigin) - .def_property_readonly("pedistal_radius", &ps::PlasmaSource::getPedestalRadius) + .def_property_readonly("pedestal_radius", &ps::PlasmaSource::getPedestalRadius) .def_property_readonly("ion_density_peaking_factor", &ps::PlasmaSource::getIonDensityPeaking) .def_property_readonly("ion_temperature_peaking_factor", &ps::PlasmaSource::getIonTemperaturePeaking) .def_property_readonly("ion_temperature_beta", &ps::PlasmaSource::getIonTemperatureBeta) diff --git a/tests/conftest.py b/tests/conftest.py index 5520a3f..2bc88fa 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -6,15 +6,15 @@ "elongation": 1.557, "ion_density_origin": 1.09e20, "ion_density_peaking_factor": 1, - "ion_density_pedistal": 1.09e20, - "ion_density_seperatrix": 3e19, + "ion_density_pedestal": 1.09e20, + "ion_density_separatrix": 3e19, "ion_temperature_origin": 45.9, "ion_temperature_peaking_factor": 8.06, - "ion_temperature_pedistal": 6.09, - "ion_temperature_seperatrix": 0.1, + "ion_temperature_pedestal": 6.09, + "ion_temperature_separatrix": 0.1, "major_radius": 906.0, "minor_radius": 292.258, - "pedistal_radius": 0.8 * 292.258, + "pedestal_radius": 0.8 * 292.258, "plasma_id": 1, "shafranov_shift": 44.789, "triangularity": 0.270, diff --git a/tests/test_plasma_source.py b/tests/test_plasma_source.py index 965d35e..c208df9 100644 --- a/tests/test_plasma_source.py +++ b/tests/test_plasma_source.py @@ -33,7 +33,7 @@ def test_ion_density_boundary(self, plasma_source): """Test the ion density at the boundary.""" ion_density = plasma_source.ion_density(plasma_source.minor_radius) - assert ion_density == pytest.approx(plasma_source.ion_density_seperatrix) + assert ion_density == pytest.approx(plasma_source.ion_density_separatrix) def test_ion_temperature_magnetic_origin(self, plasma_source): """Test the ion temperature at the magnetic origin.""" @@ -58,7 +58,7 @@ def test_ion_temperature_boundary(self, plasma_source): ion_temperature = plasma_source.ion_temperature(plasma_source.minor_radius) assert ion_temperature == pytest.approx( - plasma_source.ion_temperature_seperatrix + plasma_source.ion_temperature_separatrix ) def test_dt_cross_section(self, plasma_source): From 6a4c0d9266cec4f816cc2ad785853cef3f05f7a3 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 16:12:41 +0100 Subject: [PATCH 38/40] Update README --- README.md | 114 +++++++++++++++++++++++++++++++++++++++++------------- 1 file changed, 88 insertions(+), 26 deletions(-) diff --git a/README.md b/README.md index a1ba215..185e9d4 100644 --- a/README.md +++ b/README.md @@ -8,42 +8,104 @@ The plasma source is based on a paper by [C. Fausser et al](https://www.scienced ## Installation +### Installing from PyPI + ```pip install parametric_plasma_source``` +### Installing from source + +Installation of the parametric plasma source from source requires cmake to build the underlying C++ code. This can be obtained from +your OS's package manager by e.g. `sudo apt-get install cmake` or from cmake source. + +If you intend to develop the code then it is recommended to work in a virtual environment. + +The requirements for developing the code can be installed by running: + +```pip install -r requirements-develop.txt``` + +The package can be built and installed in editable mode by: + +```pip install -e .``` + ## Usage -The parametric plasma source can be imported an used in Python 3 in the following manner. +The parametric plasma source can be sampled either directly in Python 3 or sampled in an OpenMC simulation. + +For a better understanding of the varibles take a look at the [C. Fausser et al](https://www.sciencedirect.com/science/article/pii/S0920379612000853) paper. + +### Sampling in Python + +The parametric plasma source can be imported an used in Python 3 in the following manner: ```[python] -from parametric_plasma_source import Plasma -my_plasma = Plasma(major_radius=6, - minor_radius=1.5, - elongation = 2.0 - triangularity = 0.55) -my_plasma.export_plasma_source('custom_openmc_plasma_source.so') +from parametric_plasma_source import PlasmaSource +from random import random + +plasma_params = { + "elongation": 1.557, + "ion_density_origin": 1.09e20, + "ion_density_peaking_factor": 1, + "ion_density_pedestal": 1.09e20, + "ion_density_separatrix": 3e19, + "ion_temperature_origin": 45.9, + "ion_temperature_peaking_factor": 8.06, + "ion_temperature_pedestal": 6.09, + "ion_temperature_separatrix": 0.1, + "major_radius": 906.0, + "minor_radius": 292.258, + "pedestal_radius": 0.8 * 292.258, + "plasma_id": 1, + "shafranov_shift": 44.789, + "triangularity": 0.270, + "ion_temperature_beta": 6, +} + +my_plasma = PlasmaSource(**plasma_params) +sample = my_plasma.sample([random(), random(), random(), random(), random(), random(), random(), random()]) +particle_x, particle_y, particle_z = sample[0], sample[1], sample[2] +particle_x_dir, particle_y_dir, particle_z_dir = sample[3], sample[4], sample[5] +particle_energy_mev = sample[6] ``` -In the above example the major_radius, minor_radius, elongation and triangularity while the other varibles are kept as the default values. +### Sampling in OpenMC -There are a number of additional arguments that can be passed to the Plasma class on construction. Units are in SI (e.g. meters not cm) +The parametric plasma source also contains a plugin library for OpenMC to allow the source to be sampled in an OpenMC simulation. ```[python] -ion_density_pedestal = 1.09e+20 -ion_density_separatrix = 3e+19 -ion_density_origin = 1.09e+20 -ion_temperature_pedestal = 6.09 -ion_temperature_separatrix = 0.1 -ion_temperature_origin = 45.9 -pedestal_radius = 0.8 -ion_density_peaking_factor = 1 -ion_temperature_peaking_factor = 8.06 -minor_radius = 1.56 -major_radius = 2.5 -elongation = 2.0 -triangularity = 0.55 -shafranov_shift = 0.0 -number_of_bins = 100 -plasma_type = 1 +from parametric_plasma_source import PlasmaSource, SOURCE_SAMPLING_PATH +import openmc + +plasma_params = { + "elongation": 1.557, + "ion_density_origin": 1.09e20, + "ion_density_peaking_factor": 1, + "ion_density_pedestal": 1.09e20, + "ion_density_separatrix": 3e19, + "ion_temperature_origin": 45.9, + "ion_temperature_peaking_factor": 8.06, + "ion_temperature_pedestal": 6.09, + "ion_temperature_separatrix": 0.1, + "major_radius": 906.0, + "minor_radius": 292.258, + "pedestal_radius": 0.8 * 292.258, + "plasma_id": 1, + "shafranov_shift": 44.789, + "triangularity": 0.270, + "ion_temperature_beta": 6, +} + +my_plasma = PlasmaSource(**plasma_params) +settings = openmc.Settings() +settings.run_mode = "fixed source" +settings.batches = 10 +settings.particles = 1000 +source = openmc.Source() +source.library = SOURCE_SAMPLING_PATH +source.parameters = str(my_plasma) +settings.source = source +settings.export_to_xml() ``` -For a better understanding of the varibles take a look at the [C. Fausser et al](https://www.sciencedirect.com/science/article/pii/S0920379612000853) paper. +## Running Tests + +The tests are run by executing `pytest tests` from within your virtual environment. From c0b0148467f0d3a8d7b780f566c740fa91f05b03 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 19:57:01 +0100 Subject: [PATCH 39/40] Supply length units in m Reinstate conversion to cm when sampling in OpenMC. --- README.md | 18 ++++++++++-------- examples/build_xml.py | 8 ++++---- parametric_plasma_source/source_sampling.cpp | 5 +++++ tests/conftest.py | 8 ++++---- tests/test_data/baseline_source.txt | 2 +- tests/test_plasma_source.py | 16 ++++++++-------- 6 files changed, 32 insertions(+), 25 deletions(-) diff --git a/README.md b/README.md index 185e9d4..4cc5517 100644 --- a/README.md +++ b/README.md @@ -51,11 +51,11 @@ plasma_params = { "ion_temperature_peaking_factor": 8.06, "ion_temperature_pedestal": 6.09, "ion_temperature_separatrix": 0.1, - "major_radius": 906.0, - "minor_radius": 292.258, - "pedestal_radius": 0.8 * 292.258, + "major_radius": 9.06, + "minor_radius": 2.92258, + "pedestal_radius": 0.8 * 2.92258, "plasma_id": 1, - "shafranov_shift": 44.789, + "shafranov_shift": 0.44789, "triangularity": 0.270, "ion_temperature_beta": 6, } @@ -71,6 +71,8 @@ particle_energy_mev = sample[6] The parametric plasma source also contains a plugin library for OpenMC to allow the source to be sampled in an OpenMC simulation. +When using the OpenMC sampling the inputs must be provided in meters where applicable (the sampling will convert to cm). + ```[python] from parametric_plasma_source import PlasmaSource, SOURCE_SAMPLING_PATH import openmc @@ -85,11 +87,11 @@ plasma_params = { "ion_temperature_peaking_factor": 8.06, "ion_temperature_pedestal": 6.09, "ion_temperature_separatrix": 0.1, - "major_radius": 906.0, - "minor_radius": 292.258, - "pedestal_radius": 0.8 * 292.258, + "major_radius": 9.06, + "minor_radius": 2.92258, + "pedestal_radius": 0.8 * 2.92258, "plasma_id": 1, - "shafranov_shift": 44.789, + "shafranov_shift": 0.44789, "triangularity": 0.270, "ion_temperature_beta": 6, } diff --git a/examples/build_xml.py b/examples/build_xml.py index 8ee085f..477ca4b 100644 --- a/examples/build_xml.py +++ b/examples/build_xml.py @@ -19,11 +19,11 @@ "ion_temperature_peaking_factor": 8.06, "ion_temperature_pedestal": 6.09, "ion_temperature_separatrix": 0.1, - "major_radius": 906.0, - "minor_radius": 292.258, - "pedestal_radius": 0.8 * 292.258, + "major_radius": 9.06, + "minor_radius": 2.92258, + "pedestal_radius": 0.8 * 2.92258, "plasma_id": 1, - "shafranov_shift": 44.789, + "shafranov_shift": 0.44789, "triangularity": 0.270, "ion_temperature_beta": 6, } diff --git a/parametric_plasma_source/source_sampling.cpp b/parametric_plasma_source/source_sampling.cpp index fbd2cb3..49b71a5 100644 --- a/parametric_plasma_source/source_sampling.cpp +++ b/parametric_plasma_source/source_sampling.cpp @@ -42,6 +42,11 @@ class SampledSource : public openmc::CustomSource { particle.wgt = 1.0; particle.delayed_group = 0; + // position (note units converted m -> cm) + particle.r.x *= 100.; + particle.r.y *= 100.; + particle.r.z *= 100.; + // energy particle.E = E * 1e6; // convert from MeV -> eV diff --git a/tests/conftest.py b/tests/conftest.py index 2bc88fa..4b8a656 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -12,11 +12,11 @@ "ion_temperature_peaking_factor": 8.06, "ion_temperature_pedestal": 6.09, "ion_temperature_separatrix": 0.1, - "major_radius": 906.0, - "minor_radius": 292.258, - "pedestal_radius": 0.8 * 292.258, + "major_radius": 9.06, + "minor_radius": 2.92258, + "pedestal_radius": 0.8 * 2.92258, "plasma_id": 1, - "shafranov_shift": 44.789, + "shafranov_shift": 0.44789, "triangularity": 0.270, "ion_temperature_beta": 6, } diff --git a/tests/test_data/baseline_source.txt b/tests/test_data/baseline_source.txt index f33a26d..e249ef4 100644 --- a/tests/test_data/baseline_source.txt +++ b/tests/test_data/baseline_source.txt @@ -1 +1 @@ -major_radius=906, minor_radius=292.258, elongation=1.557, triangularity=0.27, shafranov_shift=44.789, pedestal_radius=233.806, ion_density_pedestal=1.09e+20, ion_density_separatrix=3e+19, ion_density_origin=1.09e+20, ion_density_alpha=1, ion_temperature_pedestal=6.09, ion_temperature_separatrix=0.1, ion_temperature_origin=45.9, ion_temperature_alpha=8.06, ion_temperature_beta=6, plasma_type=plasma, plasma_id=1, number_of_bins=100, minimum_toroidal_angle=0, maximum_toroidal_angle=360 \ No newline at end of file +major_radius=9.06, minor_radius=2.92258, elongation=1.557, triangularity=0.27, shafranov_shift=0.44789, pedestal_radius=2.33806, ion_density_pedestal=1.09e+20, ion_density_separatrix=3e+19, ion_density_origin=1.09e+20, ion_density_alpha=1, ion_temperature_pedestal=6.09, ion_temperature_separatrix=0.1, ion_temperature_origin=45.9, ion_temperature_alpha=8.06, ion_temperature_beta=6, plasma_type=plasma, plasma_id=1, number_of_bins=100, minimum_toroidal_angle=0, maximum_toroidal_angle=360 \ No newline at end of file diff --git a/tests/test_plasma_source.py b/tests/test_plasma_source.py index c208df9..88cb045 100644 --- a/tests/test_plasma_source.py +++ b/tests/test_plasma_source.py @@ -19,13 +19,13 @@ def test_ion_density_magnetic_origin(self, plasma_source): def test_ion_density_inside_pedestal(self, plasma_source): """Test the ion density inside the pedestal.""" - ion_density = plasma_source.ion_density(20.0) + ion_density = plasma_source.ion_density(0.2) assert ion_density == pytest.approx(1.09e20) def test_ion_density_outside_pedestal(self, plasma_source): """Test the ion density outside the pedestal.""" - ion_density = plasma_source.ion_density(240.0) + ion_density = plasma_source.ion_density(2.4) assert ion_density == pytest.approx(1.00629067e20) @@ -43,13 +43,13 @@ def test_ion_temperature_magnetic_origin(self, plasma_source): def test_ion_temperature_inside_pedestal(self, plasma_source): """Test the ion temperature inside the pedestal.""" - ion_temperature = plasma_source.ion_temperature(20.0) + ion_temperature = plasma_source.ion_temperature(0.2) assert ion_temperature == pytest.approx(45.89987429) def test_ion_temperature_outside_pedestal(self, plasma_source): """Test the ion temperature outside the pedestal.""" - ion_temperature = plasma_source.ion_temperature(240.0) + ion_temperature = plasma_source.ion_temperature(2.4) assert ion_temperature == pytest.approx(5.45529258) @@ -111,9 +111,9 @@ def test_sampling_regression(self, plasma_source): ] expected = ( - 490.1585757452634, - 785.7624748651705, - -19.32184336005464, + 4.9015857574526, + 7.8576247486517, + -0.19321843360054, -0.5702309715680232, -0.06740484811110535, 0.8187143735856279, @@ -122,4 +122,4 @@ def test_sampling_regression(self, plasma_source): sample = plasma_source.sample(randoms) - assert sample == expected + assert sample == pytest.approx(expected) From 0e34dc6c81a69a43319a8c808bf646acef745579 Mon Sep 17 00:00:00 2001 From: Dan Short Date: Wed, 9 Sep 2020 20:03:15 +0100 Subject: [PATCH 40/40] Use openmc model for examples --- examples/{build_xml.py => build_and_run_model.py} | 10 ++++++---- examples/plot_flux.py | 2 +- 2 files changed, 7 insertions(+), 5 deletions(-) rename examples/{build_xml.py => build_and_run_model.py} (94%) diff --git a/examples/build_xml.py b/examples/build_and_run_model.py similarity index 94% rename from examples/build_xml.py rename to examples/build_and_run_model.py index 477ca4b..8cf8516 100644 --- a/examples/build_xml.py +++ b/examples/build_and_run_model.py @@ -35,7 +35,6 @@ iron.set_density("g/cm3", 5.0) iron.add_element("Fe", 1.0) mats = openmc.Materials([iron]) -mats.export_to_xml() # Create a 5 cm x 5 cm box filled with iron cells = [] @@ -46,7 +45,6 @@ cells += [openmc.Cell(fill=None, region=+inner_box1 & -inner_box2)] cells += [openmc.Cell(fill=iron, region=+inner_box2 & outer_box)] geometry = openmc.Geometry(cells) -geometry.export_to_xml() # Tell OpenMC we're going to use our custom source settings = openmc.Settings() @@ -57,7 +55,6 @@ source.library = SOURCE_SAMPLING_PATH source.parameters = str(plasma) settings.source = source -settings.export_to_xml() # Finally, define a mesh tally so that we can see the resulting flux mesh = openmc.RegularMesh() @@ -69,4 +66,9 @@ tally.filters = [openmc.MeshFilter(mesh)] tally.scores = ["flux"] tallies = openmc.Tallies([tally]) -tallies.export_to_xml() + +model = openmc.model.Model( + materials=mats, geometry=geometry, settings=settings, tallies=tallies +) + +model.run() diff --git a/examples/plot_flux.py b/examples/plot_flux.py index 32ee02b..13381b8 100644 --- a/examples/plot_flux.py +++ b/examples/plot_flux.py @@ -7,7 +7,7 @@ # Get the flux from the statepoint with openmc.StatePoint("statepoint.10.h5") as sp: - flux = np.log(sp.tallies[1].mean) + flux = np.log(sp.get_tally(scores=["flux"]).mean) flux.shape = (50, 50) # Plot the flux