Completed replacement of physics_mg.F90 with C++

This commit is contained in:
Adam G Nelson 2018-10-12 10:40:35 -04:00
parent ef9908e774
commit e3e388bb6d
10 changed files with 140 additions and 140 deletions

View file

@ -345,7 +345,6 @@ add_library(libopenmc SHARED
src/photon_physics.F90
src/physics_common.F90
src/physics.F90
src/physics_mg.F90
src/plot.F90
src/plot_header.F90
src/progress_header.F90

View file

@ -26,6 +26,7 @@ class Material
public:
int32_t id; //!< Unique ID
double volume_ {-1.0}; //!< Volume in [cm^3]
bool fissionable {false}; //!< Does this material contain fissionable nuclides
//! \brief Default temperature for cells containing this material.
//!

View file

@ -14,19 +14,33 @@ namespace openmc {
// SCATTER
//==============================================================================
//! \brief Samples the scattering event
//! \brief samples particle behavior after a collision event.
extern "C" void
collision_mg(Particle* p, Bank* fission_bank, const int64_t fission_bank_size,
const double* energy_bin_avg, const MaterialMacroXS& material_xs);
//! \brief samples a reaction type.
//!
//! Note that there is special logic when suvival biasing is turned on since
//! fission and disappearance are treated implicitly.
void
sample_reaction(Particle* p, Bank* fission_bank,
const int64_t fission_bank_size, const double* energy_bin_avg,
const MaterialMacroXS& material_xs);
//! \brief Samples the scattering event
void
scatter(Particle* p, const double* energy_bin_avg);
//! \brief Determines the average total, prompt and delayed neutrons produced
//! from fission and creates the appropriate bank sites.
extern "C" void
void
create_fission_sites(Particle* p, Bank* bank_array, int64_t& size_bank,
int64_t& bank_array_size, MaterialMacroXS& material_xs);
const int64_t& bank_array_size, const MaterialMacroXS& material_xs);
//! \brief Handles an absorption event
extern "C" void
absorption(Particle* p, MaterialMacroXS& material_xs);
void
absorption(Particle* p, const MaterialMacroXS& material_xs);
} // namespace openmc
#endif // OPENMC_PHYSICS_MG_H

View file

@ -2986,7 +2986,7 @@ contains
! Check if material is fissionable
if (nuclides(materials(i) % nuclide(j)) % fissionable) then
materials(i) % fissionable = .true.
call materials(i) % set_fissionable(logical(.true., C_BOOL))
end if
end do

View file

@ -122,6 +122,13 @@ extern "C" {
material_map[id] = index - 1;
}
bool material_fissionable(Material* mat) {return mat->fissionable;}
void material_set_fissionable(Material* mat, bool fissionable, int32_t index)
{
mat->fissionable = fissionable;
}
void extend_materials_c(int32_t n)
{
materials.reserve(materials.size() + n);

View file

@ -51,6 +51,20 @@ module material_header
integer(C_INT32_T), intent(in), value :: index
end subroutine material_set_id_c
function material_fissionable_c(mat_ptr) &
bind(C, name='material_fissionable') result(fissionable)
import C_PTR, C_BOOL
type(C_PTR), intent(in), value :: mat_ptr
logical(C_BOOL) :: fissionable
end function material_fissionable_c
subroutine material_set_fissionable_c(mat_ptr, fissionable) &
bind(C, name='material_set_fissionable')
import C_PTR, C_BOOL
type(C_PTR), intent(in), value :: mat_ptr
logical(C_BOOL), intent(in), value :: fissionable
end subroutine material_set_fissionable_c
subroutine extend_materials_c(n) bind(C)
import C_INT32_T
integer(C_INT32_T), intent(in), value :: n
@ -95,7 +109,6 @@ module material_header
character(20), allocatable :: sab_names(:) ! name of S(a,b) table
! Does this material contain fissionable nuclides? Is it depletable?
logical :: fissionable = .false.
logical :: depletable = .false.
! enforce isotropic scattering in lab for specific nuclides
@ -105,6 +118,8 @@ module material_header
contains
procedure :: id => material_id
procedure :: set_id => material_set_id
procedure :: fissionable => material_fissionable
procedure :: set_fissionable => material_set_fissionable
procedure :: set_density => material_set_density
procedure :: init_nuclide_index => material_init_nuclide_index
procedure :: assign_sab_tables => material_assign_sab_tables
@ -139,6 +154,18 @@ contains
call material_set_id_c(this % ptr, id, index)
end subroutine material_set_id
function material_fissionable(this) result(fissionable)
class(Material), intent(in) :: this
logical(C_BOOL) :: fissionable
fissionable = material_fissionable_c(this % ptr)
end function material_fissionable
subroutine material_set_fissionable(this, fissionable)
class(Material), intent(in) :: this
logical(C_BOOL), intent(in) :: fissionable
call material_set_fissionable_c(this % ptr, fissionable)
end subroutine material_set_fissionable
function material_set_density(this, density) result(err)
class(Material), intent(inout) :: this
real(8), intent(in) :: density
@ -685,7 +712,7 @@ contains
integer(C_INT) :: err
if (index >= 1 .and. index <= size(materials)) then
fissionable = materials(index) % fissionable
fissionable = materials(index) % fissionable()
err = 0
else
err = E_OUT_OF_BOUNDS

View file

@ -90,7 +90,8 @@ contains
end if
end do NUCLIDE_LOOP
mat % fissionable = query_fissionable_c(mat % n_nuclides, mat % nuclide)
call mat % set_fissionable(query_fissionable_c(mat % n_nuclides, &
mat % nuclide))
end do MATERIAL_LOOP

View file

@ -1,126 +0,0 @@
module physics_mg
! This module contains the multi-group specific physics routines so as to not
! hinder performance of the CE versions with multiple if-thens.
use bank_header
use constants
use error, only: fatal_error, warning, write_message
use material_header, only: Material, materials
use math, only: rotate_angle
use mgxs_interface
use message_passing
use nuclide_header, only: MaterialMacroXS, material_xs
use particle_header
use physics_common
use random_lcg, only: prn
use settings
use simulation_header
use string, only: to_str
use tally_header
implicit none
interface
subroutine scatter(p, energy_bin_avg) bind(C)
import Particle, C_DOUBLE
type(Particle), intent(inout) :: p
real(C_DOUBLE), intent(in) :: energy_bin_avg(*)
end subroutine scatter
subroutine create_fission_sites(p, bank_array, size_bank, bank_array_size, &
material_xs) bind(C)
import Particle, Bank, C_INT64_T, MaterialMacroXS
type(Particle), intent(inout) :: p
type(Bank), intent(inout) :: bank_array(*)
integer(C_INT64_T), intent(inout) :: size_bank
integer(C_INT64_T), intent(in) :: bank_array_size
type(MaterialMacroXS), intent(in) :: material_xs
end subroutine create_fission_sites
subroutine absorption(p, material_xs) bind(C)
import Particle, MaterialMacroXS
type(Particle), intent(inout) :: p
type(MaterialMacroXS), intent(in) :: material_xs
end subroutine absorption
end interface
contains
!===============================================================================
! COLLISION_MG samples a nuclide and reaction and then calls the appropriate
! routine for that reaction
!===============================================================================
subroutine collision_mg(p)
type(Particle), intent(inout) :: p
! Add to collision counter for particle
p % n_collision = p % n_collision + 1
! Sample nuclide/reaction for the material the particle is in
call sample_reaction(p)
! Display information about collision
if (verbosity >= 10 .or. trace) then
call write_message(" " // "Energy Group = " // trim(to_str(p % g)))
end if
end subroutine collision_mg
!===============================================================================
! SAMPLE_REACTION samples a nuclide based on the macroscopic cross sections for
! each nuclide within a material and then samples a reaction for that nuclide
! and calls the appropriate routine to process the physics. Note that there is
! special logic when suvival biasing is turned on since fission and
! disappearance are treated implicitly.
!===============================================================================
subroutine sample_reaction(p)
type(Particle), intent(inout) :: p
type(Material), pointer :: mat
mat => materials(p % material)
! Create fission bank sites. Note that while a fission reaction is sampled,
! it never actually "happens", i.e. the weight of the particle does not
! change when sampling fission sites. The following block handles all
! absorption (including fission)
if (mat % fissionable) then
if (run_mode == MODE_EIGENVALUE) then
call create_fission_sites(p, fission_bank, n_bank, &
size(fission_bank, KIND=C_INT64_T), &
material_xs)
elseif (run_mode == MODE_FIXEDSOURCE .and. create_fission_neutrons) then
call create_fission_sites(p, p % secondary_bank, p % n_secondary, &
size(p % secondary_bank, KIND=C_INT64_T), &
material_xs)
end if
end if
! If survival biasing is being used, the following subroutine adjusts the
! weight of the particle. Otherwise, it checks to see if absorption occurs
if (material_xs % absorption > ZERO) then
call absorption(p, material_xs)
else
p % absorb_wgt = ZERO
end if
if (.not. p % alive) return
! Sample a scattering reaction and determine the secondary energy of the
! exiting neutron
call scatter(p, energy_bin_avg)
! Play russian roulette if survival biasing is turned on
if (survival_biasing) then
call russian_roulette(p)
if (.not. p % alive) return
end if
end subroutine sample_reaction
end module physics_mg

View file

@ -8,14 +8,76 @@
#include "openmc/constants.h"
#include "openmc/eigenvalue.h"
#include "openmc/error.h"
#include "openmc/material.h"
#include "openmc/math_functions.h"
#include "openmc/message_passing.h"
#include "openmc/mgxs_interface.h"
#include "openmc/physics_common.h"
#include "openmc/random_lcg.h"
#include "openmc/settings.h"
#include "openmc/simulation.h"
namespace openmc {
void
collision_mg(Particle* p, Bank* fission_bank, const int64_t fission_bank_size,
const double* energy_bin_avg, const MaterialMacroXS& material_xs)
{
// Add to the collision counter for the particle
p->n_collision++;
// Sample the reaction type
sample_reaction(p, fission_bank, fission_bank_size, energy_bin_avg,
material_xs);
// Display information about collision
if ((settings::verbosity >= 10) || (openmc_trace)) {
std::stringstream msg;
msg << " Energy Group = " << p->g;
write_message(msg, 1);
}
}
void
sample_reaction(Particle* p, Bank* fission_bank,
const int64_t fission_bank_size, const double* energy_bin_avg,
const MaterialMacroXS& material_xs)
{
// Create fission bank sites. Note that while a fission reaction is sampled,
// it never actually "happens", i.e. the weight of the particle does not
// change when sampling fission sites. The following block handles all
// absorption (including fission)
if (materials[p->material - 1]->fissionable) {
if (settings::run_mode == RUN_MODE_EIGENVALUE) {
create_fission_sites(p, fission_bank, n_bank, fission_bank_size,
material_xs);
} else if ((settings::run_mode == RUN_MODE_FIXEDSOURCE) &&
(settings::create_fission_neutrons)) {
create_fission_sites(p, p->secondary_bank, p->n_secondary,
MAX_SECONDARY, material_xs);
}
}
// If survival biasing is being used, the following subroutine adjusts the
// weight of the particle. Otherwise, it checks to see if absorption occurs.
if (material_xs.absorption > 0.) {
absorption(p, material_xs);
} else {
p->absorb_wgt = 0.;
}
if (!p->alive) return;
// Sample a scattering event to determine the energy of the exiting neutron
scatter(p, energy_bin_avg);
// Play Russian roulette if survival biasing is turned on
if (settings::survival_biasing) {
russian_roulette(p);
if (!p->alive) return;
}
}
void
scatter(Particle* p, const double* energy_bin_avg)
{
@ -42,7 +104,8 @@ scatter(Particle* p, const double* energy_bin_avg)
void
create_fission_sites(Particle* p, Bank* bank_array, int64_t& size_bank,
int64_t& bank_array_size, MaterialMacroXS& material_xs)
const int64_t& bank_array_size,
const MaterialMacroXS& material_xs)
{
// TODO: Heat generation from fission
@ -146,7 +209,7 @@ create_fission_sites(Particle* p, Bank* bank_array, int64_t& size_bank,
}
void
absorption(Particle* p, MaterialMacroXS& material_xs)
absorption(Particle* p, const MaterialMacroXS& material_xs)
{
if (settings::survival_biasing) {
// Determine weight absorbed in survival biasing

View file

@ -17,7 +17,7 @@ module tracking
use nuclide_header
use particle_header
use physics, only: collision
use physics_mg, only: collision_mg
! use physics_mg, only: collision_mg
use random_lcg, only: prn, prn_set_stream
use settings
use simulation_header
@ -33,6 +33,19 @@ module tracking
implicit none
interface
subroutine collision_mg(p, fission_bank, fission_bank_size, &
energy_bin_avg, material_xs) bind(C)
import Particle, Bank, C_INT64_T, C_DOUBLE, MaterialMacroXS
type(Particle), intent(inout) :: p
type(Bank), intent(inout) :: fission_bank(*)
integer(C_INT64_T), value, intent(in) :: fission_bank_size
real(C_DOUBLE), intent(in) :: energy_bin_avg(*)
type(MaterialMacroXS), intent(in) :: material_xs
end subroutine collision_mg
end interface
contains
!===============================================================================
@ -226,7 +239,8 @@ contains
if (run_CE) then
call collision(p)
else
call collision_mg(p)
call collision_mg(p, fission_bank, &
size(fission_bank, kind=C_INT64_T), energy_bin_avg, material_xs)
end if
! Score collision estimator tallies -- this is done after a collision