mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-27 05:35:49 -04:00
Remove Material p0 and element attributes on Fortran side
This commit is contained in:
Paul Romano
parent
08d55249ae
commit
f407982ff6
5 changed files with 25 additions and 122 deletions
|
|
@ -73,7 +73,6 @@ public:
|
|||
double volume_ {-1.0}; //!< Volume in [cm^3]
|
||||
bool fissionable_ {false}; //!< Does this material contain fissionable nuclides
|
||||
bool depletable_ {false}; //!< Is the material depletable?
|
||||
bool has_isotropic_nuclides_ {false};
|
||||
std::vector<bool> p0_; //!< Indicate which nuclides are to be treated with iso-in-lab scattering
|
||||
|
||||
// Thermal scattering tables
|
||||
|
|
|
|||
|
|
@ -388,22 +388,18 @@ contains
|
|||
subroutine read_materials_xml() bind(C)
|
||||
integer :: i ! loop index for materials
|
||||
integer :: j ! loop index for nuclides
|
||||
integer :: k ! loop index
|
||||
integer :: n ! number of nuclides
|
||||
integer :: index_nuclide ! index in nuclides
|
||||
integer :: index_element ! index in elements
|
||||
integer :: index_sab ! index in sab_tables
|
||||
logical :: file_exists ! does materials.xml exist?
|
||||
character(20) :: name ! name of nuclide, e.g. U235
|
||||
character(3) :: element ! name of element, e.g. Zr
|
||||
character(MAX_WORD_LEN) :: units ! units on density
|
||||
character(MAX_LINE_LEN) :: filename ! absolute path to materials.xml
|
||||
character(MAX_WORD_LEN), allocatable :: sarray(:)
|
||||
real(8) :: val ! value entered for density
|
||||
real(8) :: temp_dble ! temporary double prec. real
|
||||
logical :: sum_density ! density is sum of nuclide densities
|
||||
type(VectorChar) :: names ! temporary list of nuclide names
|
||||
type(VectorChar) :: list_iso_lab ! temporary list of isotropic lab scatterers
|
||||
type(VectorReal) :: densities ! temporary list of nuclide densities
|
||||
type(Material), pointer :: mat => null()
|
||||
type(XMLDocument) :: doc
|
||||
|
|
@ -624,19 +620,6 @@ contains
|
|||
end do INDIVIDUAL_NUCLIDES
|
||||
end if
|
||||
|
||||
! =======================================================================
|
||||
! READ AND PARSE <isotropic> element
|
||||
|
||||
if (check_for_node(node_mat, "isotropic")) then
|
||||
n = node_word_count(node_mat, "isotropic")
|
||||
allocate(sarray(n))
|
||||
call get_node_array(node_mat, "isotropic", sarray)
|
||||
do j = 1, n
|
||||
call list_iso_lab % push_back(sarray(j))
|
||||
end do
|
||||
deallocate(sarray)
|
||||
end if
|
||||
|
||||
! ========================================================================
|
||||
! COPY NUCLIDES TO ARRAYS IN MATERIAL
|
||||
|
||||
|
|
@ -645,7 +628,6 @@ contains
|
|||
mat % n_nuclides = n
|
||||
allocate(mat % names(n))
|
||||
allocate(mat % nuclide(n))
|
||||
allocate(mat % element(n))
|
||||
allocate(mat % atom_density(n))
|
||||
|
||||
ALL_NUCLIDES: do j = 1, mat % n_nuclides
|
||||
|
|
@ -667,51 +649,12 @@ contains
|
|||
mat % nuclide(j) = nuclide_dict % get(name)
|
||||
end if
|
||||
|
||||
! If the corresponding element hasn't been encountered yet and photon
|
||||
! transport will be used, we need to add its symbol to the element_dict
|
||||
if (photon_transport) then
|
||||
element = name(1:scan(name, '0123456789') - 1)
|
||||
|
||||
! Make sure photon cross section data is available
|
||||
if (.not. library_present(LIBRARY_PHOTON, element)) then
|
||||
call fatal_error("Could not find element " // trim(element) &
|
||||
// " in cross_sections data file!")
|
||||
end if
|
||||
|
||||
if (.not. element_dict % has(element)) then
|
||||
index_element = index_element + 1
|
||||
mat % element(j) = index_element
|
||||
|
||||
call element_dict % set(element, index_element)
|
||||
else
|
||||
mat % element(j) = element_dict % get(element)
|
||||
end if
|
||||
end if
|
||||
|
||||
! Copy name and atom/weight percent
|
||||
mat % names(j) = name
|
||||
mat % atom_density(j) = densities % data(j)
|
||||
|
||||
end do ALL_NUCLIDES
|
||||
|
||||
if (run_CE) then
|
||||
! By default, isotropic-in-lab is not used
|
||||
if (list_iso_lab % size() > 0) then
|
||||
mat % has_isotropic_nuclides = .true.
|
||||
allocate(mat % p0(n))
|
||||
mat % p0(:) = .false.
|
||||
|
||||
! Apply isotropic-in-lab treatment to specified nuclides
|
||||
do j = 1, list_iso_lab % size()
|
||||
do k = 1, n
|
||||
if (names % data(k) == list_iso_lab % data(j)) then
|
||||
mat % p0(k) = .true.
|
||||
end if
|
||||
end do
|
||||
end do
|
||||
end if
|
||||
end if
|
||||
|
||||
! Check to make sure either all atom percents or all weight percents are
|
||||
! given
|
||||
if (.not. (all(mat % atom_density >= ZERO) .or. &
|
||||
|
|
@ -726,7 +669,6 @@ contains
|
|||
! Clear lists
|
||||
call names % clear()
|
||||
call densities % clear()
|
||||
call list_iso_lab % clear()
|
||||
|
||||
! Add material to dictionary
|
||||
call material_dict % set(mat % id(), i)
|
||||
|
|
|
|||
|
|
@ -252,7 +252,6 @@ Material::Material(pugi::xml_node node)
|
|||
if (settings::run_CE) {
|
||||
// By default, isotropic-in-lab is not used
|
||||
if (iso_lab.size() > 0) {
|
||||
has_isotropic_nuclides_ = true;
|
||||
p0_.resize(n);
|
||||
|
||||
// Apply isotropic-in-lab treatment to specified nuclides
|
||||
|
|
@ -458,14 +457,8 @@ void Material::init_bremsstrahlung()
|
|||
auto n_k = data::ttb_k_grid.size();
|
||||
auto n_e = data::ttb_e_grid.size();
|
||||
|
||||
// Get pointers to nuclides, elements, densities
|
||||
int32_t index;
|
||||
openmc_get_material_index(id_, &index);
|
||||
int* nuclide_;
|
||||
double* atom_density_;
|
||||
int n;
|
||||
openmc_material_get_densities(index, &nuclide_, &atom_density_, &n);
|
||||
int* element_ = material_element(index);
|
||||
// Determine number of elements
|
||||
int n = element_.size();
|
||||
|
||||
for (int particle = 0; particle < 2; ++particle) {
|
||||
// Loop over logic twice, once for electron, once for positron
|
||||
|
|
@ -497,10 +490,8 @@ void Material::init_bremsstrahlung()
|
|||
// fixed in the future.
|
||||
for (int i = 0; i < n; ++i) {
|
||||
// Get pointer to current element
|
||||
// TODO: off-by-one
|
||||
const auto& elm = data::elements[element_[i] - 1];
|
||||
// TODO: off-by-one
|
||||
double awr = data::nuclides[nuclide_[i] - 1]->awr_;
|
||||
const auto& elm = data::elements[element_[i]];
|
||||
double awr = data::nuclides[nuclide_[i]]->awr_;
|
||||
|
||||
// Get atomic density and mass density of nuclide given atom/weight percent
|
||||
double atom_density = (atom_density_[0] > 0.0) ?
|
||||
|
|
|
|||
|
|
@ -70,7 +70,6 @@ module material_header
|
|||
character(len=104) :: name = "" ! User-defined name
|
||||
integer :: n_nuclides = 0 ! number of nuclides
|
||||
integer, allocatable :: nuclide(:) ! index in nuclides array
|
||||
integer, allocatable :: element(:) ! index in elements array
|
||||
real(8) :: density ! total atom density in atom/b-cm
|
||||
real(C_DOUBLE), allocatable :: atom_density(:) ! nuclide atom density in atom/b-cm
|
||||
real(8) :: density_gpcc ! total density in g/cm^3
|
||||
|
|
@ -93,11 +92,6 @@ module material_header
|
|||
|
||||
! Does this material contain fissionable nuclides? Is it depletable?
|
||||
logical :: depletable = .false.
|
||||
|
||||
! enforce isotropic scattering in lab for specific nuclides
|
||||
logical :: has_isotropic_nuclides = .false.
|
||||
logical, allocatable :: p0(:)
|
||||
|
||||
contains
|
||||
procedure :: id => material_id
|
||||
procedure :: set_id => material_set_id
|
||||
|
|
@ -233,27 +227,4 @@ contains
|
|||
err = 0
|
||||
end function openmc_extend_materials
|
||||
|
||||
!===============================================================================
|
||||
! Fortran compatibility
|
||||
!===============================================================================
|
||||
|
||||
function material_isotropic(i_material, i_nuc_mat) result(iso) bind(C)
|
||||
integer(C_INT), value :: i_material
|
||||
integer(C_INT), value :: i_nuc_mat
|
||||
logical(C_BOOL) :: iso
|
||||
|
||||
iso = .false.
|
||||
associate (mat => materials(i_material))
|
||||
if (mat % has_isotropic_nuclides) then
|
||||
iso = mat % p0(i_nuc_mat)
|
||||
end if
|
||||
end associate
|
||||
end function
|
||||
|
||||
function material_element(i_material) result(ptr) bind(C)
|
||||
integer(C_INT), value :: i_material
|
||||
type(C_PTR) :: ptr
|
||||
ptr = C_LOC(materials(i_material) % element(1))
|
||||
end function
|
||||
|
||||
end module material_header
|
||||
|
|
|
|||
|
|
@ -483,11 +483,8 @@ int sample_element(Particle* p)
|
|||
double cutoff = prn() * simulation::material_xs.total;
|
||||
|
||||
// Get pointers to elements, densities
|
||||
int* nuclide;
|
||||
double* density;
|
||||
int n;
|
||||
openmc_material_get_densities(p->material, &nuclide, &density, &n);
|
||||
int* element = material_element(p->material);
|
||||
const auto& mat {model::materials[p->material - 1]};
|
||||
int n = mat->nuclide_.size();
|
||||
|
||||
int i = 0;
|
||||
double prob = 0.0;
|
||||
|
|
@ -500,9 +497,8 @@ int sample_element(Particle* p)
|
|||
}
|
||||
|
||||
// Find atom density
|
||||
// TODO: off-by-one
|
||||
i_element = element[i] - 1;
|
||||
double atom_density = density[i];
|
||||
i_element = mat->element_[i];
|
||||
double atom_density = mat->atom_density_[i];
|
||||
|
||||
// Determine microscopic cross section
|
||||
double sigma = atom_density * simulation::micro_photon_xs[i_element].total;
|
||||
|
|
@ -710,21 +706,25 @@ void scatter(Particle* p, int i_nuclide, int i_nuc_mat)
|
|||
p->event = EVENT_SCATTER;
|
||||
|
||||
// Sample new outgoing angle for isotropic-in-lab scattering
|
||||
if (material_isotropic(p->material, i_nuc_mat)) {
|
||||
// Sample isotropic-in-lab outgoing direction
|
||||
double mu = 2.0*prn() - 1.0;
|
||||
double phi = 2.0*PI*prn();
|
||||
Direction u_new;
|
||||
u_new.x = mu;
|
||||
u_new.y = std::sqrt(1.0 - mu*mu)*std::cos(phi);
|
||||
u_new.z = std::sqrt(1.0 - mu*mu)*std::sin(phi);
|
||||
// TODO: off-by-one
|
||||
const auto& mat {model::materials[p->material - 1]};
|
||||
if (!mat->p0_.empty()) {
|
||||
if (mat->p0_[i_nuc_mat]) {
|
||||
// Sample isotropic-in-lab outgoing direction
|
||||
double mu = 2.0*prn() - 1.0;
|
||||
double phi = 2.0*PI*prn();
|
||||
Direction u_new;
|
||||
u_new.x = mu;
|
||||
u_new.y = std::sqrt(1.0 - mu*mu)*std::cos(phi);
|
||||
u_new.z = std::sqrt(1.0 - mu*mu)*std::sin(phi);
|
||||
|
||||
p->mu = u_old.dot(u_new);
|
||||
p->mu = u_old.dot(u_new);
|
||||
|
||||
// Change direction of particle
|
||||
p->coord[0].uvw[0] = u_new.x;
|
||||
p->coord[0].uvw[1] = u_new.y;
|
||||
p->coord[0].uvw[2] = u_new.z;
|
||||
// Change direction of particle
|
||||
p->coord[0].uvw[0] = u_new.x;
|
||||
p->coord[0].uvw[1] = u_new.y;
|
||||
p->coord[0].uvw[2] = u_new.z;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue