diff --git a/include/openmc/simulation.h b/include/openmc/simulation.h index 02e434aba..2002e95fe 100644 --- a/include/openmc/simulation.h +++ b/include/openmc/simulation.h @@ -11,6 +11,10 @@ namespace openmc { +constexpr int STATUS_EXIT_NORMAL {0}; +constexpr int STATUS_EXIT_MAX_BATCH {1}; +constexpr int STATUS_EXIT_ON_TRIGGER {2}; + //============================================================================== // Global variable declarations //============================================================================== @@ -55,32 +59,29 @@ extern "C" int thread_id; //!< ID of a given thread //! Determine number of particles to transport per process void calculate_work(); -//! Initialize simulation -extern "C" void openmc_simulation_init_c(); - //! Initialize a batch -extern "C" void initialize_batch(); +void initialize_batch(); //! Initialize a fission generation -extern "C" void initialize_generation(); +void initialize_generation(); -extern "C" void initialize_history(Particle* p, int64_t index_source); +void initialize_history(Particle* p, int64_t index_source); //! Finalize a batch //! //! Handles synchronization and accumulation of tallies, calculation of Shannon //! entropy, getting single-batch estimate of keff, and turning on tallies when //! appropriate -extern "C" void finalize_batch(); +void finalize_batch(); //! Finalize a fission generation -extern "C" void finalize_generation(); +void finalize_generation(); //! Determine overall generation number extern "C" int overall_generation(); #ifdef OPENMC_MPI -extern "C" void broadcast_results(); +void broadcast_results(); #endif } // namespace openmc diff --git a/include/openmc/timer.h b/include/openmc/timer.h index 699e0fcbe..518923f99 100644 --- a/include/openmc/timer.h +++ b/include/openmc/timer.h @@ -47,6 +47,7 @@ extern Timer time_inactive; extern Timer time_initialize; extern Timer time_tallies; extern Timer time_total; +extern Timer time_transport; } // namespace openmc diff --git a/src/api.F90 b/src/api.F90 index fa8b81475..3d757e138 100644 --- a/src/api.F90 +++ b/src/api.F90 @@ -73,7 +73,6 @@ module openmc_api public :: openmc_material_filter_set_bins public :: openmc_mesh_filter_set_mesh public :: openmc_meshsurface_filter_set_mesh - public :: openmc_next_batch public :: openmc_nuclide_name public :: openmc_plot_geometry public :: openmc_reset @@ -271,8 +270,6 @@ contains ! Reset timers call time_read_xs % reset() call time_unionize % reset() - call time_tallies % reset() - call time_transport % reset() call reset_timers() err = 0 diff --git a/src/input_xml.F90 b/src/input_xml.F90 index cb51672da..1c8bcc22a 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -1222,7 +1222,11 @@ contains ! Read derivative attributes. do i = 1, size(node_deriv_list) +!$omp parallel +!$omp critical (ReadTallyDeriv) call tally_derivs(i) % from_xml(node_deriv_list(i)) +!$omp end critical (ReadTallyDeriv) +!$omp end parallel ! Update tally derivative dictionary call tally_deriv_dict % set(tally_derivs(i) % id, i) diff --git a/src/output.F90 b/src/output.F90 index 8b4e2507b..ab37ceb50 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -505,7 +505,7 @@ contains write(ou,100) " Reading cross sections", time_read_xs % elapsed write(ou,100) "Total time in simulation", time_inactive_elapsed() + & time_active_elapsed() - write(ou,100) " Time in transport only", time_transport % elapsed + write(ou,100) " Time in transport only", time_transport_elapsed() if (run_mode == MODE_EIGENVALUE) then write(ou,100) " Time in inactive batches", time_inactive_elapsed() end if @@ -515,7 +515,7 @@ contains write(ou,100) " Sampling source sites", time_bank_sample_elapsed() write(ou,100) " SEND/RECV source sites", time_bank_sendrecv_elapsed() end if - write(ou,100) " Time accumulating tallies", time_tallies % elapsed + write(ou,100) " Time accumulating tallies", time_tallies_elapsed() if (cmfd_run) write(ou,100) " Time in CMFD", time_cmfd % elapsed if (cmfd_run) write(ou,100) " Building matrices", & time_cmfdbuild % elapsed @@ -531,7 +531,7 @@ contains speed_inactive = real(n_particles * (n_inactive - restart_batch) * & gen_per_batch) / time_inactive_elapsed() speed_active = real(n_particles * n_active * gen_per_batch) / & - time_active % elapsed + time_active_elapsed() else speed_inactive = ZERO speed_active = real(n_particles * (n_batches - restart_batch) * & diff --git a/src/simulation.F90 b/src/simulation.F90 index e8f56459c..01119b4e7 100644 --- a/src/simulation.F90 +++ b/src/simulation.F90 @@ -7,143 +7,22 @@ module simulation #endif use bank_header, only: source_bank - use cmfd_execute, only: cmfd_init_batch, cmfd_tally_init, execute_cmfd - use cmfd_header, only: cmfd_on use constants, only: ZERO - use eigenvalue, only: calculate_average_keff, calculate_generation_keff -#ifdef _OPENMP - use eigenvalue, only: join_bank_from_threads -#endif - use error, only: fatal_error, write_message - use geometry_header, only: n_cells + use error, only: fatal_error use material_header, only: n_materials, materials use message_passing - use mgxs_interface, only: energy_bins, energy_bin_avg use nuclide_header, only: micro_xs, n_nuclides - use output, only: print_batch_keff - use particle_header use photon_header, only: micro_photon_xs, n_elements - use random_lcg, only: set_particle_seed use settings use simulation_header - use state_point, only: openmc_statepoint_write, load_state_point - use string, only: to_str - use tally, only: accumulate_tallies, setup_active_tallies, & - init_tally_routines use tally_header use tally_filter_header, only: filter_matches, n_filters - use tally_derivative_header, only: tally_derivs - use timer_header - use trigger, only: check_triggers - use tracking, only: transport implicit none private - public :: openmc_next_batch - - integer(C_INT), parameter :: STATUS_EXIT_NORMAL = 0 - integer(C_INT), parameter :: STATUS_EXIT_MAX_BATCH = 1 - integer(C_INT), parameter :: STATUS_EXIT_ON_TRIGGER = 2 - - interface - subroutine initialize_source() bind(C) - end subroutine - - subroutine initialize_generation() bind(C) - end subroutine - - subroutine finalize_generation() bind(C) - end subroutine finalize_generation - - subroutine finalize_batch() bind(C) - end subroutine - - subroutine initialize_history(p, index_source) bind(C) - import Particle, C_INT64_T - type(Particle), intent(inout) :: p - integer(C_INT64_T), value :: index_source - end subroutine - - function sample_external_source() result(site) bind(C) - import Bank - type(Bank) :: site - end function - end interface contains -!=============================================================================== -! OPENMC_NEXT_BATCH -!=============================================================================== - - function openmc_next_batch(status) result(err) bind(C) - integer(C_INT), intent(out), optional :: status - integer(C_INT) :: err - - type(Particle) :: p - integer(8) :: i_work - interface - subroutine initialize_batch() bind(C) - end subroutine - end interface - - err = 0 - - ! Make sure simulation has been initialized - if (.not. simulation_initialized) then - err = E_ALLOCATE - call set_errmsg("Simulation has not been initialized yet.") - return - end if - - call initialize_batch() - - ! ======================================================================= - ! LOOP OVER GENERATIONS - GENERATION_LOOP: do current_gen = 1, gen_per_batch - - call initialize_generation() - - ! Start timer for transport - call time_transport % start() - - ! ==================================================================== - ! LOOP OVER PARTICLES -!$omp parallel do schedule(runtime) firstprivate(p) copyin(tally_derivs) - PARTICLE_LOOP: do i_work = 1, work - current_work = i_work - - ! grab source particle from bank - call initialize_history(p, current_work) - - ! transport particle - call transport(p) - - end do PARTICLE_LOOP -!$omp end parallel do - - ! Accumulate time for transport - call time_transport % stop() - - call finalize_generation() - - end do GENERATION_LOOP - - call finalize_batch() - - ! Check simulation ending criteria - if (present(status)) then - if (current_batch == n_max_batches) then - status = STATUS_EXIT_MAX_BATCH - elseif (satisfy_triggers) then - status = STATUS_EXIT_ON_TRIGGER - else - status = STATUS_EXIT_NORMAL - end if - end if - - end function openmc_next_batch - !=============================================================================== ! INITIALIZE_SIMULATION !=============================================================================== diff --git a/src/simulation.cpp b/src/simulation.cpp index 2f154b8ce..60b22af9d 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -6,6 +6,7 @@ #include "openmc/error.h" #include "openmc/message_passing.h" #include "openmc/output.h" +#include "openmc/particle.h" #include "openmc/random_lcg.h" #include "openmc/settings.h" #include "openmc/source.h" @@ -16,6 +17,8 @@ #include #include +namespace openmc { + // data/functions from Fortran side extern "C" bool cmfd_on; @@ -37,8 +40,11 @@ extern "C" void print_runtime(); extern "C" void setup_active_tallies(); extern "C" void simulation_init_f(); extern "C" void simulation_finalize_f(); +extern "C" void transport(Particle* p); extern "C" void write_tallies(); +} // namespace openmc + //============================================================================== // C API functions //============================================================================== @@ -161,6 +167,64 @@ int openmc_simulation_finalize() return 0; } +int openmc_next_batch(int* status) +{ + using namespace openmc; + using openmc::simulation::current_gen; + + // Make sure simulation has been initialized + if (!simulation::initialized) { + set_errmsg("Simulation has not been initialized yet."); + return OPENMC_E_ALLOCATE; + } + + initialize_batch(); + + // ======================================================================= + // LOOP OVER GENERATIONS + for (current_gen = 1; current_gen <= settings::gen_per_batch; ++current_gen) { + + initialize_generation(); + + // Start timer for transport + time_transport.start(); + + // ==================================================================== + // LOOP OVER PARTICLES + +#pragma omp parallel for schedule(runtime) + for (int64_t i_work = 1; i_work <= simulation::work; ++i_work) { + simulation::current_work = i_work; + + // grab source particle from bank + Particle p; + initialize_history(&p, simulation::current_work); + + // transport particle + transport(&p); + } + + // Accumulate time for transport + time_transport.stop(); + + finalize_generation(); + } + + finalize_batch(); + + // Check simulation ending criteria + if (status) { + if (simulation::current_batch == settings::n_max_batches) { + *status = STATUS_EXIT_MAX_BATCH; + } else if (simulation::satisfy_triggers) { + *status = STATUS_EXIT_ON_TRIGGER; + } else { + *status = STATUS_EXIT_NORMAL; + } + } + return 0; +} + namespace openmc { //============================================================================== diff --git a/src/state_point.F90 b/src/state_point.F90 index d99bc1b84..0ca89dcc2 100644 --- a/src/state_point.F90 +++ b/src/state_point.F90 @@ -397,7 +397,7 @@ contains call write_dataset(runtime_group, "simulation", & time_inactive_elapsed() + time_active_elapsed()) call write_dataset(runtime_group, "transport", & - time_transport % get_value()) + time_transport_elapsed()) if (run_mode == MODE_EIGENVALUE) then call write_dataset(runtime_group, "inactive batches", & time_inactive_elapsed()) @@ -413,7 +413,7 @@ contains time_bank_sendrecv_elapsed()) end if call write_dataset(runtime_group, "accumulating tallies", & - time_tallies % get_value()) + time_tallies_elapsed()) if (cmfd_run) then call write_dataset(runtime_group, "CMFD", time_cmfd % get_value()) call write_dataset(runtime_group, "CMFD building matrices", & diff --git a/src/timer.cpp b/src/timer.cpp index df8f4fb92..ecf7d3ae9 100644 --- a/src/timer.cpp +++ b/src/timer.cpp @@ -47,6 +47,7 @@ Timer time_inactive; Timer time_initialize; Timer time_tallies; Timer time_total; +Timer time_transport; //============================================================================== // Fortran compatibility @@ -61,6 +62,7 @@ extern "C" double time_inactive_elapsed() { return time_inactive.elapsed(); } extern "C" double time_initialize_elapsed() { return time_initialize.elapsed(); } extern "C" double time_tallies_elapsed() { return time_tallies.elapsed(); } extern "C" double time_total_elapsed() { return time_total.elapsed(); } +extern "C" double time_transport_elapsed() { return time_transport.elapsed(); } extern "C" void reset_timers() { @@ -69,8 +71,11 @@ extern "C" void reset_timers() time_bank_sample.reset(); time_bank_sendrecv.reset(); time_finalize.reset(); + time_inactive.reset(); time_initialize.reset(); + time_tallies.reset(); time_total.reset(); + time_transport.reset(); } } // namespace openmc diff --git a/src/timer_header.F90 b/src/timer_header.F90 index 00af1596f..d0fbc1a6d 100644 --- a/src/timer_header.F90 +++ b/src/timer_header.F90 @@ -43,6 +43,10 @@ module timer_header import C_DOUBLE real(C_DOUBLE) :: t end function + function time_transport_elapsed() result(t) bind(C) + import C_DOUBLE + real(C_DOUBLE) :: t + end function subroutine reset_timers() bind(C) end subroutine end interface @@ -68,14 +72,8 @@ module timer_header ! ============================================================================ ! TIMING VARIABLES - type(Timer) :: time_total ! timer for total run - type(Timer) :: time_initialize ! timer for initialization type(Timer) :: time_read_xs ! timer for reading cross sections type(Timer) :: time_unionize ! timer for material xs-energy grid union - type(Timer) :: time_tallies ! timer for accumulate tallies - type(Timer) :: time_inactive ! timer for inactive batches - type(Timer) :: time_active ! timer for active batches - type(Timer) :: time_transport ! timer for transport only type(Timer) :: time_finalize ! timer for finalization contains diff --git a/src/tracking.F90 b/src/tracking.F90 index d14de4a46..038c7fdb5 100644 --- a/src/tracking.F90 +++ b/src/tracking.F90 @@ -48,7 +48,7 @@ contains ! TRANSPORT encompasses the main logic for moving a particle through geometry. !=============================================================================== - subroutine transport(p) + subroutine transport(p) bind(C) type(Particle), intent(inout) :: p