Implementing minor edits from review on 1/11 by @paulromano

This commit is contained in:
Adam Nelson 2016-02-12 04:58:11 -05:00
parent 307617538c
commit 9f02388647
9 changed files with 130 additions and 160 deletions

View file

@ -126,7 +126,7 @@ contains
if (run_CE) then
call same_nuclide_list()
else
call same_NuclideMG_list()
call same_nuclidemg_list()
end if
! Construct information needed for nuclear data

View file

@ -216,8 +216,8 @@ contains
end if
! Make sure that either eigenvalue or fixed source was specified
if (.not.check_for_node(doc, "eigenvalue") .and. &
.not.check_for_node(doc, "fixed_source")) then
if (.not. check_for_node(doc, "eigenvalue") .and. &
.not. check_for_node(doc, "fixed_source")) then
call fatal_error("<eigenvalue> or <fixed_source> not specified.")
end if
@ -230,7 +230,7 @@ contains
call get_node_ptr(doc, "eigenvalue", node_mode)
! Check number of particles
if (.not.check_for_node(node_mode, "particles")) then
if (.not. check_for_node(node_mode, "particles")) then
call fatal_error("Need to specify number of particles per generation.")
end if
@ -300,7 +300,7 @@ contains
call get_node_ptr(doc, "fixed_source", node_mode)
! Check number of particles
if (.not.check_for_node(node_mode, "particles")) then
if (.not. check_for_node(node_mode, "particles")) then
call fatal_error("Need to specify number of particles per batch.")
end if
@ -2050,9 +2050,9 @@ contains
! READ AND PARSE <nuclide> TAGS
! Check to ensure material has at least one nuclide
if ((.not. check_for_node(node_mat, "nuclide") .and. &
.not. check_for_node(node_mat, "element")) .and. &
(.not. check_for_node(node_mat, "macroscopic"))) then
if (.not. check_for_node(node_mat, "nuclide") .and. &
.not. check_for_node(node_mat, "element") .and. &
.not. check_for_node(node_mat, "macroscopic")) then
call fatal_error("No macroscopic data, nuclides or natural elements &
&specified on material " // trim(to_str(mat % id)))
end if
@ -2061,7 +2061,7 @@ contains
! them as nuclides. This is all really a facade so the user thinks they
! are entering in macroscopic data but the code treats them the same
! as nuclides internally.
! Get pointer list of XML <nuclide>
! Get pointer list of XML <macroscopic>
call get_node_list(node_mat, "macroscopic", node_macro_list)
if (get_list_size(node_macro_list) > 1) then
call fatal_error("Only one macroscopic object permitted per material, " &
@ -2077,7 +2077,7 @@ contains
end if
! Check for cross section
if (.not.check_for_node(node_nuc, "xs")) then
if (.not. check_for_node(node_nuc, "xs")) then
if (default_xs == '') then
call fatal_error("No cross section specified for macroscopic data &
& in material " // trim(to_str(mat % id)))
@ -2130,13 +2130,13 @@ contains
call get_list_item(node_nuc_list, j, node_nuc)
! Check for empty name on nuclide
if (.not.check_for_node(node_nuc, "name")) then
if (.not. check_for_node(node_nuc, "name")) then
call fatal_error("No name specified on nuclide in material " &
// trim(to_str(mat % id)))
end if
! Check for cross section
if (.not.check_for_node(node_nuc, "xs")) then
if (.not. check_for_node(node_nuc, "xs")) then
if (default_xs == '') then
call fatal_error("No cross section specified for nuclide in &
&material " // trim(to_str(mat % id)))
@ -2160,8 +2160,8 @@ contains
if (units == 'macro') then
call list_density % append(ONE)
else
if (.not.check_for_node(node_nuc, "ao") .and. &
.not.check_for_node(node_nuc, "wo")) then
if (.not. check_for_node(node_nuc, "ao") .and. &
.not. check_for_node(node_nuc, "wo")) then
call fatal_error("No atom or weight percent specified for nuclide " &
// trim(name))
elseif (check_for_node(node_nuc, "ao") .and. &
@ -2192,7 +2192,7 @@ contains
call get_list_item(node_ele_list, j, node_ele)
! Check for empty name on natural element
if (.not.check_for_node(node_ele, "name")) then
if (.not. check_for_node(node_ele, "name")) then
call fatal_error("No name specified on nuclide in material " &
// trim(to_str(mat % id)))
end if
@ -2212,8 +2212,8 @@ contains
! Check if no atom/weight percents were specified or if both atom and
! weight percents were specified
if (.not.check_for_node(node_ele, "ao") .and. &
.not.check_for_node(node_ele, "wo")) then
if (.not. check_for_node(node_ele, "ao") .and. &
.not. check_for_node(node_ele, "wo")) then
call fatal_error("No atom or weight percent specified for element " &
// trim(name))
elseif (check_for_node(node_ele, "ao") .and. &
@ -2357,8 +2357,8 @@ contains
call get_list_item(node_sab_list, j, node_sab)
! Determine name of S(a,b) table
if (.not.check_for_node(node_sab, "name") .or. &
.not.check_for_node(node_sab, "xs")) then
if (.not. check_for_node(node_sab, "name") .or. &
.not. check_for_node(node_sab, "xs")) then
call fatal_error("Need to specify <name> and <xs> for S(a,b) &
&table.")
end if
@ -2588,8 +2588,8 @@ contains
end if
! Make sure either upper-right or width was specified
if (.not.check_for_node(node_mesh, "upper_right") .and. &
.not.check_for_node(node_mesh, "width")) then
if (.not. check_for_node(node_mesh, "upper_right") .and. &
.not. check_for_node(node_mesh, "width")) then
call fatal_error("Must specify either <upper_right> and <width> on a &
&tally mesh.")
end if
@ -3377,30 +3377,12 @@ contains
case ('n2n', '(n,2n)')
t % score_bins(j) = N_2N
! Disallow for MG mode since data not present
if (.not. run_CE) then
call fatal_error("Cannot tally (n,2n) reaction rate in &
&multi-group mode")
end if
case ('n3n', '(n,3n)')
t % score_bins(j) = N_3N
! Disallow for MG mode since data not present
if (.not. run_CE) then
call fatal_error("Cannot tally (n,3n) reaction rate in &
&multi-group mode")
end if
case ('n4n', '(n,4n)')
t % score_bins(j) = N_4N
! Disallow for MG mode since data not present
if (.not. run_CE) then
call fatal_error("Cannot tally (n,4n) reaction rate in &
&multi-group mode")
end if
case ('absorption')
t % score_bins(j) = SCORE_ABSORPTION
if (t % find_filter(FILTER_ENERGYOUT) > 0) then
@ -3590,10 +3572,10 @@ contains
! Do a check at the end (instead of for every case) to make sure
! the tallies are compatible with MG mode where we have less detailed
! nuclear data
if (.not. run_CE .and. t % score_bins(j) > 0) then
call fatal_error("Cannot tally " // trim(score_name) // &
" reaction rate in multi-group mode")
end if
if (.not. run_CE .and. t % score_bins(j) > 0) then
call fatal_error("Cannot tally " // trim(score_name) // &
" reaction rate in multi-group mode")
end if
end do
t % n_score_bins = n_scores
@ -4455,7 +4437,7 @@ contains
end if
! determine metastable state
if (.not.check_for_node(node_ace, "metastable")) then
if (.not. check_for_node(node_ace, "metastable")) then
listing % metastable = .false.
else
listing % metastable = .true.

View file

@ -27,10 +27,7 @@ module macroxs_header
subroutine macroxs_init_(this, mat, nuclides, groups, get_kfiss, get_fiss, &
max_order, scatt_type, legendre_mu_points, &
error_code, error_text)
import MacroXS
import Material
import NuclideMGContainer
import MAX_LINE_LEN
import MacroXS, Material, NuclideMGContainer, MAX_LINE_LEN
class(MacroXS), intent(inout) :: this ! The MacroXS to initialize
type(Material), pointer, intent(in) :: mat ! base material
type(NuclideMGContainer), intent(in) :: nuclides(:) ! List of nuclides to harvest from
@ -383,14 +380,14 @@ contains
select type(nuc => nuclides(mat % nuclide(i)) % obj)
type is (NuclideAngle)
if (npol == -1) then
npol = nuc % Npol
nazi = nuc % Nazi
npol = nuc % n_pol
nazi = nuc % n_azi
allocate(this % polar(npol))
this % polar = nuc % polar
allocate(this % azimuthal(nazi))
this % azimuthal = nuc % azimuthal
else
if ((npol /= nuc % Npol) .or. (nazi /= nuc % Nazi)) then
if ((npol /= nuc % n_pol) .or. (nazi /= nuc % n_azi)) then
error_code = 1
error_text = "All Angular Data Must Be Same Length!"
end if

View file

@ -560,6 +560,7 @@ contains
!===============================================================================
! EXPAND_HARMONIC expands a given series of real spherical harmonics
!===============================================================================
pure function expand_harmonic(data, order, uvw) result(val)
real(8), intent(in) :: data(:)
integer, intent(in) :: order
@ -584,6 +585,7 @@ contains
! EVALUATE_LEGENDRE Find the value of f(x) given a set of Legendre coefficients
! and the value of x
!===============================================================================
pure function evaluate_legendre(data, x) result(val)
real(8), intent(in) :: data(:)
real(8), intent(in) :: x
@ -591,9 +593,9 @@ contains
integer :: l
val = 0.5_8 * data(1)
val = HALF * data(1)
do l = 1, size(data) - 1
val = val + (real(l,8) + 0.5_8) * data(l + 1) * calc_pn(l,x)
val = val + (real(l,8) + HALF) * data(l + 1) * calc_pn(l,x)
end do
end function evaluate_legendre

View file

@ -125,7 +125,7 @@ contains
end select
! Now read in the data specific to the type we just declared
call NuclideMG_init(nuclides_MG(i_nuclide) % obj, node_xsdata, &
call nuclidemg_init(nuclides_MG(i_nuclide) % obj, node_xsdata, &
energy_groups, get_kfiss, get_fiss, error_code, &
error_text)
@ -169,13 +169,13 @@ contains
end subroutine read_mgxs
!===============================================================================
! SAME_NUCLIDE_LIST creates a linked list for each nuclide containing the
! SAME_NUCLIDEMG_LIST creates a linked list for each nuclide containing the
! indices in the nuclides array of all other instances of that nuclide. For
! example, the same nuclide may exist at multiple temperatures resulting
! in multiple entries in the nuclides array for a single zaid number.
!===============================================================================
subroutine same_NuclideMG_list()
subroutine same_nuclidemg_list()
integer :: i ! index in nuclides array
integer :: j ! index in nuclides array
@ -188,21 +188,21 @@ contains
end do
end do
end subroutine same_NuclideMG_list
end subroutine same_nuclidemg_list
!===============================================================================
! NUCLIDE_*_INIT reads in the data from the XML file, as already accessed
!===============================================================================
subroutine NuclideMG_init(this, node_xsdata, groups, get_kfiss, get_fiss, &
subroutine nuclidemg_init(this, node_xsdata, groups, get_kfiss, get_fiss, &
error_code, error_text)
class(NuclideMG), intent(inout) :: this ! Working Object
type(Node), pointer, intent(in) :: node_xsdata ! Data from data.xml
integer, intent(in) :: groups ! Number of Energy groups
logical, intent(in) :: get_kfiss ! Need Kappa-Fission?
logical, intent(in) :: get_fiss ! Should we get fiss data?
integer, intent(inout) :: error_code ! Code signifying error
character(MAX_LINE_LEN), intent(inout) :: error_text ! Error message to print
type(Node), pointer, intent(in) :: node_xsdata ! Data from data.xml
integer, intent(in) :: groups ! Number of Energy groups
logical, intent(in) :: get_kfiss ! Need Kappa-Fission?
logical, intent(in) :: get_fiss ! Should we get fiss data?
integer, intent(inout) :: error_code ! Code signifying error
character(MAX_LINE_LEN), intent(inout) :: error_text ! Error to print
type(Node), pointer :: node_legendre_mu
character(MAX_LINE_LEN) :: temp_str
@ -298,16 +298,16 @@ contains
select type(this)
type is (NuclideIso)
call NuclideIso_init(this, node_xsdata, groups, get_kfiss, get_fiss, &
call nuclideiso_init(this, node_xsdata, groups, get_kfiss, get_fiss, &
error_code, error_text)
type is (NuclideAngle)
call NuclideAngle_init(this, node_xsdata, groups, get_kfiss, get_fiss, &
call nuclideangle_init(this, node_xsdata, groups, get_kfiss, get_fiss, &
error_code, error_text)
end select
end subroutine NuclideMG_init
end subroutine nuclidemg_init
subroutine NuclideIso_init(this, node_xsdata, groups, get_kfiss, get_fiss, &
subroutine nuclideiso_init(this, node_xsdata, groups, get_kfiss, get_fiss, &
error_code, error_text)
class(NuclideIso), intent(inout) :: this ! Working Object
type(Node), pointer, intent(in) :: node_xsdata ! Data from data.xml
@ -436,10 +436,10 @@ contains
this % mult = ONE
end if
end subroutine NuclideIso_init
end subroutine nuclideiso_init
subroutine NuclideAngle_init(this, node_xsdata, groups, get_kfiss, get_fiss, &
error_code, error_text)
subroutine nuclideangle_init(this, node_xsdata, groups, get_kfiss, get_fiss, &
error_code, error_text)
class(NuclideAngle), intent(inout) :: this ! Working Object
type(Node), pointer, intent(in) :: node_xsdata ! Data from data.xml
integer, intent(in) :: groups ! Number of Energy groups
@ -463,7 +463,7 @@ contains
end if
if (check_for_node(node_xsdata, "num_polar")) then
call get_node_value(node_xsdata, "num_polar", this % Npol)
call get_node_value(node_xsdata, "num_polar", this % n_pol)
else
error_code = 1
error_text = "num_polar Must Be Provided!"
@ -471,7 +471,7 @@ contains
end if
if (check_for_node(node_xsdata, "num_azimuthal")) then
call get_node_value(node_xsdata, "num_azimuthal", this % Nazi)
call get_node_value(node_xsdata, "num_azimuthal", this % n_azi)
else
error_code = 1
error_text = "num_azimuthal Must Be Provided!"
@ -479,17 +479,17 @@ contains
end if
! Load angle data, if present (else equally spaced)
allocate(this % polar(this % Npol))
allocate(this % azimuthal(this % Nazi))
allocate(this % polar(this % n_pol))
allocate(this % azimuthal(this % n_azi))
if (check_for_node(node_xsdata, "polar")) then
error_code = 1
error_text = "User-Specified polar angle bins not yet supported!"
return
call get_node_array(node_xsdata, "polar", this % polar)
else
dangle = PI / real(this % Npol,8)
do iangle = 1, this % Npol
this % polar(iangle) = (real(iangle,8) - 0.5_8) * dangle
dangle = PI / real(this % n_pol,8)
do iangle = 1, this % n_pol
this % polar(iangle) = (real(iangle,8) - HALF) * dangle
end do
end if
if (check_for_node(node_xsdata, "azimuthal")) then
@ -498,9 +498,9 @@ contains
return
call get_node_array(node_xsdata, "azimuthal", this % azimuthal)
else
dangle = TWO * PI / real(this % Nazi,8)
do iangle = 1, this % Nazi
this % azimuthal(iangle) = -PI + (real(iangle,8) - 0.5_8) * dangle
dangle = TWO * PI / real(this % n_azi,8)
do iangle = 1, this % n_azi
this % azimuthal(iangle) = -PI + (real(iangle,8) - HALF) * dangle
end do
end if
@ -509,19 +509,19 @@ contains
if (check_for_node(node_xsdata, "chi")) then
! Get chi
allocate(temp_arr(groups * this % Nazi * this % Npol))
allocate(temp_arr(groups * this % n_azi * this % n_pol))
call get_node_array(node_xsdata, "chi", temp_arr)
allocate(this % chi(groups, this % Nazi, this % Npol))
this % chi = reshape(temp_arr, (/groups, this % Nazi, this % Npol/))
allocate(this % chi(groups, this % n_azi, this % n_pol))
this % chi = reshape(temp_arr, (/groups, this % n_azi, this % n_pol/))
deallocate(temp_arr)
! Get nu_fission (as a vector)
if (check_for_node(node_xsdata, "nu_fission")) then
allocate(temp_arr(groups * this % Nazi * this % Npol))
allocate(temp_arr(groups * this % n_azi * this % n_pol))
call get_node_array(node_xsdata, "nu_fission", temp_arr)
allocate(this % nu_fission(groups, 1, this % Nazi, this % Npol))
this % nu_fission = reshape(temp_arr, (/groups, 1, this % Nazi, &
this % Npol/))
allocate(this % nu_fission(groups, 1, this % n_azi, this % n_pol))
this % nu_fission = reshape(temp_arr, (/groups, 1, this % n_azi, &
this % n_pol/))
deallocate(temp_arr)
else
error_code = 1
@ -533,11 +533,11 @@ contains
! Get nu_fission (as a matrix)
if (check_for_node(node_xsdata, "nu_fission")) then
allocate(temp_arr(groups * this % Nazi * this % Npol))
allocate(temp_arr(groups * this % n_azi * this % n_pol))
call get_node_array(node_xsdata, "nu_fission", temp_arr)
allocate(this % nu_fission(groups, groups, this % Nazi, this % Npol))
allocate(this % nu_fission(groups, groups, this % n_azi, this % n_pol))
this % nu_fission = reshape(temp_arr, (/groups, groups, &
this % Nazi, this % Npol/))
this % n_azi, this % n_pol/))
deallocate(temp_arr)
else
error_code = 1
@ -547,10 +547,10 @@ contains
end if
if (get_fiss) then
if (check_for_node(node_xsdata, "fission")) then
allocate(temp_arr(groups * this % Nazi * this % Npol))
allocate(temp_arr(groups * this % n_azi * this % n_pol))
call get_node_array(node_xsdata, "fission", temp_arr)
allocate(this % fission(groups, this % Nazi, this % Npol))
this % fission = reshape(temp_arr, (/groups, this % Nazi, this % Npol/))
allocate(this % fission(groups, this % n_azi, this % n_pol))
this % fission = reshape(temp_arr, (/groups, this % n_azi, this % n_pol/))
deallocate(temp_arr)
else
error_code = 1
@ -561,10 +561,10 @@ contains
end if
if (get_kfiss) then
if (check_for_node(node_xsdata, "kappa_fission")) then
allocate(temp_arr(groups * this % Nazi * this % Npol))
allocate(temp_arr(groups * this % n_azi * this % n_pol))
call get_node_array(node_xsdata, "kappa_fission", temp_arr)
allocate(this % k_fission(groups, this % Nazi, this % Npol))
this % k_fission = reshape(temp_arr, (/groups, this % Nazi, this % Npol/))
allocate(this % k_fission(groups, this % n_azi, this % n_pol))
this % k_fission = reshape(temp_arr, (/groups, this % n_azi, this % n_pol/))
deallocate(temp_arr)
else
error_code = 1
@ -576,10 +576,10 @@ contains
end if
if (check_for_node(node_xsdata, "absorption")) then
allocate(temp_arr(groups * this % Nazi * this % Npol))
allocate(temp_arr(groups * this % n_azi * this % n_pol))
call get_node_array(node_xsdata, "absorption", temp_arr)
allocate(this % absorption(groups, this % Nazi, this % Npol))
this % absorption = reshape(temp_arr, (/groups, this % Nazi, this % Npol/))
allocate(this % absorption(groups, this % n_azi, this % n_pol))
this % absorption = reshape(temp_arr, (/groups, this % n_azi, this % n_pol/))
deallocate(temp_arr)
else
error_code = 1
@ -587,12 +587,12 @@ contains
return
end if
allocate(this % scatter(groups, groups, order_dim, this % Nazi, this % Npol))
allocate(this % scatter(groups, groups, order_dim, this % n_azi, this % n_pol))
if (check_for_node(node_xsdata, "scatter")) then
allocate(temp_arr(groups * groups * order_dim * this % Nazi * this%Npol))
allocate(temp_arr(groups * groups * order_dim * this % n_azi * this%n_pol))
call get_node_array(node_xsdata, "scatter", temp_arr)
this % scatter = reshape(temp_arr, (/groups, groups, order_dim, &
this%Nazi,this%Npol/))
this%n_azi,this%n_pol/))
deallocate(temp_arr)
else
error_code = 1
@ -601,23 +601,23 @@ contains
end if
if (check_for_node(node_xsdata, "total")) then
allocate(temp_arr(groups * this % Nazi * this % Npol))
allocate(temp_arr(groups * this % n_azi * this % n_pol))
call get_node_array(node_xsdata, "total", temp_arr)
allocate(this % total(groups, this % Nazi, this % Npol))
this % total = reshape(temp_arr, (/groups, this % Nazi, this % Npol/))
allocate(this % total(groups, this % n_azi, this % n_pol))
this % total = reshape(temp_arr, (/groups, this % n_azi, this % n_pol/))
deallocate(temp_arr)
else
this % total = this % absorption + sum(this%scatter(:,:,1,:,:),dim=1)
end if
! Get Mult Data
allocate(this % mult(groups, groups, this % Nazi, this % Npol))
allocate(this % mult(groups, groups, this % n_azi, this % n_pol))
if (check_for_node(node_xsdata, "multiplicity")) then
arr_len = get_arraysize_double(node_xsdata, "multiplicity")
if (arr_len == groups * groups * this % Nazi * this % Npol) then
if (arr_len == groups * groups * this % n_azi * this % n_pol) then
allocate(temp_arr(arr_len))
call get_node_array(node_xsdata, "multiplicity", temp_arr)
this % mult = reshape(temp_arr, (/groups, groups, this % Nazi, this % Npol/))
this % mult = reshape(temp_arr, (/groups, groups, this % n_azi, this % n_pol/))
deallocate(temp_arr)
else
error_code = 1
@ -628,7 +628,7 @@ contains
this % mult = ONE
end if
end subroutine NuclideAngle_init
end subroutine nuclideangle_init
!===============================================================================
@ -643,7 +643,6 @@ contains
logical :: get_kfiss, get_fiss
integer :: error_code
character(MAX_LINE_LEN) :: error_text
integer :: representation
integer :: scatt_type
integer :: legendre_mu_points
@ -675,19 +674,11 @@ contains
! At the same time, we will find the scattering type, as that will dictate
! how we allocate the scatter object within macroxs
legendre_mu_points = nuclides_MG(mat % nuclide(1)) % obj % legendre_mu_points
scatt_type = nuclides_MG(mat % nuclide(1)) % obj % scatt_type
select type(nuc => nuclides_MG(mat % nuclide(1)) % obj)
type is (NuclideIso)
representation = MGXS_ISOTROPIC
type is (NuclideAngle)
representation = MGXS_ANGLE
end select
scatt_type = nuclides_MG(mat % nuclide(1)) % obj % scatt_type
! Now allocate accordingly
select case(representation)
case(MGXS_ISOTROPIC)
allocate(MacroXSIso :: macro_xs(i_mat) % obj)
case(MGXS_ANGLE)
type is (NuclideAngle)
allocate(MacroXSAngle :: macro_xs(i_mat) % obj)
end select
@ -696,9 +687,7 @@ contains
scatt_type, legendre_mu_points, &
error_code, error_text)
! Handle any errors
if (error_code /= 0) then
call fatal_error(trim(error_text))
end if
if (error_code /= 0) call fatal_error(trim(error_text))
end do
end subroutine create_macro_xs

View file

@ -177,7 +177,7 @@ module nuclide_header
type, extends(NuclideMG) :: NuclideAngle
! Microscopic cross sections. Dimensions are: (Npol, Nazi, Nl, Ng, Ng)
! Microscopic cross sections. Dimensions are: (n_pol, n_azi, Nl, Ng, Ng)
real(8), allocatable :: total(:,:,:) ! total cross section
real(8), allocatable :: absorption(:,:,:) ! absorption cross section
real(8), allocatable :: scatter(:,:,:,:,:) ! scattering information
@ -188,8 +188,8 @@ module nuclide_header
real(8), allocatable :: mult(:,:,:,:) ! Scatter multiplicity (Gout x Gin)
! In all cases, right-most indices are theta, phi
integer :: Npol ! Number of polar angles
integer :: Nazi ! Number of azimuthal angles
integer :: n_pol ! Number of polar angles
integer :: n_azi ! Number of azimuthal angles
real(8), allocatable :: polar(:) ! polar angles
real(8), allocatable :: azimuthal(:) ! azimuthal angles
@ -476,8 +476,8 @@ module nuclide_header
! Write Basic Nuclide Information
call nuclidemg_print(this, unit_)
write(unit_,*) ' # of Polar Angles = ' // trim(to_str(this % Npol))
write(unit_,*) ' # of Azimuthal Angles = ' // trim(to_str(this % Nazi))
write(unit_,*) ' # of Polar Angles = ' // trim(to_str(this % n_pol))
write(unit_,*) ' # of Azimuthal Angles = ' // trim(to_str(this % n_azi))
! Determine size of mgxs and scattering matrices
size_scattmat = (size(this % scatter) + size(this % mult)) * 8

View file

@ -96,8 +96,8 @@ module particle_header
contains
procedure :: initialize => initialize_particle
procedure :: clear => clear_particle
procedure :: initialize_from_source => initialize_from_source
procedure :: create_secondary => create_secondary
procedure :: initialize_from_source
procedure :: create_secondary
end type Particle
contains

View file

@ -79,7 +79,7 @@ contains
! SCATTDATA_INIT builds the scattdata object
!===============================================================================
subroutine scattdatabase_init(this, order, energy, mult)
subroutine scattdata_init(this, order, energy, mult)
class(ScattData), intent(inout) :: this ! Object to work on
integer, intent(in) :: order ! Data Order
real(8), intent(in) :: energy(:,:) ! Energy Transfer Matrix
@ -96,7 +96,7 @@ contains
allocate(this % data(order, groups, groups))
this % data = ZERO
end subroutine scattdatabase_init
end subroutine scattdata_init
subroutine scattdatalegendre_init(this, order, energy, mult, coeffs)
class(ScattDataLegendre), intent(inout) :: this ! Object to work on
@ -105,7 +105,7 @@ contains
real(8), intent(in) :: mult(:,:) ! Scatter Prod'n Matrix
real(8), intent(in) :: coeffs(:,:,:) ! Coefficients to use
call scattdatabase_init(this, order, energy, mult)
call scattdata_init(this, order, energy, mult)
this % data = coeffs
@ -123,7 +123,7 @@ contains
groups = size(energy,dim=1)
call scattdatabase_init(this, order, energy, mult)
call scattdata_init(this, order, energy, mult)
allocate(this % mu(order))
this % dmu = TWO / real(order,8)
@ -175,7 +175,7 @@ contains
groups = size(energy,dim=1)
call scattdatabase_init(this, this_order, energy, mult)
call scattdata_init(this, this_order, energy, mult)
allocate(this % mu(this_order))
this % dmu = TWO / real(this_order - 1)

View file

@ -110,7 +110,7 @@ contains
integer, save :: num_resamples = 0 ! Number of resamples encountered
! Set weight to one by default
site%wgt = ONE
site % wgt = ONE
! Set the random number generator to the source stream.
call prn_set_stream(STREAM_SOURCE)
@ -118,10 +118,10 @@ contains
! Sample from among multiple source distributions
n_source = size(external_source)
if (n_source > 1) then
r(1) = prn()*sum(external_source(:)%strength)
r(1) = prn()*sum(external_source(:) % strength)
c = ZERO
do i = 1, n_source
c = c + external_source(i)%strength
c = c + external_source(i) % strength
if (r(1) < c) exit
end do
else
@ -132,14 +132,14 @@ contains
found = .false.
do while (.not.found)
! Set particle defaults
call p%initialize()
call p % initialize()
! Sample spatial distribution
site%xyz(:) = external_source(i)%space%sample()
site % xyz(:) = external_source(i) % space % sample()
! Fill p with needed data
p%coord(1)%xyz(:) = site%xyz
p%coord(1)%uvw(:) = [ ONE, ZERO, ZERO ]
p % coord(1) % xyz(:) = site % xyz
p % coord(1) % uvw(:) = [ ONE, ZERO, ZERO ]
! Now search to see if location exists in geometry
call find_cell(p, found)
@ -152,27 +152,27 @@ contains
end if
! Check if spatial site is in fissionable material
select type (space => external_source(i)%space)
select type (space => external_source(i) % space)
type is (SpatialBox)
if (space%only_fissionable) then
if (p%material == MATERIAL_VOID) then
if (space % only_fissionable) then
if (p % material == MATERIAL_VOID) then
found = .false.
elseif (.not. materials(p%material)%fissionable) then
elseif (.not. materials(p % material) % fissionable) then
found = .false.
end if
end if
end select
end do
call p%clear()
call p % clear()
! Sample angle
site%uvw(:) = external_source(i)%angle%sample()
site % uvw(:) = external_source(i) % angle % sample()
! Check for monoenergetic source above maximum neutron energy
select type (energy => external_source(i)%energy)
select type (energy => external_source(i) % energy)
type is (Discrete)
if (any(energy%x >= energy_max_neutron)) then
if (any(energy % x >= energy_max_neutron)) then
call fatal_error("Source energy above range of energies of at least &
&one cross section table")
end if
@ -180,23 +180,23 @@ contains
do
! Sample energy spectrum
site%E = external_source(i)%energy%sample()
site % E = external_source(i) % energy % sample()
! resample if energy is greater than maximum neutron energy
if (site%E < energy_max_neutron) exit
if (site % E < energy_max_neutron) exit
end do
! Set delayed group
site%delayed_group = 0
site % delayed_group = 0
! If running in MG, convert site%E to group
! If running in MG, convert site % E to group
if (.not. run_CE) then
if (site%E <= energy_bins(1)) then
site%g = 1
else if (site%E > energy_bins(energy_groups + 1)) then
site%g = energy_groups
if (site % E <= energy_bins(1)) then
site % g = 1
else if (site % E > energy_bins(energy_groups + 1)) then
site % g = energy_groups
else
site%g = binary_search(energy_bins, energy_groups + 1, site%E)
site % g = binary_search(energy_bins, energy_groups + 1, site % E)
end if
end if