From debc085cf692cb4f2a8a9585d0acfe4a7b06bfb9 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 22 Jun 2022 22:50:07 +0200 Subject: [PATCH 01/62] add placeholder skeleton for MCPL-input --- include/openmc/source.h | 19 +++++++++++++++++++ src/source.cpp | 21 +++++++++++++++++++++ 2 files changed, 40 insertions(+) diff --git a/include/openmc/source.h b/include/openmc/source.h index 18fddc6899..cadbfe93bb 100644 --- a/include/openmc/source.h +++ b/include/openmc/source.h @@ -100,6 +100,25 @@ private: vector sites_; //!< Source sites from a file }; +//============================================================================== +// MCPL-file input source +//============================================================================== +class MCPLFileSource : public Source { +public: + // Constructors, destructors + MCPLFileSource(std::string path); + ~MCPLFileSource(); + + // Defer implementation to custom source library + SourceSite sample(uint64_t* seed) const override + +private: + vector sites_; //! Date: Thu, 23 Jun 2022 11:34:53 +0200 Subject: [PATCH 02/62] include mcpl header --- src/source.cpp | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/src/source.cpp b/src/source.cpp index cd8b5e5743..5e26fa0a8d 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -31,6 +31,10 @@ #include "openmc/state_point.h" #include "openmc/xml_interface.h" +#ifdef HAS_MCPL +#include +#endif + namespace openmc { //============================================================================== From 5d74af4bf766ab732d320081be47b9c4f78303fe Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Thu, 23 Jun 2022 11:35:46 +0200 Subject: [PATCH 03/62] declarations of MCPL-internals --- include/openmc/source.h | 13 +++++++++++-- 1 file changed, 11 insertions(+), 2 deletions(-) diff --git a/include/openmc/source.h b/include/openmc/source.h index cadbfe93bb..69dd47eae9 100644 --- a/include/openmc/source.h +++ b/include/openmc/source.h @@ -109,11 +109,20 @@ public: MCPLFileSource(std::string path); ~MCPLFileSource(); - // Defer implementation to custom source library - SourceSite sample(uint64_t* seed) const override + // Properties + ParticleType particle_type() const { return particle_; } + + //! Sample from the external source distribution + //! \param[inout] seed Pseudorandom seed pointer + //! \return Site read from MCPL-file + SourceSite sample(uint64_t* seed) const override; private: + ParticleType particle_ {ParticleType::neutron}; //!< Type of particle emitted vector sites_; //! Date: Thu, 23 Jun 2022 11:37:05 +0200 Subject: [PATCH 04/62] internals of sample and constructor --- src/source.cpp | 51 ++++++++++++++++++++++++++++++++++++++++++++++---- 1 file changed, 47 insertions(+), 4 deletions(-) diff --git a/src/source.cpp b/src/source.cpp index 5e26fa0a8d..231f355fc2 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -337,7 +337,7 @@ CustomSourceWrapper::~CustomSourceWrapper() } //=========================================================================== -// Read an MCPL-file +// Read particles from an MCPL-file //=========================================================================== MCPLFileSource::MCPLFileSource(std::string path) { @@ -348,12 +348,55 @@ MCPLFileSource::MCPLFileSource(std::string path) // Read the source from a binary file instead of sampling from some // assumed source distribution - write_message(6, "Reading source file from {}...", path); + write_message(6, "Reading mcpl source file from {}",path); // Open the mcpl file - id_t file_id = mcpl_open(path.c_str) - + mcpl_file = mcpl_open(path.c_str); + //do checks on the mcpl_file to see if particles are many enough. + // should model this on the example source shown in the docs. + uint64_t nparticles=mcpl_hdr_nparticles(mcpl_file); + +} + +MCPLFileSource::~MCPLFilsSource(){ + mcpl_close_file(mcpl_file); +} + +SourceSite MCPLFileSource::sample(uint64_t *seed){ + SourceSite omc_particle; + mcpl_particle *mcpl_particle; + //extract particle from mcpl-file + mcpl_particle=mcpl_reazd(mcplfile); + // check if it is a neutron or a photon. otherwise skip + while ( mcpl_particle->pdgcode!=2112 && mcpl_particle->pdgcode!=22 ) + { + mcpl_particle(mcpl_read(mcplfile); + //check for file exhaustion + } + + if(mcpl_particle->pdgcode=2112) { + omc_particle_=ParticleType::neutron; + } else { + omc_particle_=Particletype::photon; + } + + //particle is good, convert to openmc-formalism + omc_particle.r.x=mcpl_particle.x; + omc_particle.r.y=mcpl_particle.y; + omc_particle.r.z=mcpl_particle.z; + + omc_particle.u.x=mcpl_particle->direction[0]; + omc_particle.u.y=mcpl_particle->direction[1]; + omc_particle.u.z=mcpl_particle->direction[2]; + + //mcpl stores particles in MeV + omc_particle.E=mcpl_particle->ekin*1e6; + + omc_particle.t=mcpl_particle->time*1e-3; + omc_particle.wgt=mcpl_particle->weight; + omc_particle.delayed_group=0; + return omc_particle; } From 909b96f31e90a8a8ce40e890b3f4e73ee516466b Mon Sep 17 00:00:00 2001 From: erkn Date: Mon, 27 Jun 2022 23:25:35 +0200 Subject: [PATCH 05/62] (wip): skeleton for input only test --- .../regression_tests/source_mcpl_file/test.py | 52 +++++++++++++++++++ 1 file changed, 52 insertions(+) create mode 100644 tests/regression_tests/source_mcpl_file/test.py diff --git a/tests/regression_tests/source_mcpl_file/test.py b/tests/regression_tests/source_mcpl_file/test.py new file mode 100644 index 0000000000..5a6a87d878 --- /dev/null +++ b/tests/regression_tests/source_mcpl_file/test.py @@ -0,0 +1,52 @@ +import pathlib as pl +import os +import shutil +import subprocess +import textwrap + +import openmc +import pytest + +from tests.testing_harness import TestHarness + +def model(): + subprocess.run(['gcc','-o gen_dummt_mcpl.out','-lm','-lmcpl','gen_dummy_mcpl.c'] ) + subprocess.run(['./gen_dummy_mcpl.out']) + +class SourceMCPLFileTestHarness(TestHarness): + def execute_test(self): + """Run OpenMC with the appropriate arguments and check the outputs.""" + try: + self._create_input() + self._run_openmc() + self._test_output_created() + results = self._get_results() + self._write_results(results) + self._compare_results() + finally: + self._cleanup() + + def update_results(self): + """Update the results_true using the current version of OpenMC.""" + try: + self._create_input() + self._run_openmc() + self._test_output_created() + results = self._get_results() + self._write_results(results) + self._overwrite_results() + finally: + self._cleanup() + + def _test_output_created(self): + """Check that the output files were created""" + stat epoint = glob.glob(os.path.join(os.getcwd(), self._sp_name)) + assert len(statepoint) == 1, 'Either multiple or no statepoint files ' \ + 'exist.' + assert statepoint[0].endswith('h5'), \ + 'statepoint file does not end with h5.' + +def test_mcpl_source_file(): + harness = SourceMCPLFileTestHarness('statepoint.10.h5') + harness.main() + From 3c30dbba5b643e8464eb9b45a8667e08bf6ad44b Mon Sep 17 00:00:00 2001 From: erkn Date: Mon, 27 Jun 2022 23:26:01 +0200 Subject: [PATCH 06/62] c-prog to generate a dummy mcpl --- .../source_mcpl_file/gen_dummy_mcpl.c | 97 +++++++++++++++++++ 1 file changed, 97 insertions(+) create mode 100644 tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c diff --git a/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c b/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c new file mode 100644 index 0000000000..7b5a933f18 --- /dev/null +++ b/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c @@ -0,0 +1,97 @@ +#include +#include +#include +#include + +const int batches=10; +const int particles=10000; + +double rand01(){ + double r=rand()/((double) RAND_MAX); + return r; +} + + +int main(int argc, char **argv){ + char outfilename[128]; + snprintf(outfilename,127,"source.%d.mcpl",batches); + + /*generate an mcpl_header_of sorts*/ + mcpl_outfile_t outfile=mcpl_create_outfile(outfilename); + mcpl_particle_t *particle, Particle; + + char line[256]; + snprintf(line,255,"gen_dummy_mcpl.c"); + mcpl_hdr_set_srcname(outfile,line); + mcpl_enable_universal_pdgcode(outfile,2112);/*all particles are neutrons*/ + snprintf(line,255,"Dummy MCPL-file for testing OpenMC mcpl input"); + mcpl_hdr_add_comment(outfile,line); + + /*also add the instrument file and the command line as blobs*/ + FILE *fp; + if( (fp=fopen("gen_dummy_mcpl.c","rb"))!=NULL){ + unsigned char *buffer; + int size,status; + /*find the file size by seeking to end, "tell" the position, and then go back again*/ + fseek(fp, 0L, SEEK_END); + size = ftell(fp); // get current file pointer + fseek(fp, 0L, SEEK_SET); // seek back to beginning of file + if ( size && (buffer=malloc(size))!=NULL){ + if (size!=(fread(buffer,1,size,fp))){ + fprintf(stderr,"Warning: c source generator file not read cleanly\n"); + } + mcpl_hdr_add_data(outfile, "mcpl_point_source_file_generator", size, buffer); + free(buffer); + } + fclose(fp); + } else { + fprintf(stderr,"Warning: could not open c source generator file, hence not embedded.\n"); + } + + + /*the main particle loop*/ + particle=&Particle; + int i; + for (i=0;iposition[0]=0; + particle->position[1]=0; + particle->position[2]=0; + + /*generate a random direction on unit sphere*/ + double nrm=2.0; + double vx,vy,vz; + int iter=0; + do { + if(iter>100){ + printf("warning: exceeed max random vector attempts: at loop %d -> %g %d\n",i,nrm,iter); + } + vx=rand01()*2.0-1.0; + vy=rand01()*2.0-1.0; + vz=rand01()*2.0-1.0; + nrm=(vx*vx + vy*vy + vz*vz); + iter++; + } while (nrm>1.0); + vx/=sqrt(nrm); + vy/=sqrt(nrm); + vz/=sqrt(nrm); + particle->direction[0]=vx; + particle->direction[1]=vy; + particle->direction[2]=vz; + + /*generate the kinetic energy (in MeV) of the particle*/ + particle->ekin=1e-3; + + /*set time=0*/ + particle->time=0; + /*set weight*/ + particle->weight=1; + + particle->userflags=0; + mcpl_add_particle(outfile,particle); + } + mcpl_closeandgzip_outfile(outfile); +} From bc797c4df83b906bebe454cd1c9d7a17bd52c49c Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 5 Jul 2022 01:30:07 +0200 Subject: [PATCH 07/62] don't zip the output-file --- .../source_mcpl_file/gen_dummy_mcpl.c | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c b/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c index 7b5a933f18..f35a07fc44 100644 --- a/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c +++ b/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c @@ -54,20 +54,20 @@ int main(int argc, char **argv){ int i; for (i=0;iposition[0]=0; particle->position[1]=0; particle->position[2]=0; - + /*generate a random direction on unit sphere*/ double nrm=2.0; double vx,vy,vz; int iter=0; do { if(iter>100){ - printf("warning: exceeed max random vector attempts: at loop %d -> %g %d\n",i,nrm,iter); + printf("Warning: Exceeded max random vector attempts: at loop %d -> %g %d\n",i,nrm,iter); } vx=rand01()*2.0-1.0; vy=rand01()*2.0-1.0; @@ -82,9 +82,9 @@ int main(int argc, char **argv){ particle->direction[1]=vy; particle->direction[2]=vz; - /*generate the kinetic energy (in MeV) of the particle*/ + /*generate the kinetic energy (in MeV) of the particle*/ particle->ekin=1e-3; - + /*set time=0*/ particle->time=0; /*set weight*/ @@ -93,5 +93,5 @@ int main(int argc, char **argv){ particle->userflags=0; mcpl_add_particle(outfile,particle); } - mcpl_closeandgzip_outfile(outfile); + mcpl_close_outfile(outfile); } From 80e69b2e0856d5af025087283964fd5a03e20339 Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 5 Jul 2022 01:53:14 +0200 Subject: [PATCH 08/62] now compiles fine if mcpl is indeed installed --- .../regression_tests/source_mcpl_file/test.py | 21 ++++++++++++------- 1 file changed, 14 insertions(+), 7 deletions(-) diff --git a/tests/regression_tests/source_mcpl_file/test.py b/tests/regression_tests/source_mcpl_file/test.py index 5a6a87d878..0adf67f1aa 100644 --- a/tests/regression_tests/source_mcpl_file/test.py +++ b/tests/regression_tests/source_mcpl_file/test.py @@ -3,16 +3,13 @@ import os import shutil import subprocess import textwrap - +import glob import openmc import pytest from tests.testing_harness import TestHarness -def model(): - subprocess.run(['gcc','-o gen_dummt_mcpl.out','-lm','-lmcpl','gen_dummy_mcpl.c'] ) - subprocess.run(['./gen_dummy_mcpl.out']) - + class SourceMCPLFileTestHarness(TestHarness): def execute_test(self): """Run OpenMC with the appropriate arguments and check the outputs.""" @@ -26,6 +23,10 @@ class SourceMCPLFileTestHarness(TestHarness): finally: self._cleanup() + def _create_input(self): + subprocess.run(['gcc','-o','gen_dummy_mcpl.out','gen_dummy_mcpl.c','-lm','-lmcpl']) + subprocess.run(['./gen_dummy_mcpl.out']) + def update_results(self): """Update the results_true using the current version of OpenMC.""" try: @@ -40,13 +41,19 @@ class SourceMCPLFileTestHarness(TestHarness): def _test_output_created(self): """Check that the output files were created""" - stat epoint = glob.glob(os.path.join(os.getcwd(), self._sp_name)) + statepoint = glob.glob(os.path.join(os.getcwd(), self._sp_name)) assert len(statepoint) == 1, 'Either multiple or no statepoint files ' \ 'exist.' assert statepoint[0].endswith('h5'), \ 'statepoint file does not end with h5.' + def _cleanup(self): + super()._cleanup() + source_mcpl=glob.glob(os.path.join(os.getcwd(),'source*.mcpl')) + for f in source_mcpl: + if (os.path.exists(f)): + os.remove(f) + def test_mcpl_source_file(): harness = SourceMCPLFileTestHarness('statepoint.10.h5') harness.main() - From 6d352aa55c8fab8bf427bceeb165c6b870d94e53 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 22 Jun 2022 22:50:07 +0200 Subject: [PATCH 09/62] add placeholder skeleton for MCPL-input --- include/openmc/source.h | 19 +++++++++++++++++++ src/source.cpp | 21 +++++++++++++++++++++ 2 files changed, 40 insertions(+) diff --git a/include/openmc/source.h b/include/openmc/source.h index 18fddc6899..cadbfe93bb 100644 --- a/include/openmc/source.h +++ b/include/openmc/source.h @@ -100,6 +100,25 @@ private: vector sites_; //!< Source sites from a file }; +//============================================================================== +// MCPL-file input source +//============================================================================== +class MCPLFileSource : public Source { +public: + // Constructors, destructors + MCPLFileSource(std::string path); + ~MCPLFileSource(); + + // Defer implementation to custom source library + SourceSite sample(uint64_t* seed) const override + +private: + vector sites_; //! Date: Thu, 23 Jun 2022 11:34:53 +0200 Subject: [PATCH 10/62] include mcpl header --- src/source.cpp | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/src/source.cpp b/src/source.cpp index cd8b5e5743..5e26fa0a8d 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -31,6 +31,10 @@ #include "openmc/state_point.h" #include "openmc/xml_interface.h" +#ifdef HAS_MCPL +#include +#endif + namespace openmc { //============================================================================== From c43989af9c222b20413ee3b892feaf20919b8707 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Thu, 23 Jun 2022 11:35:46 +0200 Subject: [PATCH 11/62] declarations of MCPL-internals --- include/openmc/source.h | 13 +++++++++++-- 1 file changed, 11 insertions(+), 2 deletions(-) diff --git a/include/openmc/source.h b/include/openmc/source.h index cadbfe93bb..69dd47eae9 100644 --- a/include/openmc/source.h +++ b/include/openmc/source.h @@ -109,11 +109,20 @@ public: MCPLFileSource(std::string path); ~MCPLFileSource(); - // Defer implementation to custom source library - SourceSite sample(uint64_t* seed) const override + // Properties + ParticleType particle_type() const { return particle_; } + + //! Sample from the external source distribution + //! \param[inout] seed Pseudorandom seed pointer + //! \return Site read from MCPL-file + SourceSite sample(uint64_t* seed) const override; private: + ParticleType particle_ {ParticleType::neutron}; //!< Type of particle emitted vector sites_; //! Date: Thu, 23 Jun 2022 11:37:05 +0200 Subject: [PATCH 12/62] internals of sample and constructor --- src/source.cpp | 51 ++++++++++++++++++++++++++++++++++++++++++++++---- 1 file changed, 47 insertions(+), 4 deletions(-) diff --git a/src/source.cpp b/src/source.cpp index 5e26fa0a8d..231f355fc2 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -337,7 +337,7 @@ CustomSourceWrapper::~CustomSourceWrapper() } //=========================================================================== -// Read an MCPL-file +// Read particles from an MCPL-file //=========================================================================== MCPLFileSource::MCPLFileSource(std::string path) { @@ -348,12 +348,55 @@ MCPLFileSource::MCPLFileSource(std::string path) // Read the source from a binary file instead of sampling from some // assumed source distribution - write_message(6, "Reading source file from {}...", path); + write_message(6, "Reading mcpl source file from {}",path); // Open the mcpl file - id_t file_id = mcpl_open(path.c_str) - + mcpl_file = mcpl_open(path.c_str); + //do checks on the mcpl_file to see if particles are many enough. + // should model this on the example source shown in the docs. + uint64_t nparticles=mcpl_hdr_nparticles(mcpl_file); + +} + +MCPLFileSource::~MCPLFilsSource(){ + mcpl_close_file(mcpl_file); +} + +SourceSite MCPLFileSource::sample(uint64_t *seed){ + SourceSite omc_particle; + mcpl_particle *mcpl_particle; + //extract particle from mcpl-file + mcpl_particle=mcpl_reazd(mcplfile); + // check if it is a neutron or a photon. otherwise skip + while ( mcpl_particle->pdgcode!=2112 && mcpl_particle->pdgcode!=22 ) + { + mcpl_particle(mcpl_read(mcplfile); + //check for file exhaustion + } + + if(mcpl_particle->pdgcode=2112) { + omc_particle_=ParticleType::neutron; + } else { + omc_particle_=Particletype::photon; + } + + //particle is good, convert to openmc-formalism + omc_particle.r.x=mcpl_particle.x; + omc_particle.r.y=mcpl_particle.y; + omc_particle.r.z=mcpl_particle.z; + + omc_particle.u.x=mcpl_particle->direction[0]; + omc_particle.u.y=mcpl_particle->direction[1]; + omc_particle.u.z=mcpl_particle->direction[2]; + + //mcpl stores particles in MeV + omc_particle.E=mcpl_particle->ekin*1e6; + + omc_particle.t=mcpl_particle->time*1e-3; + omc_particle.wgt=mcpl_particle->weight; + omc_particle.delayed_group=0; + return omc_particle; } From 103c6ffa0330ab1313489a27b2d9504cfb30f7ba Mon Sep 17 00:00:00 2001 From: erkn Date: Mon, 27 Jun 2022 23:25:35 +0200 Subject: [PATCH 13/62] (wip): skeleton for input only test --- .../regression_tests/source_mcpl_file/test.py | 52 +++++++++++++++++++ 1 file changed, 52 insertions(+) create mode 100644 tests/regression_tests/source_mcpl_file/test.py diff --git a/tests/regression_tests/source_mcpl_file/test.py b/tests/regression_tests/source_mcpl_file/test.py new file mode 100644 index 0000000000..5a6a87d878 --- /dev/null +++ b/tests/regression_tests/source_mcpl_file/test.py @@ -0,0 +1,52 @@ +import pathlib as pl +import os +import shutil +import subprocess +import textwrap + +import openmc +import pytest + +from tests.testing_harness import TestHarness + +def model(): + subprocess.run(['gcc','-o gen_dummt_mcpl.out','-lm','-lmcpl','gen_dummy_mcpl.c'] ) + subprocess.run(['./gen_dummy_mcpl.out']) + +class SourceMCPLFileTestHarness(TestHarness): + def execute_test(self): + """Run OpenMC with the appropriate arguments and check the outputs.""" + try: + self._create_input() + self._run_openmc() + self._test_output_created() + results = self._get_results() + self._write_results(results) + self._compare_results() + finally: + self._cleanup() + + def update_results(self): + """Update the results_true using the current version of OpenMC.""" + try: + self._create_input() + self._run_openmc() + self._test_output_created() + results = self._get_results() + self._write_results(results) + self._overwrite_results() + finally: + self._cleanup() + + def _test_output_created(self): + """Check that the output files were created""" + stat epoint = glob.glob(os.path.join(os.getcwd(), self._sp_name)) + assert len(statepoint) == 1, 'Either multiple or no statepoint files ' \ + 'exist.' + assert statepoint[0].endswith('h5'), \ + 'statepoint file does not end with h5.' + +def test_mcpl_source_file(): + harness = SourceMCPLFileTestHarness('statepoint.10.h5') + harness.main() + From 73fb335350bb3b6ed0199e9c3cc4ee19d73b163e Mon Sep 17 00:00:00 2001 From: erkn Date: Mon, 27 Jun 2022 23:26:01 +0200 Subject: [PATCH 14/62] c-prog to generate a dummy mcpl --- .../source_mcpl_file/gen_dummy_mcpl.c | 97 +++++++++++++++++++ 1 file changed, 97 insertions(+) create mode 100644 tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c diff --git a/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c b/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c new file mode 100644 index 0000000000..7b5a933f18 --- /dev/null +++ b/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c @@ -0,0 +1,97 @@ +#include +#include +#include +#include + +const int batches=10; +const int particles=10000; + +double rand01(){ + double r=rand()/((double) RAND_MAX); + return r; +} + + +int main(int argc, char **argv){ + char outfilename[128]; + snprintf(outfilename,127,"source.%d.mcpl",batches); + + /*generate an mcpl_header_of sorts*/ + mcpl_outfile_t outfile=mcpl_create_outfile(outfilename); + mcpl_particle_t *particle, Particle; + + char line[256]; + snprintf(line,255,"gen_dummy_mcpl.c"); + mcpl_hdr_set_srcname(outfile,line); + mcpl_enable_universal_pdgcode(outfile,2112);/*all particles are neutrons*/ + snprintf(line,255,"Dummy MCPL-file for testing OpenMC mcpl input"); + mcpl_hdr_add_comment(outfile,line); + + /*also add the instrument file and the command line as blobs*/ + FILE *fp; + if( (fp=fopen("gen_dummy_mcpl.c","rb"))!=NULL){ + unsigned char *buffer; + int size,status; + /*find the file size by seeking to end, "tell" the position, and then go back again*/ + fseek(fp, 0L, SEEK_END); + size = ftell(fp); // get current file pointer + fseek(fp, 0L, SEEK_SET); // seek back to beginning of file + if ( size && (buffer=malloc(size))!=NULL){ + if (size!=(fread(buffer,1,size,fp))){ + fprintf(stderr,"Warning: c source generator file not read cleanly\n"); + } + mcpl_hdr_add_data(outfile, "mcpl_point_source_file_generator", size, buffer); + free(buffer); + } + fclose(fp); + } else { + fprintf(stderr,"Warning: could not open c source generator file, hence not embedded.\n"); + } + + + /*the main particle loop*/ + particle=&Particle; + int i; + for (i=0;iposition[0]=0; + particle->position[1]=0; + particle->position[2]=0; + + /*generate a random direction on unit sphere*/ + double nrm=2.0; + double vx,vy,vz; + int iter=0; + do { + if(iter>100){ + printf("warning: exceeed max random vector attempts: at loop %d -> %g %d\n",i,nrm,iter); + } + vx=rand01()*2.0-1.0; + vy=rand01()*2.0-1.0; + vz=rand01()*2.0-1.0; + nrm=(vx*vx + vy*vy + vz*vz); + iter++; + } while (nrm>1.0); + vx/=sqrt(nrm); + vy/=sqrt(nrm); + vz/=sqrt(nrm); + particle->direction[0]=vx; + particle->direction[1]=vy; + particle->direction[2]=vz; + + /*generate the kinetic energy (in MeV) of the particle*/ + particle->ekin=1e-3; + + /*set time=0*/ + particle->time=0; + /*set weight*/ + particle->weight=1; + + particle->userflags=0; + mcpl_add_particle(outfile,particle); + } + mcpl_closeandgzip_outfile(outfile); +} From c6c876d835fe2177426ee64f4ae8dc15d5e616f2 Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 5 Jul 2022 01:30:07 +0200 Subject: [PATCH 15/62] don't zip the output-file --- .../source_mcpl_file/gen_dummy_mcpl.c | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c b/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c index 7b5a933f18..f35a07fc44 100644 --- a/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c +++ b/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c @@ -54,20 +54,20 @@ int main(int argc, char **argv){ int i; for (i=0;iposition[0]=0; particle->position[1]=0; particle->position[2]=0; - + /*generate a random direction on unit sphere*/ double nrm=2.0; double vx,vy,vz; int iter=0; do { if(iter>100){ - printf("warning: exceeed max random vector attempts: at loop %d -> %g %d\n",i,nrm,iter); + printf("Warning: Exceeded max random vector attempts: at loop %d -> %g %d\n",i,nrm,iter); } vx=rand01()*2.0-1.0; vy=rand01()*2.0-1.0; @@ -82,9 +82,9 @@ int main(int argc, char **argv){ particle->direction[1]=vy; particle->direction[2]=vz; - /*generate the kinetic energy (in MeV) of the particle*/ + /*generate the kinetic energy (in MeV) of the particle*/ particle->ekin=1e-3; - + /*set time=0*/ particle->time=0; /*set weight*/ @@ -93,5 +93,5 @@ int main(int argc, char **argv){ particle->userflags=0; mcpl_add_particle(outfile,particle); } - mcpl_closeandgzip_outfile(outfile); + mcpl_close_outfile(outfile); } From 8719a539f9f8c9134ded079b2eec64c6cdf3aa82 Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 5 Jul 2022 01:53:14 +0200 Subject: [PATCH 16/62] now compiles fine if mcpl is indeed installed --- .../regression_tests/source_mcpl_file/test.py | 21 ++++++++++++------- 1 file changed, 14 insertions(+), 7 deletions(-) diff --git a/tests/regression_tests/source_mcpl_file/test.py b/tests/regression_tests/source_mcpl_file/test.py index 5a6a87d878..0adf67f1aa 100644 --- a/tests/regression_tests/source_mcpl_file/test.py +++ b/tests/regression_tests/source_mcpl_file/test.py @@ -3,16 +3,13 @@ import os import shutil import subprocess import textwrap - +import glob import openmc import pytest from tests.testing_harness import TestHarness -def model(): - subprocess.run(['gcc','-o gen_dummt_mcpl.out','-lm','-lmcpl','gen_dummy_mcpl.c'] ) - subprocess.run(['./gen_dummy_mcpl.out']) - + class SourceMCPLFileTestHarness(TestHarness): def execute_test(self): """Run OpenMC with the appropriate arguments and check the outputs.""" @@ -26,6 +23,10 @@ class SourceMCPLFileTestHarness(TestHarness): finally: self._cleanup() + def _create_input(self): + subprocess.run(['gcc','-o','gen_dummy_mcpl.out','gen_dummy_mcpl.c','-lm','-lmcpl']) + subprocess.run(['./gen_dummy_mcpl.out']) + def update_results(self): """Update the results_true using the current version of OpenMC.""" try: @@ -40,13 +41,19 @@ class SourceMCPLFileTestHarness(TestHarness): def _test_output_created(self): """Check that the output files were created""" - stat epoint = glob.glob(os.path.join(os.getcwd(), self._sp_name)) + statepoint = glob.glob(os.path.join(os.getcwd(), self._sp_name)) assert len(statepoint) == 1, 'Either multiple or no statepoint files ' \ 'exist.' assert statepoint[0].endswith('h5'), \ 'statepoint file does not end with h5.' + def _cleanup(self): + super()._cleanup() + source_mcpl=glob.glob(os.path.join(os.getcwd(),'source*.mcpl')) + for f in source_mcpl: + if (os.path.exists(f)): + os.remove(f) + def test_mcpl_source_file(): harness = SourceMCPLFileTestHarness('statepoint.10.h5') harness.main() - From 0d1c346b5a12d2a452450b435a80cf4a2fbbe608 Mon Sep 17 00:00:00 2001 From: erkn Date: Sat, 9 Jul 2022 16:31:00 +0200 Subject: [PATCH 17/62] fix misprints --- src/source.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/source.cpp b/src/source.cpp index 231f355fc2..37e580b18d 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -345,7 +345,7 @@ MCPLFileSource::MCPLFileSource(std::string path) if (!file_exists(path)) { fatal_error(fmt::format("Source file '{}' does not exist.", path)); } - + // Read the source from a binary file instead of sampling from some // assumed source distribution write_message(6, "Reading mcpl source file from {}",path); @@ -367,7 +367,7 @@ SourceSite MCPLFileSource::sample(uint64_t *seed){ SourceSite omc_particle; mcpl_particle *mcpl_particle; //extract particle from mcpl-file - mcpl_particle=mcpl_reazd(mcplfile); + mcpl_particle=mcpl_read(mcplfile); // check if it is a neutron or a photon. otherwise skip while ( mcpl_particle->pdgcode!=2112 && mcpl_particle->pdgcode!=22 ) { @@ -375,9 +375,9 @@ SourceSite MCPLFileSource::sample(uint64_t *seed){ //check for file exhaustion } - if(mcpl_particle->pdgcode=2112) { + if(mcpl_particle->pdgcode==2112) { omc_particle_=ParticleType::neutron; - } else { + } else if (mcpl_particle->pdgcode==22){ omc_particle_=Particletype::photon; } From 7a6c69fbd5daefcca4b612affdbac6ecd46e9d32 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Mon, 11 Jul 2022 09:49:29 +0200 Subject: [PATCH 18/62] fix misprints --- src/source.cpp | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/source.cpp b/src/source.cpp index 231f355fc2..9ec1ff798c 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -367,7 +367,7 @@ SourceSite MCPLFileSource::sample(uint64_t *seed){ SourceSite omc_particle; mcpl_particle *mcpl_particle; //extract particle from mcpl-file - mcpl_particle=mcpl_reazd(mcplfile); + mcpl_particle=mcpl_read(mcplfile); // check if it is a neutron or a photon. otherwise skip while ( mcpl_particle->pdgcode!=2112 && mcpl_particle->pdgcode!=22 ) { @@ -375,9 +375,9 @@ SourceSite MCPLFileSource::sample(uint64_t *seed){ //check for file exhaustion } - if(mcpl_particle->pdgcode=2112) { + if(mcpl_particle->pdgcode==2112) { omc_particle_=ParticleType::neutron; - } else { + } else if (mcpl_particle->pdgcode==22){ omc_particle_=Particletype::photon; } @@ -390,9 +390,9 @@ SourceSite MCPLFileSource::sample(uint64_t *seed){ omc_particle.u.y=mcpl_particle->direction[1]; omc_particle.u.z=mcpl_particle->direction[2]; - //mcpl stores particles in MeV + //mcpl stores kinetic energy in MeV omc_particle.E=mcpl_particle->ekin*1e6; - + //mcpl stores time in ms omc_particle.t=mcpl_particle->time*1e-3; omc_particle.wgt=mcpl_particle->weight; omc_particle.delayed_group=0; From 34746926d349d0e04134893aa987c49d9508e121 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 12 Jul 2022 13:48:55 +0200 Subject: [PATCH 19/62] add necessary xmls --- tests/regression_tests/source_mcpl_file/geometry.xml | 8 ++++++++ tests/regression_tests/source_mcpl_file/materials.xml | 9 +++++++++ 2 files changed, 17 insertions(+) create mode 100644 tests/regression_tests/source_mcpl_file/geometry.xml create mode 100644 tests/regression_tests/source_mcpl_file/materials.xml diff --git a/tests/regression_tests/source_mcpl_file/geometry.xml b/tests/regression_tests/source_mcpl_file/geometry.xml new file mode 100644 index 0000000000..bc56030e18 --- /dev/null +++ b/tests/regression_tests/source_mcpl_file/geometry.xml @@ -0,0 +1,8 @@ + + + + + + + + diff --git a/tests/regression_tests/source_mcpl_file/materials.xml b/tests/regression_tests/source_mcpl_file/materials.xml new file mode 100644 index 0000000000..2472a74717 --- /dev/null +++ b/tests/regression_tests/source_mcpl_file/materials.xml @@ -0,0 +1,9 @@ + + + + + + + + + From 25cfc328d04c930a8d75ba41eb46ed7e26c64cae Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 12 Jul 2022 13:49:56 +0200 Subject: [PATCH 20/62] set seed to make test deteministic --- tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c | 1 + 1 file changed, 1 insertion(+) diff --git a/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c b/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c index f35a07fc44..61eea95c6f 100644 --- a/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c +++ b/tests/regression_tests/source_mcpl_file/gen_dummy_mcpl.c @@ -52,6 +52,7 @@ int main(int argc, char **argv){ /*the main particle loop*/ particle=&Particle; int i; + srand(1234); for (i=0;i Date: Tue, 12 Jul 2022 13:50:35 +0200 Subject: [PATCH 21/62] test target result --- tests/regression_tests/source_mcpl_file/results_true.dat | 2 ++ 1 file changed, 2 insertions(+) create mode 100644 tests/regression_tests/source_mcpl_file/results_true.dat diff --git a/tests/regression_tests/source_mcpl_file/results_true.dat b/tests/regression_tests/source_mcpl_file/results_true.dat new file mode 100644 index 0000000000..f2209006a1 --- /dev/null +++ b/tests/regression_tests/source_mcpl_file/results_true.dat @@ -0,0 +1,2 @@ +k-combined: +3.017557E-01 3.398770E-03 From b0bb63cb3508d0dff59fd950a60dc35db2188b19 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 12 Jul 2022 13:50:53 +0200 Subject: [PATCH 22/62] yet another necessary file --- tests/regression_tests/source_mcpl_file/settings.xml | 11 +++++++++++ 1 file changed, 11 insertions(+) create mode 100644 tests/regression_tests/source_mcpl_file/settings.xml diff --git a/tests/regression_tests/source_mcpl_file/settings.xml b/tests/regression_tests/source_mcpl_file/settings.xml new file mode 100644 index 0000000000..8b1510e0ce --- /dev/null +++ b/tests/regression_tests/source_mcpl_file/settings.xml @@ -0,0 +1,11 @@ + + + + 10 + 5 + 1000 + + + source.10.mcpl + + From e480b5c851c542e53ff6e83236994f23c2bfdba3 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 12 Jul 2022 14:40:01 +0200 Subject: [PATCH 23/62] add a check that the test input got created as intended --- tests/regression_tests/source_mcpl_file/test.py | 12 +++++++++++- 1 file changed, 11 insertions(+), 1 deletion(-) diff --git a/tests/regression_tests/source_mcpl_file/test.py b/tests/regression_tests/source_mcpl_file/test.py index 0adf67f1aa..922ddddd31 100644 --- a/tests/regression_tests/source_mcpl_file/test.py +++ b/tests/regression_tests/source_mcpl_file/test.py @@ -15,6 +15,7 @@ class SourceMCPLFileTestHarness(TestHarness): """Run OpenMC with the appropriate arguments and check the outputs.""" try: self._create_input() + self._test_input_created() self._run_openmc() self._test_output_created() results = self._get_results() @@ -31,6 +32,7 @@ class SourceMCPLFileTestHarness(TestHarness): """Update the results_true using the current version of OpenMC.""" try: self._create_input() + self._test_input_created() self._run_openmc() self._test_output_created() results = self._get_results() @@ -39,6 +41,14 @@ class SourceMCPLFileTestHarness(TestHarness): finally: self._cleanup() + def _test_input_created(self): + """Check that the input mcpl.file was generated as it should""" + mcplfile=glob.glob(os.path.join(os.getcwd(),'source.10.mcpl')) + assert len(mcplfile) == 1, 'Either multiple or no mcpl files ' \ + 'exist.' + assert mcplfile[0].endswith('mcpl'), \ + 'output file does not end with mcpl.' + def _test_output_created(self): """Check that the output files were created""" statepoint = glob.glob(os.path.join(os.getcwd(), self._sp_name)) @@ -55,5 +65,5 @@ class SourceMCPLFileTestHarness(TestHarness): os.remove(f) def test_mcpl_source_file(): - harness = SourceMCPLFileTestHarness('statepoint.10.h5') + harness = SourceMCPLFileTestHarness('source.10.mcpl') harness.main() From 8c9bb761e5d3838a33822f7fa283acab106cc736 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 13 Jul 2022 10:25:55 +0200 Subject: [PATCH 24/62] fix bugs and finalize interface to MCPL --- include/openmc/source.h | 22 +++++++------ src/source.cpp | 73 +++++++++++++++++++++++++---------------- 2 files changed, 57 insertions(+), 38 deletions(-) diff --git a/include/openmc/source.h b/include/openmc/source.h index 69dd47eae9..4d7ec104aa 100644 --- a/include/openmc/source.h +++ b/include/openmc/source.h @@ -12,6 +12,10 @@ #include "openmc/particle.h" #include "openmc/vector.h" +#ifdef OPENMC_MCPL +#include +#endif + namespace openmc { //============================================================================== @@ -100,6 +104,7 @@ private: vector sites_; //!< Source sites from a file }; +#ifdef OPENMC_MCPL //============================================================================== // MCPL-file input source //============================================================================== @@ -109,24 +114,21 @@ public: MCPLFileSource(std::string path); ~MCPLFileSource(); - // Properties - ParticleType particle_type() const { return particle_; } - + // Methods //! Sample from the external source distribution //! \param[inout] seed Pseudorandom seed pointer //! \return Site read from MCPL-file - SourceSite sample(uint64_t* seed) const override; + SourceSite sample(uint64_t* seed) const; private: - ParticleType particle_ {ParticleType::neutron}; //!< Type of particle emitted + SourceSite read_single_particle() const; + void read_source_bank(vector &sites_); + vector sites_; //! #endif @@ -336,6 +336,7 @@ CustomSourceWrapper::~CustomSourceWrapper() #endif } +#ifdef OPENMC_MCPL //=========================================================================== // Read particles from an MCPL-file //=========================================================================== @@ -351,54 +352,70 @@ MCPLFileSource::MCPLFileSource(std::string path) write_message(6, "Reading mcpl source file from {}",path); // Open the mcpl file - mcpl_file = mcpl_open(path.c_str); + mcpl_file = mcpl_open_file(path.c_str()); //do checks on the mcpl_file to see if particles are many enough. // should model this on the example source shown in the docs. - uint64_t nparticles=mcpl_hdr_nparticles(mcpl_file); + n_sites=mcpl_hdr_nparticles(mcpl_file); + read_source_bank(sites_); } -MCPLFileSource::~MCPLFilsSource(){ +MCPLFileSource::~MCPLFileSource(){ mcpl_close_file(mcpl_file); } -SourceSite MCPLFileSource::sample(uint64_t *seed){ - SourceSite omc_particle; - mcpl_particle *mcpl_particle; +SourceSite MCPLFileSource::sample(uint64_t* seed) const +{ + size_t i_site = sites_.size() * prn(seed); + return sites_[i_site]; +} + +void MCPLFileSource::read_source_bank(vector &sites_) +{ + sites_.resize(n_sites); + for (int i=0;ipdgcode!=2112 && mcpl_particle->pdgcode!=22 ) - { - mcpl_particle(mcpl_read(mcplfile); - //check for file exhaustion + mcpl_particle=mcpl_read(mcpl_file); + // check if it is a neutron or a photon. otherwise skip + while ( mcpl_particle->pdgcode!=2112 && mcpl_particle->pdgcode!=22 ) { + mcpl_particle=mcpl_read(mcpl_file); + //should check for file exhaustion This could happen if particles are other than + //neutrons or photons } if(mcpl_particle->pdgcode==2112) { - omc_particle_=ParticleType::neutron; - } else if (mcpl_particle->pdgcode==22){ - omc_particle_=Particletype::photon; + omc_particle_.particle=ParticleType::neutron; + } else if (mcpl_particle->pdgcode==22) { + omc_particle_.particle=ParticleType::photon; } //particle is good, convert to openmc-formalism - omc_particle.r.x=mcpl_particle.x; - omc_particle.r.y=mcpl_particle.y; - omc_particle.r.z=mcpl_particle.z; + omc_particle_.r.x=mcpl_particle->position[0]; + omc_particle_.r.y=mcpl_particle->position[1]; + omc_particle_.r.z=mcpl_particle->position[2]; - omc_particle.u.x=mcpl_particle->direction[0]; - omc_particle.u.y=mcpl_particle->direction[1]; - omc_particle.u.z=mcpl_particle->direction[2]; + omc_particle_.u.x=mcpl_particle->direction[0]; + omc_particle_.u.y=mcpl_particle->direction[1]; + omc_particle_.u.z=mcpl_particle->direction[2]; //mcpl stores kinetic energy in MeV - omc_particle.E=mcpl_particle->ekin*1e6; + omc_particle_.E=mcpl_particle->ekin*1e6; //mcpl stores time in ms - omc_particle.t=mcpl_particle->time*1e-3; - omc_particle.wgt=mcpl_particle->weight; - omc_particle.delayed_group=0; - return omc_particle; + omc_particle_.time=mcpl_particle->time*1e-3; + omc_particle_.wgt=mcpl_particle->weight; + omc_particle_.delayed_group=0; + return omc_particle_; } - +#endif //OPENMC_MCPL //============================================================================== // Non-member functions From 915d69fb3d217cf705158dba680b8f02a99208c7 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 13 Jul 2022 10:26:24 +0200 Subject: [PATCH 25/62] add an option to turn on mcpl-support --- CMakeLists.txt | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/CMakeLists.txt b/CMakeLists.txt index 85dcd7084d..6af9643cdb 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -36,6 +36,7 @@ option(OPENMC_ENABLE_COVERAGE "Compile with coverage analysis flags" option(OPENMC_USE_DAGMC "Enable support for DAGMC (CAD) geometry" OFF) option(OPENMC_USE_LIBMESH "Enable support for libMesh unstructured mesh tallies" OFF) option(OPENMC_USE_MPI "Enable MPI" OFF) +option(OPENMC_USE_MCPL "Enable MPCPL" OFF) #=============================================================================== # Set a default build configuration if not explicitly specified @@ -154,6 +155,10 @@ if(OPENMC_ENABLE_COVERAGE) list(APPEND ldflags --coverage) endif() +if(OPENMC_USE_MCPL) + list(APPEND ldflags -lmcpl) +endif() + # Show flags being used message(STATUS "OpenMC C++ flags: ${cxxflags}") message(STATUS "OpenMC Linker flags: ${ldflags}") @@ -408,6 +413,10 @@ endif() if (OPENMC_USE_MPI) target_compile_definitions(libopenmc PUBLIC -DOPENMC_MPI) endif() +if(OPENMC_USE_MCPL) + target_compile_definitions(libopenmc PUBLIC -DOPENMC_MCPL) +endif() + # Set git SHA1 hash as a compile definition if(GIT_FOUND) From b140eeb7cb4c7a1ee76172cd21c32c2e9951a24b Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 13 Jul 2022 10:27:18 +0200 Subject: [PATCH 26/62] add a source type to the settings interface --- src/settings.cpp | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/settings.cpp b/src/settings.cpp index eb0d4936c7..12e0a80509 100644 --- a/src/settings.cpp +++ b/src/settings.cpp @@ -430,6 +430,9 @@ void read_settings_xml() if (check_for_node(node, "file")) { auto path = get_node_value(node, "file", false, true); model::external_sources.push_back(make_unique(path)); + } else if (check_for_node(node, "mcpl")) { + auto path = get_node_value(node, "mcpl", false, true); + model::external_sources.push_back(make_unique(path)); } else if (check_for_node(node, "library")) { // Get shared library path and parameters auto path = get_node_value(node, "library", false, true); From cd163ef0b91fb4acc30630528b0728d8bab3098e Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 13 Jul 2022 12:13:10 +0200 Subject: [PATCH 27/62] only add mcpl-support if OPENMC_MCPL is defined --- src/settings.cpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/settings.cpp b/src/settings.cpp index 12e0a80509..d419cab06f 100644 --- a/src/settings.cpp +++ b/src/settings.cpp @@ -430,9 +430,11 @@ void read_settings_xml() if (check_for_node(node, "file")) { auto path = get_node_value(node, "file", false, true); model::external_sources.push_back(make_unique(path)); +#ifdef OPENMC_MCPL } else if (check_for_node(node, "mcpl")) { auto path = get_node_value(node, "mcpl", false, true); model::external_sources.push_back(make_unique(path)); +#endif } else if (check_for_node(node, "library")) { // Get shared library path and parameters auto path = get_node_value(node, "library", false, true); From a79233b3ec85f59647fb800da7db0973ce9365e8 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 13 Jul 2022 14:13:33 +0200 Subject: [PATCH 28/62] needed to avoid name clashes with other tests --- tests/regression_tests/source_mcpl_file/__init__.py | 0 1 file changed, 0 insertions(+), 0 deletions(-) create mode 100644 tests/regression_tests/source_mcpl_file/__init__.py diff --git a/tests/regression_tests/source_mcpl_file/__init__.py b/tests/regression_tests/source_mcpl_file/__init__.py new file mode 100644 index 0000000000..e69de29bb2 From eff690fd90032abd0d4a5e68f5babcc38c476de1 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 13 Jul 2022 20:55:53 +0200 Subject: [PATCH 29/62] fail without tripping the entire test-process --- tests/regression_tests/source_mcpl_file/test.py | 8 ++------ 1 file changed, 2 insertions(+), 6 deletions(-) diff --git a/tests/regression_tests/source_mcpl_file/test.py b/tests/regression_tests/source_mcpl_file/test.py index 922ddddd31..f55a769047 100644 --- a/tests/regression_tests/source_mcpl_file/test.py +++ b/tests/regression_tests/source_mcpl_file/test.py @@ -1,11 +1,6 @@ -import pathlib as pl import os -import shutil import subprocess -import textwrap import glob -import openmc -import pytest from tests.testing_harness import TestHarness @@ -25,7 +20,8 @@ class SourceMCPLFileTestHarness(TestHarness): self._cleanup() def _create_input(self): - subprocess.run(['gcc','-o','gen_dummy_mcpl.out','gen_dummy_mcpl.c','-lm','-lmcpl']) + compiled=subprocess.run(['gcc','-o','gen_dummy_mcpl.out','gen_dummy_mcpl.c','-lm','-lmcpl']) + assert compiled==0, 'Could not compile mcpl-file generator code' subprocess.run(['./gen_dummy_mcpl.out']) def update_results(self): From 1ba696741412b17ac676cc2a0b4a78c606906d35 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Fri, 15 Jul 2022 21:35:59 +0200 Subject: [PATCH 30/62] fix typo Co-authored-by: Paul Romano --- CMakeLists.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 6af9643cdb..be56f1c12f 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -36,7 +36,7 @@ option(OPENMC_ENABLE_COVERAGE "Compile with coverage analysis flags" option(OPENMC_USE_DAGMC "Enable support for DAGMC (CAD) geometry" OFF) option(OPENMC_USE_LIBMESH "Enable support for libMesh unstructured mesh tallies" OFF) option(OPENMC_USE_MPI "Enable MPI" OFF) -option(OPENMC_USE_MCPL "Enable MPCPL" OFF) +option(OPENMC_USE_MCPL "Enable MCPL" OFF) #=============================================================================== # Set a default build configuration if not explicitly specified From 380cfc150b449ef56debb5261caa010a52c42bdd Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 9 Aug 2022 08:12:19 +0200 Subject: [PATCH 31/62] refactor MCPL file read part of FileSource A new constructor is added to FileSource that reads from a given MCPL-file. --- include/openmc/source.h | 30 +--------- src/settings.cpp | 2 +- src/source.cpp | 128 +++++++++++++++------------------------- 3 files changed, 52 insertions(+), 108 deletions(-) diff --git a/include/openmc/source.h b/include/openmc/source.h index 4d7ec104aa..5c2dae9311 100644 --- a/include/openmc/source.h +++ b/include/openmc/source.h @@ -96,7 +96,9 @@ class FileSource : public Source { public: // Constructors explicit FileSource(std::string path); - +#ifdef OPENMC_MCPL + explicit FileSource(mcpl_file_t mcpl_file); +#endif // Methods SourceSite sample(uint64_t* seed) const override; @@ -104,32 +106,6 @@ private: vector sites_; //!< Source sites from a file }; -#ifdef OPENMC_MCPL -//============================================================================== -// MCPL-file input source -//============================================================================== -class MCPLFileSource : public Source { -public: - // Constructors, destructors - MCPLFileSource(std::string path); - ~MCPLFileSource(); - - // Methods - //! Sample from the external source distribution - //! \param[inout] seed Pseudorandom seed pointer - //! \return Site read from MCPL-file - SourceSite sample(uint64_t* seed) const; - -private: - SourceSite read_single_particle() const; - void read_source_bank(vector &sites_); - - vector sites_; //!(path)); + model::external_sources.push_back(make_unique(mcpl_open_file(path.c_str()))); #endif } else if (check_for_node(node, "library")) { // Get shared library path and parameters diff --git a/src/source.cpp b/src/source.cpp index 82978f1030..49e51c4b16 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -277,6 +277,54 @@ FileSource::FileSource(std::string path) file_close(file_id); } +#ifdef OPENMC_MCPL +FileSource::FileSource(mcpl_file_t mcpl_file) +{ + //do checks on the mcpl_file to see if particles are many enough. + // should model this on the example source shown in the docs. + size_t n_sites=mcpl_hdr_nparticles(mcpl_file); + + sites_.resize(n_sites); + for (int i=0;ipdgcode!=2112 && mcpl_particle->pdgcode!=22 ) { + mcpl_particle=mcpl_read(mcpl_file); + //should check for file exhaustion This could happen if particles are other than + //neutrons or photons + } + + if(mcpl_particle->pdgcode==2112) { + site_.particle=ParticleType::neutron; + } else if (mcpl_particle->pdgcode==22) { + site_.particle=ParticleType::photon; + } + + //particle is good, convert to openmc-formalism + site_.r.x=mcpl_particle->position[0]; + site_.r.y=mcpl_particle->position[1]; + site_.r.z=mcpl_particle->position[2]; + + site_.u.x=mcpl_particle->direction[0]; + site_.u.y=mcpl_particle->direction[1]; + site_.u.z=mcpl_particle->direction[2]; + + //mcpl stores kinetic energy in MeV + site_.E=mcpl_particle->ekin*1e6; + //mcpl stores time in ms + site_.time=mcpl_particle->time*1e-3; + site_.wgt=mcpl_particle->weight; + site_.delayed_group=0; + sites_[i]=site_; + } + mcpl_close_file(mcpl_file); +} +#endif //OPENMC_MCPL + SourceSite FileSource::sample(uint64_t* seed) const { size_t i_site = sites_.size() * prn(seed); @@ -336,86 +384,6 @@ CustomSourceWrapper::~CustomSourceWrapper() #endif } -#ifdef OPENMC_MCPL -//=========================================================================== -// Read particles from an MCPL-file -//=========================================================================== -MCPLFileSource::MCPLFileSource(std::string path) -{ - // Check if source file exists - if (!file_exists(path)) { - fatal_error(fmt::format("Source file '{}' does not exist.", path)); - } - - // Read the source from a binary file instead of sampling from some - // assumed source distribution - write_message(6, "Reading mcpl source file from {}",path); - - // Open the mcpl file - mcpl_file = mcpl_open_file(path.c_str()); - - //do checks on the mcpl_file to see if particles are many enough. - // should model this on the example source shown in the docs. - n_sites=mcpl_hdr_nparticles(mcpl_file); - - read_source_bank(sites_); -} - -MCPLFileSource::~MCPLFileSource(){ - mcpl_close_file(mcpl_file); -} - -SourceSite MCPLFileSource::sample(uint64_t* seed) const -{ - size_t i_site = sites_.size() * prn(seed); - return sites_[i_site]; -} - -void MCPLFileSource::read_source_bank(vector &sites_) -{ - sites_.resize(n_sites); - for (int i=0;ipdgcode!=2112 && mcpl_particle->pdgcode!=22 ) { - mcpl_particle=mcpl_read(mcpl_file); - //should check for file exhaustion This could happen if particles are other than - //neutrons or photons - } - - if(mcpl_particle->pdgcode==2112) { - omc_particle_.particle=ParticleType::neutron; - } else if (mcpl_particle->pdgcode==22) { - omc_particle_.particle=ParticleType::photon; - } - - //particle is good, convert to openmc-formalism - omc_particle_.r.x=mcpl_particle->position[0]; - omc_particle_.r.y=mcpl_particle->position[1]; - omc_particle_.r.z=mcpl_particle->position[2]; - - omc_particle_.u.x=mcpl_particle->direction[0]; - omc_particle_.u.y=mcpl_particle->direction[1]; - omc_particle_.u.z=mcpl_particle->direction[2]; - - //mcpl stores kinetic energy in MeV - omc_particle_.E=mcpl_particle->ekin*1e6; - //mcpl stores time in ms - omc_particle_.time=mcpl_particle->time*1e-3; - omc_particle_.wgt=mcpl_particle->weight; - omc_particle_.delayed_group=0; - return omc_particle_; -} -#endif //OPENMC_MCPL //============================================================================== // Non-member functions From aa571c2c33b421e2524e249ca00f9584e891dd0a Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 6 Sep 2022 02:22:21 +0200 Subject: [PATCH 32/62] wip add MCPL-out to openmc code for writing a source_bank/point compiles but is untested --- include/openmc/state_point.h | 9 +++ src/state_point.cpp | 125 +++++++++++++++++++++++++++++++++++ 2 files changed, 134 insertions(+) diff --git a/include/openmc/state_point.h b/include/openmc/state_point.h index 5ada2bf88d..2d92c935a6 100644 --- a/include/openmc/state_point.h +++ b/include/openmc/state_point.h @@ -9,6 +9,10 @@ #include "openmc/particle.h" #include "openmc/vector.h" +#ifdef OPENMC_MCPL +#include +#endif + namespace openmc { void load_state_point(); @@ -21,5 +25,10 @@ void write_tally_results_nr(hid_t file_id); void restart_set_keff(); void write_unstructured_mesh_results(); +#ifdef OPENMC_MCPL +void write_mcpl_source_point(const char *filename_, bool surf_source_bank = false); +void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank); +#endif + } // namespace openmc #endif // OPENMC_STATE_POINT_H diff --git a/src/state_point.cpp b/src/state_point.cpp index 9ea0459153..babd670a2d 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -28,6 +28,10 @@ #include "openmc/timer.h" #include "openmc/vector.h" +#ifdef OPENMC_MCPL +#include +#endif + namespace openmc { extern "C" int openmc_statepoint_write(const char* filename, bool* write_source) @@ -598,6 +602,45 @@ void write_source_point(const char* filename, bool surf_source_bank) file_close(file_id); } +void write_mcpl_point(const char *filename, bool surf_source_bank) +{ + std::string filename_; + if (filename) { + filename_ = filename; + } else { + // Determine width for zero padding + int w = std::to_string(settings::n_max_batches).size(); + + filename_ = fmt::format("{0}source.{1:0{2}}.mcpl", settings::path_output, + simulation::current_batch, w); + } + + mcpl_outfile_t file_id; + + std::string line; + if (mpi::master) { + // this must be rewritten: file_open is h5-specific + file_id = mcpl_create_outfile(filename_.c_str()); + //write_attribute(file_id, "filetype", "source"); + //write header stuff (oopy in xml-files as binary blobs for instance)) + if (VERSION_DEV){ + line=fmt::format("OpenMC {VERSION_MAJOR}.{VERSION_MINOR}.{VERSION_RELEASE}-development"); + } else { + line=fmt::format("OpenMC {VERSION_MAJOR}.{VERSION_MINOR}.{VERSION_RELEASE}"); + } + mcpl_hdr_set_srcname(file_id,line.c_str()); + } + + write_mcpl_source_bank(file_id, surf_source_bank); + + if (mpi::master) { + //change this - this is h5 specific + mcpl_closeandgzip_outfile(file_id); + } + +} + + void write_source_bank(hid_t group_id, bool surf_source_bank) { hid_t banktype = h5banktype(); @@ -714,6 +757,88 @@ void write_source_bank(hid_t group_id, bool surf_source_bank) H5Tclose(banktype); } +void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) +{ + int64_t dims_size = settings::n_particles; + int64_t count_size = simulation::work_per_rank; + + // Set vectors for source bank and starting bank index of each process + vector* bank_index = &simulation::work_index; + vector* source_bank = &simulation::source_bank; + vector surf_source_index_vector; + vector surf_source_bank_vector; + + if(surf_source_bank) { + surf_source_index_vector = calculate_surf_source_size(); + dims_size = surf_source_index_vector[mpi::n_procs]; + count_size = simulation::surf_source_bank.size(); + + bank_index = &surf_source_index_vector; + + // Copy data in a SharedArray into a vector. + surf_source_bank_vector.resize(count_size); + surf_source_bank_vector.assign(simulation::surf_source_bank.data(), + simulation::surf_source_bank.data() + count_size); + source_bank = &surf_source_bank_vector; + } + + if (mpi::master) { + //write particles from the master node + //loop over the other nodes and receive data - then write those. + for (int i = 0; i < mpi::n_procs; ++i) { +#ifdef OPENMC_MPI + if (i>0) + MPI_Recv(source_bank->data(), count[0], mpi::source_site, i, i, + mpi::intracomm, MPI_STATUS_IGNORE); +#endif + //now write the source_banke data again. + for (vector::iterator _site=source_bank->begin(); _site!=source_bank->end();_site++){ + //particle is now at the iterator + //write it to the mcpl-file + mcpl_particle_t p; + p.position[0]=_site->r.x; + p.position[1]=_site->r.y; + p.position[2]=_site->r.z; + + p.direction[0]=_site->u.x; + p.direction[1]=_site->u.y; + p.direction[2]=_site->u.x; + + //mcpl stores kinetic energy in MeV + p.ekin=_site->E*1e-6; + + p.time=_site->time*1e3; + + p.weight=_site->wgt; + + switch(_site->particle){ + case ParticleType::neutron: + p.pdgcode=2112; + break; + case ParticleType::photon: + p.pdgcode=22; + break; + case ParticleType::electron: + p.pdgcode=11; + break; + case ParticleType::positron: + p.pdgcode=-11; + break; + } + + mcpl_add_particle(file_id,&p); + } + } + } else { +#ifdef OPENMC_MPI + MPI_Send(source_bank->data(), count_size, mpi::source_site, 0, mpi::rank, + mpi::intracomm); +#endif + } + +} + + // Determine member names of a compound HDF5 datatype std::string dtype_member_names(hid_t dtype_id) { From be09b442fbfc41d029ed0eccef87ea25477f40aa Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 6 Sep 2022 09:11:41 +0200 Subject: [PATCH 33/62] allow electrons and positrons as well --- src/source.cpp | 29 ++++++++++++++++++++--------- 1 file changed, 20 insertions(+), 9 deletions(-) diff --git a/src/source.cpp b/src/source.cpp index 82978f1030..e7c944cb58 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -383,19 +383,30 @@ SourceSite MCPLFileSource::read_single_particle() const { SourceSite omc_particle_; const mcpl_particle_t *mcpl_particle; - //extract particle from mcpl-file + // extract particle from mcpl-file mcpl_particle=mcpl_read(mcpl_file); - // check if it is a neutron or a photon. otherwise skip - while ( mcpl_particle->pdgcode!=2112 && mcpl_particle->pdgcode!=22 ) { + // check if it is a neutron, photon, electron, or positron. Otherwise skip. + int pdg=mcpl_particle->pdgcode; + while ( pdg!=2112 && pdg!=22 && pdg!=11 && pdg!=-11) { mcpl_particle=mcpl_read(mcpl_file); - //should check for file exhaustion This could happen if particles are other than - //neutrons or photons + pdg=mcpl_particle->pdgcode; + //should check for file exhaustion. This could happen if particles are other than + //neutrons, photons, electrons, or positrons. } - if(mcpl_particle->pdgcode==2112) { - omc_particle_.particle=ParticleType::neutron; - } else if (mcpl_particle->pdgcode==22) { - omc_particle_.particle=ParticleType::photon; + switch(pdg){ + case 2112: + omc_particle_.particle=ParticleType::neutron; + break; + case 22: + omc_particle_.particle=ParticleType::photon; + break; + case 11: + omc_particle_.particle=ParticleType::electron; + break; + case -11: + omc_particle_.particle=ParticleType::positron; + break; } //particle is good, convert to openmc-formalism From 4e5ba4b4601209f389c394cd3747891e76b43f1c Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 6 Sep 2022 09:12:01 +0200 Subject: [PATCH 34/62] expose a python function to check for mcpl --- openmc/lib/__init__.py | 3 +++ src/source.cpp | 6 ++++++ 2 files changed, 9 insertions(+) diff --git a/openmc/lib/__init__.py b/openmc/lib/__init__.py index c14b0d9c28..1e337fc073 100644 --- a/openmc/lib/__init__.py +++ b/openmc/lib/__init__.py @@ -48,6 +48,9 @@ def _coord_levels(): def _libmesh_enabled(): return c_bool.in_dll(_dll, "LIBMESH_ENABLED").value +def _mcpl_enabled(): + return c_bool.in_dll(_dll, "MCPL_ENABLED").value + from .error import * from .core import * from .nuclide import * diff --git a/src/source.cpp b/src/source.cpp index e7c944cb58..9b426c3cfa 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -41,6 +41,12 @@ namespace openmc { // Global variables //============================================================================== +#ifdef OPENMC_MCPL +const bool MCPL_ENABLED = true; +#else +const bool MCPL_ENABLED = false; +#endif + namespace model { vector> external_sources; From 32e71395b3ad44ecc7080ba39d4b920adfb88702 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 6 Sep 2022 23:06:19 +0200 Subject: [PATCH 35/62] add a trigger to simulation control --- src/simulation.cpp | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/src/simulation.cpp b/src/simulation.cpp index b70535c330..05c8c04508 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -406,6 +406,13 @@ void finalize_batch() auto filename = settings::path_output + "surface_source.h5"; write_source_point(filename.c_str(), true); } + + if(settings::surf_mcpl_write && simulation::current_batch == settings::n_batches){ + auto filename = settings::path_output + ".mcpl"; + write_mcplt_source_point(filename.c_str(), true); + } + + } void initialize_generation() From ab665fb4563be28b6cfd365634049c90f133f485 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Fri, 15 Jul 2022 21:35:59 +0200 Subject: [PATCH 36/62] fix typo Co-authored-by: Paul Romano --- CMakeLists.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 6af9643cdb..be56f1c12f 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -36,7 +36,7 @@ option(OPENMC_ENABLE_COVERAGE "Compile with coverage analysis flags" option(OPENMC_USE_DAGMC "Enable support for DAGMC (CAD) geometry" OFF) option(OPENMC_USE_LIBMESH "Enable support for libMesh unstructured mesh tallies" OFF) option(OPENMC_USE_MPI "Enable MPI" OFF) -option(OPENMC_USE_MCPL "Enable MPCPL" OFF) +option(OPENMC_USE_MCPL "Enable MCPL" OFF) #=============================================================================== # Set a default build configuration if not explicitly specified From 71568d4f4bffa21a2a9c26725619cd180a9a59ff Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 9 Aug 2022 08:12:19 +0200 Subject: [PATCH 37/62] refactor MCPL file read part of FileSource A new constructor is added to FileSource that reads from a given MCPL-file. --- include/openmc/source.h | 30 ++--------- src/settings.cpp | 2 +- src/source.cpp | 116 +++++++++++++++++----------------------- 3 files changed, 53 insertions(+), 95 deletions(-) diff --git a/include/openmc/source.h b/include/openmc/source.h index 4d7ec104aa..5c2dae9311 100644 --- a/include/openmc/source.h +++ b/include/openmc/source.h @@ -96,7 +96,9 @@ class FileSource : public Source { public: // Constructors explicit FileSource(std::string path); - +#ifdef OPENMC_MCPL + explicit FileSource(mcpl_file_t mcpl_file); +#endif // Methods SourceSite sample(uint64_t* seed) const override; @@ -104,32 +106,6 @@ private: vector sites_; //!< Source sites from a file }; -#ifdef OPENMC_MCPL -//============================================================================== -// MCPL-file input source -//============================================================================== -class MCPLFileSource : public Source { -public: - // Constructors, destructors - MCPLFileSource(std::string path); - ~MCPLFileSource(); - - // Methods - //! Sample from the external source distribution - //! \param[inout] seed Pseudorandom seed pointer - //! \return Site read from MCPL-file - SourceSite sample(uint64_t* seed) const; - -private: - SourceSite read_single_particle() const; - void read_source_bank(vector &sites_); - - vector sites_; //!(path)); + model::external_sources.push_back(make_unique(mcpl_open_file(path.c_str()))); #endif } else if (check_for_node(node, "library")) { // Get shared library path and parameters diff --git a/src/source.cpp b/src/source.cpp index 9b426c3cfa..9cdf1632bf 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -283,6 +283,54 @@ FileSource::FileSource(std::string path) file_close(file_id); } +#ifdef OPENMC_MCPL +FileSource::FileSource(mcpl_file_t mcpl_file) +{ + //do checks on the mcpl_file to see if particles are many enough. + // should model this on the example source shown in the docs. + size_t n_sites=mcpl_hdr_nparticles(mcpl_file); + + sites_.resize(n_sites); + for (int i=0;ipdgcode!=2112 && mcpl_particle->pdgcode!=22 ) { + mcpl_particle=mcpl_read(mcpl_file); + //should check for file exhaustion This could happen if particles are other than + //neutrons or photons + } + + if(mcpl_particle->pdgcode==2112) { + site_.particle=ParticleType::neutron; + } else if (mcpl_particle->pdgcode==22) { + site_.particle=ParticleType::photon; + } + + //particle is good, convert to openmc-formalism + site_.r.x=mcpl_particle->position[0]; + site_.r.y=mcpl_particle->position[1]; + site_.r.z=mcpl_particle->position[2]; + + site_.u.x=mcpl_particle->direction[0]; + site_.u.y=mcpl_particle->direction[1]; + site_.u.z=mcpl_particle->direction[2]; + + //mcpl stores kinetic energy in MeV + site_.E=mcpl_particle->ekin*1e6; + //mcpl stores time in ms + site_.time=mcpl_particle->time*1e-3; + site_.wgt=mcpl_particle->weight; + site_.delayed_group=0; + sites_[i]=site_; + } + mcpl_close_file(mcpl_file); +} +#endif //OPENMC_MCPL + SourceSite FileSource::sample(uint64_t* seed) const { size_t i_site = sites_.size() * prn(seed); @@ -366,73 +414,7 @@ MCPLFileSource::MCPLFileSource(std::string path) read_source_bank(sites_); } - -MCPLFileSource::~MCPLFileSource(){ - mcpl_close_file(mcpl_file); -} - -SourceSite MCPLFileSource::sample(uint64_t* seed) const -{ - size_t i_site = sites_.size() * prn(seed); - return sites_[i_site]; -} - -void MCPLFileSource::read_source_bank(vector &sites_) -{ - sites_.resize(n_sites); - for (int i=0;ipdgcode; - while ( pdg!=2112 && pdg!=22 && pdg!=11 && pdg!=-11) { - mcpl_particle=mcpl_read(mcpl_file); - pdg=mcpl_particle->pdgcode; - //should check for file exhaustion. This could happen if particles are other than - //neutrons, photons, electrons, or positrons. - } - - switch(pdg){ - case 2112: - omc_particle_.particle=ParticleType::neutron; - break; - case 22: - omc_particle_.particle=ParticleType::photon; - break; - case 11: - omc_particle_.particle=ParticleType::electron; - break; - case -11: - omc_particle_.particle=ParticleType::positron; - break; - } - - //particle is good, convert to openmc-formalism - omc_particle_.r.x=mcpl_particle->position[0]; - omc_particle_.r.y=mcpl_particle->position[1]; - omc_particle_.r.z=mcpl_particle->position[2]; - - omc_particle_.u.x=mcpl_particle->direction[0]; - omc_particle_.u.y=mcpl_particle->direction[1]; - omc_particle_.u.z=mcpl_particle->direction[2]; - - //mcpl stores kinetic energy in MeV - omc_particle_.E=mcpl_particle->ekin*1e6; - //mcpl stores time in ms - omc_particle_.time=mcpl_particle->time*1e-3; - omc_particle_.wgt=mcpl_particle->weight; - omc_particle_.delayed_group=0; - return omc_particle_; -} -#endif //OPENMC_MCPL +#endif //============================================================================== // Non-member functions From 73a391daef86e3a93c0c7180d9ff7f670afff240 Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 6 Sep 2022 02:22:21 +0200 Subject: [PATCH 38/62] wip add MCPL-out to openmc code for writing a source_bank/point compiles but is untested --- include/openmc/state_point.h | 9 +++ src/state_point.cpp | 125 +++++++++++++++++++++++++++++++++++ 2 files changed, 134 insertions(+) diff --git a/include/openmc/state_point.h b/include/openmc/state_point.h index 5ada2bf88d..2d92c935a6 100644 --- a/include/openmc/state_point.h +++ b/include/openmc/state_point.h @@ -9,6 +9,10 @@ #include "openmc/particle.h" #include "openmc/vector.h" +#ifdef OPENMC_MCPL +#include +#endif + namespace openmc { void load_state_point(); @@ -21,5 +25,10 @@ void write_tally_results_nr(hid_t file_id); void restart_set_keff(); void write_unstructured_mesh_results(); +#ifdef OPENMC_MCPL +void write_mcpl_source_point(const char *filename_, bool surf_source_bank = false); +void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank); +#endif + } // namespace openmc #endif // OPENMC_STATE_POINT_H diff --git a/src/state_point.cpp b/src/state_point.cpp index 9ea0459153..babd670a2d 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -28,6 +28,10 @@ #include "openmc/timer.h" #include "openmc/vector.h" +#ifdef OPENMC_MCPL +#include +#endif + namespace openmc { extern "C" int openmc_statepoint_write(const char* filename, bool* write_source) @@ -598,6 +602,45 @@ void write_source_point(const char* filename, bool surf_source_bank) file_close(file_id); } +void write_mcpl_point(const char *filename, bool surf_source_bank) +{ + std::string filename_; + if (filename) { + filename_ = filename; + } else { + // Determine width for zero padding + int w = std::to_string(settings::n_max_batches).size(); + + filename_ = fmt::format("{0}source.{1:0{2}}.mcpl", settings::path_output, + simulation::current_batch, w); + } + + mcpl_outfile_t file_id; + + std::string line; + if (mpi::master) { + // this must be rewritten: file_open is h5-specific + file_id = mcpl_create_outfile(filename_.c_str()); + //write_attribute(file_id, "filetype", "source"); + //write header stuff (oopy in xml-files as binary blobs for instance)) + if (VERSION_DEV){ + line=fmt::format("OpenMC {VERSION_MAJOR}.{VERSION_MINOR}.{VERSION_RELEASE}-development"); + } else { + line=fmt::format("OpenMC {VERSION_MAJOR}.{VERSION_MINOR}.{VERSION_RELEASE}"); + } + mcpl_hdr_set_srcname(file_id,line.c_str()); + } + + write_mcpl_source_bank(file_id, surf_source_bank); + + if (mpi::master) { + //change this - this is h5 specific + mcpl_closeandgzip_outfile(file_id); + } + +} + + void write_source_bank(hid_t group_id, bool surf_source_bank) { hid_t banktype = h5banktype(); @@ -714,6 +757,88 @@ void write_source_bank(hid_t group_id, bool surf_source_bank) H5Tclose(banktype); } +void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) +{ + int64_t dims_size = settings::n_particles; + int64_t count_size = simulation::work_per_rank; + + // Set vectors for source bank and starting bank index of each process + vector* bank_index = &simulation::work_index; + vector* source_bank = &simulation::source_bank; + vector surf_source_index_vector; + vector surf_source_bank_vector; + + if(surf_source_bank) { + surf_source_index_vector = calculate_surf_source_size(); + dims_size = surf_source_index_vector[mpi::n_procs]; + count_size = simulation::surf_source_bank.size(); + + bank_index = &surf_source_index_vector; + + // Copy data in a SharedArray into a vector. + surf_source_bank_vector.resize(count_size); + surf_source_bank_vector.assign(simulation::surf_source_bank.data(), + simulation::surf_source_bank.data() + count_size); + source_bank = &surf_source_bank_vector; + } + + if (mpi::master) { + //write particles from the master node + //loop over the other nodes and receive data - then write those. + for (int i = 0; i < mpi::n_procs; ++i) { +#ifdef OPENMC_MPI + if (i>0) + MPI_Recv(source_bank->data(), count[0], mpi::source_site, i, i, + mpi::intracomm, MPI_STATUS_IGNORE); +#endif + //now write the source_banke data again. + for (vector::iterator _site=source_bank->begin(); _site!=source_bank->end();_site++){ + //particle is now at the iterator + //write it to the mcpl-file + mcpl_particle_t p; + p.position[0]=_site->r.x; + p.position[1]=_site->r.y; + p.position[2]=_site->r.z; + + p.direction[0]=_site->u.x; + p.direction[1]=_site->u.y; + p.direction[2]=_site->u.x; + + //mcpl stores kinetic energy in MeV + p.ekin=_site->E*1e-6; + + p.time=_site->time*1e3; + + p.weight=_site->wgt; + + switch(_site->particle){ + case ParticleType::neutron: + p.pdgcode=2112; + break; + case ParticleType::photon: + p.pdgcode=22; + break; + case ParticleType::electron: + p.pdgcode=11; + break; + case ParticleType::positron: + p.pdgcode=-11; + break; + } + + mcpl_add_particle(file_id,&p); + } + } + } else { +#ifdef OPENMC_MPI + MPI_Send(source_bank->data(), count_size, mpi::source_site, 0, mpi::rank, + mpi::intracomm); +#endif + } + +} + + // Determine member names of a compound HDF5 datatype std::string dtype_member_names(hid_t dtype_id) { From d07223bf681e5d6e1d7cbde0dac6f97059a12767 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 6 Sep 2022 23:06:19 +0200 Subject: [PATCH 39/62] add a trigger to simulation control --- src/simulation.cpp | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/src/simulation.cpp b/src/simulation.cpp index b70535c330..05c8c04508 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -406,6 +406,13 @@ void finalize_batch() auto filename = settings::path_output + "surface_source.h5"; write_source_point(filename.c_str(), true); } + + if(settings::surf_mcpl_write && simulation::current_batch == settings::n_batches){ + auto filename = settings::path_output + ".mcpl"; + write_mcplt_source_point(filename.c_str(), true); + } + + } void initialize_generation() From 55f2be19e6c0ef51d7a7946ac15a9f2df51da6da Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Fri, 16 Sep 2022 16:07:22 +0200 Subject: [PATCH 40/62] refactor MCPL file read part of FileSource A new constructor is added to FileSource that reads from a given MCPL-file. --- src/source.cpp | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/src/source.cpp b/src/source.cpp index 626d02f9d7..21522438cf 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -426,6 +426,7 @@ MCPLFileSource::MCPLFileSource(std::string path) read_source_bank(sites_); } +<<<<<<< HEAD MCPLFileSource::~MCPLFileSource(){ mcpl_close_file(mcpl_file); @@ -495,6 +496,9 @@ SourceSite MCPLFileSource::read_single_particle() const #endif //OPENMC_MCPL ======= >>>>>>> 380cfc150b449ef56debb5261caa010a52c42bdd +======= +#endif +>>>>>>> 71568d4f4 (refactor MCPL file read part of FileSource) //============================================================================== // Non-member functions From f0af56b181555f52d05f245240846f4ce8eadcaa Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 6 Sep 2022 02:22:21 +0200 Subject: [PATCH 41/62] wip add MCPL-out to openmc code for writing a source_bank/point compiles but is untested --- include/openmc/state_point.h | 9 +++ src/state_point.cpp | 125 +++++++++++++++++++++++++++++++++++ 2 files changed, 134 insertions(+) diff --git a/include/openmc/state_point.h b/include/openmc/state_point.h index 5ada2bf88d..2d92c935a6 100644 --- a/include/openmc/state_point.h +++ b/include/openmc/state_point.h @@ -9,6 +9,10 @@ #include "openmc/particle.h" #include "openmc/vector.h" +#ifdef OPENMC_MCPL +#include +#endif + namespace openmc { void load_state_point(); @@ -21,5 +25,10 @@ void write_tally_results_nr(hid_t file_id); void restart_set_keff(); void write_unstructured_mesh_results(); +#ifdef OPENMC_MCPL +void write_mcpl_source_point(const char *filename_, bool surf_source_bank = false); +void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank); +#endif + } // namespace openmc #endif // OPENMC_STATE_POINT_H diff --git a/src/state_point.cpp b/src/state_point.cpp index 9ea0459153..babd670a2d 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -28,6 +28,10 @@ #include "openmc/timer.h" #include "openmc/vector.h" +#ifdef OPENMC_MCPL +#include +#endif + namespace openmc { extern "C" int openmc_statepoint_write(const char* filename, bool* write_source) @@ -598,6 +602,45 @@ void write_source_point(const char* filename, bool surf_source_bank) file_close(file_id); } +void write_mcpl_point(const char *filename, bool surf_source_bank) +{ + std::string filename_; + if (filename) { + filename_ = filename; + } else { + // Determine width for zero padding + int w = std::to_string(settings::n_max_batches).size(); + + filename_ = fmt::format("{0}source.{1:0{2}}.mcpl", settings::path_output, + simulation::current_batch, w); + } + + mcpl_outfile_t file_id; + + std::string line; + if (mpi::master) { + // this must be rewritten: file_open is h5-specific + file_id = mcpl_create_outfile(filename_.c_str()); + //write_attribute(file_id, "filetype", "source"); + //write header stuff (oopy in xml-files as binary blobs for instance)) + if (VERSION_DEV){ + line=fmt::format("OpenMC {VERSION_MAJOR}.{VERSION_MINOR}.{VERSION_RELEASE}-development"); + } else { + line=fmt::format("OpenMC {VERSION_MAJOR}.{VERSION_MINOR}.{VERSION_RELEASE}"); + } + mcpl_hdr_set_srcname(file_id,line.c_str()); + } + + write_mcpl_source_bank(file_id, surf_source_bank); + + if (mpi::master) { + //change this - this is h5 specific + mcpl_closeandgzip_outfile(file_id); + } + +} + + void write_source_bank(hid_t group_id, bool surf_source_bank) { hid_t banktype = h5banktype(); @@ -714,6 +757,88 @@ void write_source_bank(hid_t group_id, bool surf_source_bank) H5Tclose(banktype); } +void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) +{ + int64_t dims_size = settings::n_particles; + int64_t count_size = simulation::work_per_rank; + + // Set vectors for source bank and starting bank index of each process + vector* bank_index = &simulation::work_index; + vector* source_bank = &simulation::source_bank; + vector surf_source_index_vector; + vector surf_source_bank_vector; + + if(surf_source_bank) { + surf_source_index_vector = calculate_surf_source_size(); + dims_size = surf_source_index_vector[mpi::n_procs]; + count_size = simulation::surf_source_bank.size(); + + bank_index = &surf_source_index_vector; + + // Copy data in a SharedArray into a vector. + surf_source_bank_vector.resize(count_size); + surf_source_bank_vector.assign(simulation::surf_source_bank.data(), + simulation::surf_source_bank.data() + count_size); + source_bank = &surf_source_bank_vector; + } + + if (mpi::master) { + //write particles from the master node + //loop over the other nodes and receive data - then write those. + for (int i = 0; i < mpi::n_procs; ++i) { +#ifdef OPENMC_MPI + if (i>0) + MPI_Recv(source_bank->data(), count[0], mpi::source_site, i, i, + mpi::intracomm, MPI_STATUS_IGNORE); +#endif + //now write the source_banke data again. + for (vector::iterator _site=source_bank->begin(); _site!=source_bank->end();_site++){ + //particle is now at the iterator + //write it to the mcpl-file + mcpl_particle_t p; + p.position[0]=_site->r.x; + p.position[1]=_site->r.y; + p.position[2]=_site->r.z; + + p.direction[0]=_site->u.x; + p.direction[1]=_site->u.y; + p.direction[2]=_site->u.x; + + //mcpl stores kinetic energy in MeV + p.ekin=_site->E*1e-6; + + p.time=_site->time*1e3; + + p.weight=_site->wgt; + + switch(_site->particle){ + case ParticleType::neutron: + p.pdgcode=2112; + break; + case ParticleType::photon: + p.pdgcode=22; + break; + case ParticleType::electron: + p.pdgcode=11; + break; + case ParticleType::positron: + p.pdgcode=-11; + break; + } + + mcpl_add_particle(file_id,&p); + } + } + } else { +#ifdef OPENMC_MPI + MPI_Send(source_bank->data(), count_size, mpi::source_site, 0, mpi::rank, + mpi::intracomm); +#endif + } + +} + + // Determine member names of a compound HDF5 datatype std::string dtype_member_names(hid_t dtype_id) { From 04c8ffc81c7ed7ced0c3525b5ec13f15f8fc079c Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 6 Sep 2022 23:06:19 +0200 Subject: [PATCH 42/62] add a trigger to simulation control --- src/simulation.cpp | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/src/simulation.cpp b/src/simulation.cpp index b70535c330..05c8c04508 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -406,6 +406,13 @@ void finalize_batch() auto filename = settings::path_output + "surface_source.h5"; write_source_point(filename.c_str(), true); } + + if(settings::surf_mcpl_write && simulation::current_batch == settings::n_batches){ + auto filename = settings::path_output + ".mcpl"; + write_mcplt_source_point(filename.c_str(), true); + } + + } void initialize_generation() From 314f3dbbab75104f217d1e20e2146ced28108707 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 4 Oct 2022 10:11:06 +0200 Subject: [PATCH 43/62] do not define these functions if no MCPL present --- src/state_point.cpp | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/src/state_point.cpp b/src/state_point.cpp index babd670a2d..5b00f302dc 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -602,6 +602,7 @@ void write_source_point(const char* filename, bool surf_source_bank) file_close(file_id); } +#ifdef OPENMC_MCPL void write_mcpl_point(const char *filename, bool surf_source_bank) { std::string filename_; @@ -639,7 +640,7 @@ void write_mcpl_point(const char *filename, bool surf_source_bank) } } - +#endif void write_source_bank(hid_t group_id, bool surf_source_bank) { @@ -757,6 +758,7 @@ void write_source_bank(hid_t group_id, bool surf_source_bank) H5Tclose(banktype); } +#ifdef OPENMC_MCPL void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) { int64_t dims_size = settings::n_particles; @@ -837,7 +839,7 @@ void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) } } - +#endif // Determine member names of a compound HDF5 datatype std::string dtype_member_names(hid_t dtype_id) From 2de20231ad584fa7a77a5b135c82f6e995671bcc Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 5 Oct 2022 10:56:36 +0200 Subject: [PATCH 44/62] use find_package instead - MCPL now supports this --- CMakeLists.txt | 20 ++++++++++++-------- 1 file changed, 12 insertions(+), 8 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index be56f1c12f..98dc4fcaf0 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -128,6 +128,13 @@ if(${HDF5_VERSION} VERSION_GREATER_EQUAL 1.12.0) list(APPEND cxxflags -DH5Oget_info_by_idx_vers=1 -DH5O_info_t_vers=1) endif() +#=============================================================================== +# MCPL +#=============================================================================== +if (OPENMC_USE_MCPL) + find_package(MCPL REQUIRED) +endif() + #=============================================================================== # Set compile/link flags based on which compiler is being used #=============================================================================== @@ -155,10 +162,6 @@ if(OPENMC_ENABLE_COVERAGE) list(APPEND ldflags --coverage) endif() -if(OPENMC_USE_MCPL) - list(APPEND ldflags -lmcpl) -endif() - # Show flags being used message(STATUS "OpenMC C++ flags: ${cxxflags}") message(STATUS "OpenMC Linker flags: ${ldflags}") @@ -413,10 +416,6 @@ endif() if (OPENMC_USE_MPI) target_compile_definitions(libopenmc PUBLIC -DOPENMC_MPI) endif() -if(OPENMC_USE_MCPL) - target_compile_definitions(libopenmc PUBLIC -DOPENMC_MCPL) -endif() - # Set git SHA1 hash as a compile definition if(GIT_FOUND) @@ -460,6 +459,11 @@ if (OPENMC_USE_MPI) target_link_libraries(libopenmc MPI::MPI_CXX) endif() +if (OPENMC_USE_MCPL) + target_compile_definitions(libopenmc PUBLIC OPENMC_MCPL) + target_link_libraries(libopenmc MCPL::mcpl) +endif() + #=============================================================================== # openmc executable #=============================================================================== From 7418cb99dd7c22e80532dcc46a9d2c4c04ceb63d Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 5 Oct 2022 10:59:02 +0200 Subject: [PATCH 45/62] add mcpl-output settings --- include/openmc/settings.h | 1 + src/settings.cpp | 8 +++++++- 2 files changed, 8 insertions(+), 1 deletion(-) diff --git a/include/openmc/settings.h b/include/openmc/settings.h index 1e061b235b..66825867d8 100644 --- a/include/openmc/settings.h +++ b/include/openmc/settings.h @@ -48,6 +48,7 @@ extern bool source_latest; //!< write latest source at each batch? extern bool source_separate; //!< write source to separate file? extern bool source_write; //!< write source in HDF5 files? extern bool surf_source_write; //!< write surface source file? +extern bool surf_mcpl_write; //!< write surface mcpl file? extern bool surf_source_read; //!< read surface source file? extern bool survival_biasing; //!< use survival biasing? extern bool temperature_multipole; //!< use multipole data? diff --git a/src/settings.cpp b/src/settings.cpp index 11dd96e025..d1762596c6 100644 --- a/src/settings.cpp +++ b/src/settings.cpp @@ -61,6 +61,7 @@ bool source_latest {false}; bool source_separate {false}; bool source_write {true}; bool surf_source_write {false}; +bool surf_mcpl_write {false}; bool surf_source_read {false}; bool survival_biasing {false}; bool temperature_multipole {false}; @@ -682,9 +683,14 @@ void read_settings_xml() max_surface_particles = std::stoll(get_node_value(node_ssw, "max_particles")); } +#ifdef OPENMC_MCPL + if(check_for_node(node_ssw, "mcpl")){ + surf_mcpl_write=true; + } +#endif } - // If source is not seperate and is to be written out in the statepoint file, + // If source is not separate and is to be written out in the statepoint file, // make sure that the sourcepoint batch numbers are contained in the // statepoint list if (!source_separate) { From 172c518f0e161d0176385c67d3717ad9fb229038 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 5 Oct 2022 10:59:28 +0200 Subject: [PATCH 46/62] fix typo --- src/simulation.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/simulation.cpp b/src/simulation.cpp index 05c8c04508..353ef99265 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -409,7 +409,7 @@ void finalize_batch() if(settings::surf_mcpl_write && simulation::current_batch == settings::n_batches){ auto filename = settings::path_output + ".mcpl"; - write_mcplt_source_point(filename.c_str(), true); + write_mcpl_source_point(filename.c_str(), true); } From 83e2921034f5c448d500b93634999cdb3f350e0b Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 5 Oct 2022 11:00:39 +0200 Subject: [PATCH 47/62] add cond. compilation guards for MCPL functions --- src/state_point.cpp | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/src/state_point.cpp b/src/state_point.cpp index babd670a2d..7c1a454476 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -602,7 +602,8 @@ void write_source_point(const char* filename, bool surf_source_bank) file_close(file_id); } -void write_mcpl_point(const char *filename, bool surf_source_bank) +#ifdef OPENMC_MCPL +void write_mcpl_source_point(const char *filename, bool surf_source_bank) { std::string filename_; if (filename) { @@ -639,7 +640,7 @@ void write_mcpl_point(const char *filename, bool surf_source_bank) } } - +#endif void write_source_bank(hid_t group_id, bool surf_source_bank) { @@ -757,6 +758,7 @@ void write_source_bank(hid_t group_id, bool surf_source_bank) H5Tclose(banktype); } +#ifdef OPENMC_MCPL void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) { int64_t dims_size = settings::n_particles; @@ -837,7 +839,7 @@ void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) } } - +#endif // Determine member names of a compound HDF5 datatype std::string dtype_member_names(hid_t dtype_id) From 5e82bb2e1f1b2857375566f98669be6c93b22638 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 6 Sep 2022 09:12:01 +0200 Subject: [PATCH 48/62] expose a python function to check for mcpl --- openmc/lib/__init__.py | 3 +++ src/source.cpp | 6 ++++++ 2 files changed, 9 insertions(+) diff --git a/openmc/lib/__init__.py b/openmc/lib/__init__.py index c14b0d9c28..1e337fc073 100644 --- a/openmc/lib/__init__.py +++ b/openmc/lib/__init__.py @@ -48,6 +48,9 @@ def _coord_levels(): def _libmesh_enabled(): return c_bool.in_dll(_dll, "LIBMESH_ENABLED").value +def _mcpl_enabled(): + return c_bool.in_dll(_dll, "MCPL_ENABLED").value + from .error import * from .core import * from .nuclide import * diff --git a/src/source.cpp b/src/source.cpp index 49e51c4b16..eace931d63 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -41,6 +41,12 @@ namespace openmc { // Global variables //============================================================================== +#ifdef OPENMC_MCPL +const bool MCPL_ENABLED = true; +#else +const bool MCPL_ENABLED = false; +#endif + namespace model { vector> external_sources; From 6038555def28547a0f392f0eeded8d8d49269d51 Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 6 Sep 2022 02:22:21 +0200 Subject: [PATCH 49/62] wip add MCPL-out to openmc code for writing a source_bank/point compiles but is untested --- include/openmc/state_point.h | 9 +++ src/state_point.cpp | 125 +++++++++++++++++++++++++++++++++++ 2 files changed, 134 insertions(+) diff --git a/include/openmc/state_point.h b/include/openmc/state_point.h index 5ada2bf88d..2d92c935a6 100644 --- a/include/openmc/state_point.h +++ b/include/openmc/state_point.h @@ -9,6 +9,10 @@ #include "openmc/particle.h" #include "openmc/vector.h" +#ifdef OPENMC_MCPL +#include +#endif + namespace openmc { void load_state_point(); @@ -21,5 +25,10 @@ void write_tally_results_nr(hid_t file_id); void restart_set_keff(); void write_unstructured_mesh_results(); +#ifdef OPENMC_MCPL +void write_mcpl_source_point(const char *filename_, bool surf_source_bank = false); +void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank); +#endif + } // namespace openmc #endif // OPENMC_STATE_POINT_H diff --git a/src/state_point.cpp b/src/state_point.cpp index 9ea0459153..babd670a2d 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -28,6 +28,10 @@ #include "openmc/timer.h" #include "openmc/vector.h" +#ifdef OPENMC_MCPL +#include +#endif + namespace openmc { extern "C" int openmc_statepoint_write(const char* filename, bool* write_source) @@ -598,6 +602,45 @@ void write_source_point(const char* filename, bool surf_source_bank) file_close(file_id); } +void write_mcpl_point(const char *filename, bool surf_source_bank) +{ + std::string filename_; + if (filename) { + filename_ = filename; + } else { + // Determine width for zero padding + int w = std::to_string(settings::n_max_batches).size(); + + filename_ = fmt::format("{0}source.{1:0{2}}.mcpl", settings::path_output, + simulation::current_batch, w); + } + + mcpl_outfile_t file_id; + + std::string line; + if (mpi::master) { + // this must be rewritten: file_open is h5-specific + file_id = mcpl_create_outfile(filename_.c_str()); + //write_attribute(file_id, "filetype", "source"); + //write header stuff (oopy in xml-files as binary blobs for instance)) + if (VERSION_DEV){ + line=fmt::format("OpenMC {VERSION_MAJOR}.{VERSION_MINOR}.{VERSION_RELEASE}-development"); + } else { + line=fmt::format("OpenMC {VERSION_MAJOR}.{VERSION_MINOR}.{VERSION_RELEASE}"); + } + mcpl_hdr_set_srcname(file_id,line.c_str()); + } + + write_mcpl_source_bank(file_id, surf_source_bank); + + if (mpi::master) { + //change this - this is h5 specific + mcpl_closeandgzip_outfile(file_id); + } + +} + + void write_source_bank(hid_t group_id, bool surf_source_bank) { hid_t banktype = h5banktype(); @@ -714,6 +757,88 @@ void write_source_bank(hid_t group_id, bool surf_source_bank) H5Tclose(banktype); } +void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) +{ + int64_t dims_size = settings::n_particles; + int64_t count_size = simulation::work_per_rank; + + // Set vectors for source bank and starting bank index of each process + vector* bank_index = &simulation::work_index; + vector* source_bank = &simulation::source_bank; + vector surf_source_index_vector; + vector surf_source_bank_vector; + + if(surf_source_bank) { + surf_source_index_vector = calculate_surf_source_size(); + dims_size = surf_source_index_vector[mpi::n_procs]; + count_size = simulation::surf_source_bank.size(); + + bank_index = &surf_source_index_vector; + + // Copy data in a SharedArray into a vector. + surf_source_bank_vector.resize(count_size); + surf_source_bank_vector.assign(simulation::surf_source_bank.data(), + simulation::surf_source_bank.data() + count_size); + source_bank = &surf_source_bank_vector; + } + + if (mpi::master) { + //write particles from the master node + //loop over the other nodes and receive data - then write those. + for (int i = 0; i < mpi::n_procs; ++i) { +#ifdef OPENMC_MPI + if (i>0) + MPI_Recv(source_bank->data(), count[0], mpi::source_site, i, i, + mpi::intracomm, MPI_STATUS_IGNORE); +#endif + //now write the source_banke data again. + for (vector::iterator _site=source_bank->begin(); _site!=source_bank->end();_site++){ + //particle is now at the iterator + //write it to the mcpl-file + mcpl_particle_t p; + p.position[0]=_site->r.x; + p.position[1]=_site->r.y; + p.position[2]=_site->r.z; + + p.direction[0]=_site->u.x; + p.direction[1]=_site->u.y; + p.direction[2]=_site->u.x; + + //mcpl stores kinetic energy in MeV + p.ekin=_site->E*1e-6; + + p.time=_site->time*1e3; + + p.weight=_site->wgt; + + switch(_site->particle){ + case ParticleType::neutron: + p.pdgcode=2112; + break; + case ParticleType::photon: + p.pdgcode=22; + break; + case ParticleType::electron: + p.pdgcode=11; + break; + case ParticleType::positron: + p.pdgcode=-11; + break; + } + + mcpl_add_particle(file_id,&p); + } + } + } else { +#ifdef OPENMC_MPI + MPI_Send(source_bank->data(), count_size, mpi::source_site, 0, mpi::rank, + mpi::intracomm); +#endif + } + +} + + // Determine member names of a compound HDF5 datatype std::string dtype_member_names(hid_t dtype_id) { From 2ba7f8d1882cf5e2562b473e233b4ca84ab81080 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 6 Sep 2022 23:06:19 +0200 Subject: [PATCH 50/62] add a trigger to simulation control --- src/simulation.cpp | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/src/simulation.cpp b/src/simulation.cpp index b70535c330..05c8c04508 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -406,6 +406,13 @@ void finalize_batch() auto filename = settings::path_output + "surface_source.h5"; write_source_point(filename.c_str(), true); } + + if(settings::surf_mcpl_write && simulation::current_batch == settings::n_batches){ + auto filename = settings::path_output + ".mcpl"; + write_mcplt_source_point(filename.c_str(), true); + } + + } void initialize_generation() From 91cd304f0e4fd32eb6df0bdb75dd3f7d169013ce Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 5 Oct 2022 10:56:36 +0200 Subject: [PATCH 51/62] use find_package instead - MCPL now supports this --- CMakeLists.txt | 20 ++++++++++++-------- 1 file changed, 12 insertions(+), 8 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index be56f1c12f..98dc4fcaf0 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -128,6 +128,13 @@ if(${HDF5_VERSION} VERSION_GREATER_EQUAL 1.12.0) list(APPEND cxxflags -DH5Oget_info_by_idx_vers=1 -DH5O_info_t_vers=1) endif() +#=============================================================================== +# MCPL +#=============================================================================== +if (OPENMC_USE_MCPL) + find_package(MCPL REQUIRED) +endif() + #=============================================================================== # Set compile/link flags based on which compiler is being used #=============================================================================== @@ -155,10 +162,6 @@ if(OPENMC_ENABLE_COVERAGE) list(APPEND ldflags --coverage) endif() -if(OPENMC_USE_MCPL) - list(APPEND ldflags -lmcpl) -endif() - # Show flags being used message(STATUS "OpenMC C++ flags: ${cxxflags}") message(STATUS "OpenMC Linker flags: ${ldflags}") @@ -413,10 +416,6 @@ endif() if (OPENMC_USE_MPI) target_compile_definitions(libopenmc PUBLIC -DOPENMC_MPI) endif() -if(OPENMC_USE_MCPL) - target_compile_definitions(libopenmc PUBLIC -DOPENMC_MCPL) -endif() - # Set git SHA1 hash as a compile definition if(GIT_FOUND) @@ -460,6 +459,11 @@ if (OPENMC_USE_MPI) target_link_libraries(libopenmc MPI::MPI_CXX) endif() +if (OPENMC_USE_MCPL) + target_compile_definitions(libopenmc PUBLIC OPENMC_MCPL) + target_link_libraries(libopenmc MCPL::mcpl) +endif() + #=============================================================================== # openmc executable #=============================================================================== From 9b58c2f6dc93c4c51ce6ee7a26e7ffc645434020 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 5 Oct 2022 10:59:02 +0200 Subject: [PATCH 52/62] add mcpl-output settings --- include/openmc/settings.h | 1 + src/settings.cpp | 8 +++++++- 2 files changed, 8 insertions(+), 1 deletion(-) diff --git a/include/openmc/settings.h b/include/openmc/settings.h index 1e061b235b..66825867d8 100644 --- a/include/openmc/settings.h +++ b/include/openmc/settings.h @@ -48,6 +48,7 @@ extern bool source_latest; //!< write latest source at each batch? extern bool source_separate; //!< write source to separate file? extern bool source_write; //!< write source in HDF5 files? extern bool surf_source_write; //!< write surface source file? +extern bool surf_mcpl_write; //!< write surface mcpl file? extern bool surf_source_read; //!< read surface source file? extern bool survival_biasing; //!< use survival biasing? extern bool temperature_multipole; //!< use multipole data? diff --git a/src/settings.cpp b/src/settings.cpp index 11dd96e025..d1762596c6 100644 --- a/src/settings.cpp +++ b/src/settings.cpp @@ -61,6 +61,7 @@ bool source_latest {false}; bool source_separate {false}; bool source_write {true}; bool surf_source_write {false}; +bool surf_mcpl_write {false}; bool surf_source_read {false}; bool survival_biasing {false}; bool temperature_multipole {false}; @@ -682,9 +683,14 @@ void read_settings_xml() max_surface_particles = std::stoll(get_node_value(node_ssw, "max_particles")); } +#ifdef OPENMC_MCPL + if(check_for_node(node_ssw, "mcpl")){ + surf_mcpl_write=true; + } +#endif } - // If source is not seperate and is to be written out in the statepoint file, + // If source is not separate and is to be written out in the statepoint file, // make sure that the sourcepoint batch numbers are contained in the // statepoint list if (!source_separate) { From 43538ac2021f36caf0d72e002a119a806aa14b71 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 5 Oct 2022 10:59:28 +0200 Subject: [PATCH 53/62] fix typo --- src/simulation.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/simulation.cpp b/src/simulation.cpp index 05c8c04508..353ef99265 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -409,7 +409,7 @@ void finalize_batch() if(settings::surf_mcpl_write && simulation::current_batch == settings::n_batches){ auto filename = settings::path_output + ".mcpl"; - write_mcplt_source_point(filename.c_str(), true); + write_mcpl_source_point(filename.c_str(), true); } From cba745b28b6e384311393c177456b518b60b4f55 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Wed, 5 Oct 2022 11:00:39 +0200 Subject: [PATCH 54/62] add cond. compilation guards for MCPL functions --- src/state_point.cpp | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/src/state_point.cpp b/src/state_point.cpp index babd670a2d..7c1a454476 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -602,7 +602,8 @@ void write_source_point(const char* filename, bool surf_source_bank) file_close(file_id); } -void write_mcpl_point(const char *filename, bool surf_source_bank) +#ifdef OPENMC_MCPL +void write_mcpl_source_point(const char *filename, bool surf_source_bank) { std::string filename_; if (filename) { @@ -639,7 +640,7 @@ void write_mcpl_point(const char *filename, bool surf_source_bank) } } - +#endif void write_source_bank(hid_t group_id, bool surf_source_bank) { @@ -757,6 +758,7 @@ void write_source_bank(hid_t group_id, bool surf_source_bank) H5Tclose(banktype); } +#ifdef OPENMC_MCPL void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) { int64_t dims_size = settings::n_particles; @@ -837,7 +839,7 @@ void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) } } - +#endif // Determine member names of a compound HDF5 datatype std::string dtype_member_names(hid_t dtype_id) From 179df586853b65600518a0ea62526cc2193e0206 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Fri, 16 Sep 2022 16:07:22 +0200 Subject: [PATCH 55/62] refactor MCPL file read part of FileSource A new constructor is added to FileSource that reads from a given MCPL-file. --- src/source.cpp | 73 ++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 73 insertions(+) diff --git a/src/source.cpp b/src/source.cpp index eace931d63..bd76969d17 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -391,6 +391,79 @@ CustomSourceWrapper::~CustomSourceWrapper() } +MCPLFileSource::~MCPLFileSource(){ + mcpl_close_file(mcpl_file); +} + +SourceSite MCPLFileSource::sample(uint64_t* seed) const +{ + size_t i_site = sites_.size() * prn(seed); + return sites_[i_site]; +} + +void MCPLFileSource::read_source_bank(vector &sites_) +{ + sites_.resize(n_sites); + for (int i=0;ipdgcode; + while ( pdg!=2112 && pdg!=22 && pdg!=11 && pdg!=-11) { + mcpl_particle=mcpl_read(mcpl_file); + pdg=mcpl_particle->pdgcode; + //should check for file exhaustion. This could happen if particles are other than + //neutrons, photons, electrons, or positrons. + } + + switch(pdg){ + case 2112: + omc_particle_.particle=ParticleType::neutron; + break; + case 22: + omc_particle_.particle=ParticleType::photon; + break; + case 11: + omc_particle_.particle=ParticleType::electron; + break; + case -11: + omc_particle_.particle=ParticleType::positron; + break; + } + + //particle is good, convert to openmc-formalism + omc_particle_.r.x=mcpl_particle->position[0]; + omc_particle_.r.y=mcpl_particle->position[1]; + omc_particle_.r.z=mcpl_particle->position[2]; + + omc_particle_.u.x=mcpl_particle->direction[0]; + omc_particle_.u.y=mcpl_particle->direction[1]; + omc_particle_.u.z=mcpl_particle->direction[2]; + + //mcpl stores kinetic energy in MeV + omc_particle_.E=mcpl_particle->ekin*1e6; + //mcpl stores time in ms + omc_particle_.time=mcpl_particle->time*1e-3; + omc_particle_.wgt=mcpl_particle->weight; + omc_particle_.delayed_group=0; + return omc_particle_; +} +#endif //OPENMC_MCPL +======= +>>>>>>> 380cfc150b449ef56debb5261caa010a52c42bdd +======= +#endif +>>>>>>> 71568d4f4 (refactor MCPL file read part of FileSource) +>>>>>>> 55f2be19e (refactor MCPL file read part of FileSource) + //============================================================================== // Non-member functions //============================================================================== From e3fbeac636c2e9a47599c60e04b2a1ba4d9ba46d Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 6 Sep 2022 02:22:21 +0200 Subject: [PATCH 56/62] wip add MCPL-out to openmc code for writing a source_bank/point compiles but is untested --- src/source.cpp | 6 ------ src/state_point.cpp | 1 + 2 files changed, 1 insertion(+), 6 deletions(-) diff --git a/src/source.cpp b/src/source.cpp index bd76969d17..f205a0a256 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -457,12 +457,6 @@ SourceSite MCPLFileSource::read_single_particle() const return omc_particle_; } #endif //OPENMC_MCPL -======= ->>>>>>> 380cfc150b449ef56debb5261caa010a52c42bdd -======= -#endif ->>>>>>> 71568d4f4 (refactor MCPL file read part of FileSource) ->>>>>>> 55f2be19e (refactor MCPL file read part of FileSource) //============================================================================== // Non-member functions diff --git a/src/state_point.cpp b/src/state_point.cpp index 7c1a454476..31979c0c56 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -642,6 +642,7 @@ void write_mcpl_source_point(const char *filename, bool surf_source_bank) } #endif + void write_source_bank(hid_t group_id, bool surf_source_bank) { hid_t banktype = h5banktype(); From 3b7c4f6d2f141311cba5ad9e5f773b3bc3381852 Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Tue, 4 Oct 2022 10:11:06 +0200 Subject: [PATCH 57/62] do not define these functions if no MCPL present --- src/state_point.cpp | 1 - 1 file changed, 1 deletion(-) diff --git a/src/state_point.cpp b/src/state_point.cpp index 31979c0c56..7c1a454476 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -642,7 +642,6 @@ void write_mcpl_source_point(const char *filename, bool surf_source_bank) } #endif - void write_source_bank(hid_t group_id, bool surf_source_bank) { hid_t banktype = h5banktype(); From 982377b0cdfc573469aad799852a0f2dd50b2e68 Mon Sep 17 00:00:00 2001 From: erkn Date: Wed, 5 Oct 2022 14:31:05 +0200 Subject: [PATCH 58/62] fix misnamed variable and remove accientally committed code --- src/source.cpp | 169 ++----------------------------------------------- 1 file changed, 4 insertions(+), 165 deletions(-) diff --git a/src/source.cpp b/src/source.cpp index e1891c5e94..e35734747f 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -308,16 +308,16 @@ FileSource::FileSource(mcpl_file_t mcpl_file) switch(pdg){ case 2112: - omc_particle_.particle=ParticleType::neutron; + site_.particle=ParticleType::neutron; break; case 22: - omc_particle_.particle=ParticleType::photon; + site_.particle=ParticleType::photon; break; case 11: - omc_particle_.particle=ParticleType::electron; + site_.particle=ParticleType::electron; break; case -11: - omc_particle_.particle=ParticleType::positron; + site_.particle=ParticleType::positron; break; } @@ -401,167 +401,6 @@ CustomSourceWrapper::~CustomSourceWrapper() #endif } -<<<<<<< HEAD -<<<<<<< HEAD -#ifdef OPENMC_MCPL -//=========================================================================== -// Read particles from an MCPL-file -//=========================================================================== -MCPLFileSource::MCPLFileSource(std::string path) -{ - // Check if source file exists - if (!file_exists(path)) { - fatal_error(fmt::format("Source file '{}' does not exist.", path)); - } - - // Read the source from a binary file instead of sampling from some - // assumed source distribution - write_message(6, "Reading mcpl source file from {}",path); - - // Open the mcpl file - mcpl_file = mcpl_open_file(path.c_str()); - - //do checks on the mcpl_file to see if particles are many enough. - // should model this on the example source shown in the docs. - n_sites=mcpl_hdr_nparticles(mcpl_file); - - read_source_bank(sites_); -} -<<<<<<< HEAD - -MCPLFileSource::~MCPLFileSource(){ - mcpl_close_file(mcpl_file); -} - -SourceSite MCPLFileSource::sample(uint64_t* seed) const -{ - size_t i_site = sites_.size() * prn(seed); - return sites_[i_site]; -} - -void MCPLFileSource::read_source_bank(vector &sites_) -{ - sites_.resize(n_sites); - for (int i=0;ipdgcode; - while ( pdg!=2112 && pdg!=22 && pdg!=11 && pdg!=-11) { - mcpl_particle=mcpl_read(mcpl_file); - pdg=mcpl_particle->pdgcode; - //should check for file exhaustion. This could happen if particles are other than - //neutrons, photons, electrons, or positrons. - } - - switch(pdg){ - case 2112: - omc_particle_.particle=ParticleType::neutron; - break; - case 22: - omc_particle_.particle=ParticleType::photon; - break; - case 11: - omc_particle_.particle=ParticleType::electron; - break; - case -11: - omc_particle_.particle=ParticleType::positron; - break; - } - - //particle is good, convert to openmc-formalism - omc_particle_.r.x=mcpl_particle->position[0]; - omc_particle_.r.y=mcpl_particle->position[1]; - omc_particle_.r.z=mcpl_particle->position[2]; - - omc_particle_.u.x=mcpl_particle->direction[0]; - omc_particle_.u.y=mcpl_particle->direction[1]; - omc_particle_.u.z=mcpl_particle->direction[2]; - - //mcpl stores kinetic energy in MeV - omc_particle_.E=mcpl_particle->ekin*1e6; - //mcpl stores time in ms - omc_particle_.time=mcpl_particle->time*1e-3; - omc_particle_.wgt=mcpl_particle->weight; - omc_particle_.delayed_group=0; - return omc_particle_; -} -#endif //OPENMC_MCPL - -MCPLFileSource::~MCPLFileSource(){ - mcpl_close_file(mcpl_file); -} - -SourceSite MCPLFileSource::sample(uint64_t* seed) const -{ - size_t i_site = sites_.size() * prn(seed); - return sites_[i_site]; -} - -void MCPLFileSource::read_source_bank(vector &sites_) -{ - sites_.resize(n_sites); - for (int i=0;ipdgcode; - while ( pdg!=2112 && pdg!=22 && pdg!=11 && pdg!=-11) { - mcpl_particle=mcpl_read(mcpl_file); - pdg=mcpl_particle->pdgcode; - //should check for file exhaustion. This could happen if particles are other than - //neutrons, photons, electrons, or positrons. - } - - switch(pdg){ - case 2112: - omc_particle_.particle=ParticleType::neutron; - break; - case 22: - omc_particle_.particle=ParticleType::photon; - break; - case 11: - omc_particle_.particle=ParticleType::electron; - break; - case -11: - omc_particle_.particle=ParticleType::positron; - break; - } - - //particle is good, convert to openmc-formalism - omc_particle_.r.x=mcpl_particle->position[0]; - omc_particle_.r.y=mcpl_particle->position[1]; - omc_particle_.r.z=mcpl_particle->position[2]; - - omc_particle_.u.x=mcpl_particle->direction[0]; - omc_particle_.u.y=mcpl_particle->direction[1]; - omc_particle_.u.z=mcpl_particle->direction[2]; - - //mcpl stores kinetic energy in MeV - omc_particle_.E=mcpl_particle->ekin*1e6; - //mcpl stores time in ms - omc_particle_.time=mcpl_particle->time*1e-3; - omc_particle_.wgt=mcpl_particle->weight; - omc_particle_.delayed_group=0; - return omc_particle_; -} -#endif //OPENMC_MCPL //============================================================================== // Non-member functions From 12f088a5ce68c639c9772574eb69c4be5de452fe Mon Sep 17 00:00:00 2001 From: Erik B Knudsen Date: Thu, 6 Oct 2022 01:33:21 +0200 Subject: [PATCH 59/62] must have a full filename --- src/simulation.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/simulation.cpp b/src/simulation.cpp index 353ef99265..2c9e7e8566 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -408,7 +408,7 @@ void finalize_batch() } if(settings::surf_mcpl_write && simulation::current_batch == settings::n_batches){ - auto filename = settings::path_output + ".mcpl"; + auto filename = settings::path_output + "surface_source.mcpl"; write_mcpl_source_point(filename.c_str(), true); } From 0706a0eca305c8fc4a0a5cbfb0c5b6e46099753e Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 11 Oct 2022 00:46:32 +0200 Subject: [PATCH 60/62] fix malformed format string --- src/state_point.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/state_point.cpp b/src/state_point.cpp index 7c1a454476..e197ee162a 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -625,9 +625,9 @@ void write_mcpl_source_point(const char *filename, bool surf_source_bank) //write_attribute(file_id, "filetype", "source"); //write header stuff (oopy in xml-files as binary blobs for instance)) if (VERSION_DEV){ - line=fmt::format("OpenMC {VERSION_MAJOR}.{VERSION_MINOR}.{VERSION_RELEASE}-development"); + line=fmt::format("OpenMC {0}.{1}.{2}-development",VERSION_MAJOR,VERSION_MINOR,VERSION_RELEASE); } else { - line=fmt::format("OpenMC {VERSION_MAJOR}.{VERSION_MINOR}.{VERSION_RELEASE}"); + line=fmt::format("OpenMC {0}.{1}.{2}",VERSION_MAJOR,VERSION_MINOR,VERSION_RELEASE); } mcpl_hdr_set_srcname(file_id,line.c_str()); } From aeadaf0b2f1a59fcd2c3f8cb3ad8fa7dab36e342 Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 11 Oct 2022 00:48:44 +0200 Subject: [PATCH 61/62] fix typo causing unit length errors --- src/state_point.cpp | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/state_point.cpp b/src/state_point.cpp index e197ee162a..a8a4be034e 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -802,9 +802,11 @@ void write_mcpl_source_bank(mcpl_outfile_t file_id, bool surf_source_bank) p.position[1]=_site->r.y; p.position[2]=_site->r.z; + //mcpl requires that the direction vector is unit length + //which is also the case in openmc p.direction[0]=_site->u.x; p.direction[1]=_site->u.y; - p.direction[2]=_site->u.x; + p.direction[2]=_site->u.z; //mcpl stores kinetic energy in MeV p.ekin=_site->E*1e-6; From 72251eb30a996e6aaa9672292a1483c69d53a652 Mon Sep 17 00:00:00 2001 From: erkn Date: Tue, 11 Oct 2022 02:00:33 +0200 Subject: [PATCH 62/62] enable mcpl surf_source_write in the python layer --- openmc/settings.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/openmc/settings.py b/openmc/settings.py index 2d81216599..45ff0edc68 100644 --- a/openmc/settings.py +++ b/openmc/settings.py @@ -165,6 +165,7 @@ class Settings: banked (int) :max_particles: Maximum number of particles to be banked on surfaces per process (int) + :mcpl: Output in the form of an MCPL-file (bool) survival_biasing : bool Indicate whether survival biasing is to be used tabular_legendre : dict @@ -643,7 +644,7 @@ class Settings: cv.check_type('surface source writing options', surf_source_write, Mapping) for key, value in surf_source_write.items(): cv.check_value('surface source writing key', key, - ('surface_ids', 'max_particles')) + ('surface_ids', 'max_particles', 'mcpl')) if key == 'surface_ids': cv.check_type('surface ids for source banking', value, Iterable, Integral) @@ -655,6 +656,9 @@ class Settings: value, Integral) cv.check_greater_than('maximum particle banks on surfaces per process', value, 0) + elif key == 'mcpl': + cv.check_type('write to an MCPL-format file', value, bool) + self._surf_source_write = surf_source_write @confidence_intervals.setter @@ -1023,6 +1027,9 @@ class Settings: if 'max_particles' in self._surf_source_write: subelement = ET.SubElement(element, "max_particles") subelement.text = str(self._surf_source_write['max_particles']) + if 'mcpl' in self._surf_source_write: + subelement = ET.SubElement(element, "mcpl") + subelement.text = str(self._surf_source_write['mcpl']).lower() def _create_confidence_intervals(self, root): if self._confidence_intervals is not None: