Changed attributes of Material to better reflect what they actually represent.

This commit is contained in:
Paul Romano 2011-08-01 22:45:17 -04:00
parent f42fd72144
commit 8b3a8111bf
9 changed files with 272 additions and 272 deletions

View file

@ -18,7 +18,7 @@ USE_MPI = no
ifeq ($(USE_MPI),yes)
MPI = /opt/mpich2/1.4-intel
F90 = $(MPI)/bin/mpif90
F90FLAGS += -DMPI -f90=ifort
F90FLAGS += -DMPI
endif
#===============================================================================

View file

@ -161,8 +161,8 @@ module constants
! Versioning numbers
integer, parameter :: VERSION_MAJOR = 0
integer, parameter :: VERSION_MINOR = 2
integer, parameter :: VERSION_RELEASE = 2
integer, parameter :: VERSION_MINOR = 3
integer, parameter :: VERSION_RELEASE = 0
! Maximum number of words in a single line, length of line, and length of
! single word

View file

@ -35,7 +35,7 @@ contains
subroutine read_xs()
integer :: i ! index in materials array
integer :: j ! index over isotopes in material
integer :: j ! index over nuclides in material
integer :: index ! index in xsdatas array
integer :: n ! length of key
integer :: index_nuclides ! index in nuclides
@ -56,9 +56,9 @@ contains
index_sab = 0
do i = 1, n_materials
mat => materials(i)
do j = 1, mat%n_isotopes
index = mat%isotopes(j)
key = xsdatas(index)%id
do j = 1, mat % n_nuclides
index = mat % xsdata(j)
key = xsdatas(index) % id
n = len_trim(key)
call lower_case(key)
select case (key(n:n))
@ -66,9 +66,9 @@ contains
if (.not. dict_has_key(nuclide_dict, key)) then
index_nuclides = index_nuclides + 1
call dict_add_key(nuclide_dict, key, index_nuclides)
mat%table(j) = index_nuclides
mat % nuclide(j) = index_nuclides
else
mat%table(j) = dict_get_key(nuclide_dict, key)
mat % nuclide(j) = dict_get_key(nuclide_dict, key)
end if
case ('t')
! Need same logic as for continuous tables
@ -80,11 +80,11 @@ contains
end do
end do
n_nuclides = index_nuclides
n_sab_tables = index_sab
n_nuclides_total = index_nuclides
n_sab_tables = index_sab
! allocate arrays for ACE table storage
allocate(nuclides(n_nuclides))
allocate(nuclides(n_nuclides_total))
allocate(sab_tables(n_sab_tables))
! loop over all nuclides in xsdata
@ -94,9 +94,9 @@ contains
index_sab = 0
do i = 1, n_materials
mat => materials(i)
do j = 1, mat%n_isotopes
index = mat%isotopes(j)
key = xsdatas(index)%id
do j = 1, mat % n_nuclides
index = mat % xsdata(j)
key = xsdatas(index) % id
n = len_trim(key)
call lower_case(key)
select case (key(n:n))
@ -141,12 +141,12 @@ contains
character(MAX_WORD_LEN) :: words(MAX_WORDS) ! words on a line
character(MAX_WORD_LEN) :: filename ! name of ACE library file
character(10) :: tablename ! name of cross section table
type(Nuclide), pointer :: table => null()
type(Nuclide), pointer :: nuc => null()
filename = xsdatas(index)%path
tablename = xsdatas(index)%id
table => nuclides(index_table)
nuc => nuclides(index_table)
! Check if input file exists and is readable
inquire(FILE=filename, EXIST=file_exists, READ=readable)
@ -181,10 +181,10 @@ contains
call split_string(line, words, n)
if (trim(words(1)) == trim(tablename)) then
found_xs = .true.
table%name = words(1)
table%awr = str_to_real(words(2))
nuc%name = words(1)
nuc%awr = str_to_real(words(2))
kT = str_to_real(words(3))
table%temp = kT / K_BOLTZMANN
nuc%temp = kT / K_BOLTZMANN
end if
! Skip 5 lines
@ -216,16 +216,16 @@ contains
end do
call read_esz(table)
call read_nu_data(table)
call read_reactions(table)
call read_angular_dist(table)
call read_energy_dist(table)
call read_unr_res(table)
call read_esz(nuc)
call read_nu_data(nuc)
call read_reactions(nuc)
call read_angular_dist(nuc)
call read_energy_dist(nuc)
call read_unr_res(nuc)
! Free memory from XSS array
if(allocated(XSS)) deallocate(XSS)
if(associated(table)) nullify(table)
if(associated(nuc)) nullify(nuc)
close(unit=in)
@ -236,33 +236,33 @@ contains
! total xs, absorption xs, elastic scattering xs, and heating numbers.
!===============================================================================
subroutine read_esz(table)
subroutine read_esz(nuc)
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
integer :: NE ! number of energy points for total and elastic cross sections
! determine number of energy points
NE = NXS(3)
table%n_grid = NE
nuc%n_grid = NE
! allocate storage for arrays
allocate(table%energy(NE))
allocate(table%sigma_t(NE))
allocate(table%sigma_a(NE))
allocate(table%sigma_el(NE))
allocate(table%heating(NE))
allocate(nuc%energy(NE))
allocate(nuc%sigma_t(NE))
allocate(nuc%sigma_a(NE))
allocate(nuc%sigma_el(NE))
allocate(nuc%heating(NE))
! read data from XSS -- right now the total, absorption and elastic
! scattering are read in to these special arrays, but in reality, it should
! be necessary to only store elastic scattering and possibly total
! cross-section for total material xs generation.
XSS_index = 1
table%energy = get_real(NE)
table%sigma_t = get_real(NE)
table%sigma_a = get_real(NE)
table%sigma_el = get_real(NE)
table%heating = get_real(NE)
nuc%energy = get_real(NE)
nuc%sigma_t = get_real(NE)
nuc%sigma_a = get_real(NE)
nuc%sigma_el = get_real(NE)
nuc%heating = get_real(NE)
end subroutine read_esz
@ -272,9 +272,9 @@ contains
! down into prompt and delayed neutrons emitted as well.
!===============================================================================
subroutine read_nu_data(table)
subroutine read_nu_data(nuc)
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
integer :: i, j ! loop indices
integer :: JXS2 ! location for fission nu data
@ -301,8 +301,8 @@ contains
if (JXS2 == 0) then
! =======================================================================
! NO PROMPT/TOTAL NU DATA
table % nu_t_type = NU_NONE
table % nu_p_type = NU_NONE
nuc % nu_t_type = NU_NONE
nuc % nu_p_type = NU_NONE
elseif (XSS(JXS2) > 0) then
! =======================================================================
@ -311,16 +311,16 @@ contains
LNU = int(XSS(KNU))
if (LNU == 1) then
! Polynomial data
table % nu_t_type = NU_POLYNOMIAL
table % nu_p_type = NU_NONE
nuc % nu_t_type = NU_POLYNOMIAL
nuc % nu_p_type = NU_NONE
! allocate determine how many coefficients for polynomial
NC = int(XSS(KNU+1))
length = NC + 1
elseif (LNU == 2) then
! Tabular data
table % nu_t_type = NU_TABULAR
table % nu_p_type = NU_NONE
nuc % nu_t_type = NU_TABULAR
nuc % nu_p_type = NU_NONE
! determine number of interpolation regions and number of energies
NR = int(XSS(KNU+1))
@ -329,13 +329,13 @@ contains
end if
! allocate space for nu data storage
allocate(table % nu_t_data(length))
allocate(nuc % nu_t_data(length))
! read data -- for polynomial, this is the number of coefficients and the
! coefficients themselves, and for tabular, this is interpolation data
! and tabular E/nu
XSS_index = KNU + 1
table % nu_t_data = get_real(length)
nuc % nu_t_data = get_real(length)
elseif (XSS(JXS2) < 0) then
! =======================================================================
@ -344,14 +344,14 @@ contains
LNU = XSS(KNU)
if (LNU == 1) then
! Polynomial data
table % nu_p_type = NU_POLYNOMIAL
nuc % nu_p_type = NU_POLYNOMIAL
! allocate determine how many coefficients for polynomial
NC = XSS(KNU+1)
length = NC + 1
elseif (LNU == 2) then
! Tabular data
table % nu_p_type = NU_TABULAR
nuc % nu_p_type = NU_TABULAR
! determine number of interpolation regions and number of energies
NR = XSS(KNU+1)
@ -360,25 +360,25 @@ contains
end if
! allocate space for nu data storage
allocate(table % nu_p_data(length))
allocate(nuc % nu_p_data(length))
! read data
XSS_index = KNU + 1
table % nu_p_data = get_real(length)
nuc % nu_p_data = get_real(length)
! Now read total nu data
KNU = JXS2 + abs(XSS(JXS2)) + 1
LNU = XSS(KNU)
if (LNU == 1) then
! Polynomial data
table % nu_t_type = NU_POLYNOMIAL
nuc % nu_t_type = NU_POLYNOMIAL
! allocate determine how many coefficients for polynomial
NC = int(XSS(KNU+1))
length = NC + 1
elseif (LNU == 2) then
! Tabular data
table % nu_t_type = NU_TABULAR
nuc % nu_t_type = NU_TABULAR
! determine number of interpolation regions and number of energies
NR = int(XSS(KNU+1))
@ -387,18 +387,18 @@ contains
end if
! allocate space for nu data storage
allocate(table % nu_t_data(length))
allocate(nuc % nu_t_data(length))
! read data
XSS_index = KNU + 1
table % nu_t_data = get_real(length)
nuc % nu_t_data = get_real(length)
end if
if (JXS24 > 0) then
! =======================================================================
! DELAYED NU DATA
table % nu_d_type = NU_TABULAR
nuc % nu_d_type = NU_TABULAR
KNU = JXS24
! determine size of tabular delayed nu data
@ -407,19 +407,19 @@ contains
length = 2 + 2*NR + 2*NE
! allocate space for delayed nu data
allocate(table % nu_d_data(length))
allocate(nuc % nu_d_data(length))
! read delayed nu data
XSS_index = KNU + 1
table % nu_d_data = get_real(length)
nuc % nu_d_data = get_real(length)
! =======================================================================
! DELAYED NEUTRON ENERGY DISTRIBUTION
! Allocate space for secondary energy distribution
NPCR = NXS(8)
table % n_precursor = NPCR
allocate(table % nu_d_edist(NPCR))
nuc % n_precursor = NPCR
allocate(nuc % nu_d_edist(NPCR))
LED = JXS(26)
LDIS = JXS(27)
@ -433,30 +433,30 @@ contains
LAW = XSS(LDIS + LOCC)
IDAT = XSS(LDIS + LOCC + 1)
NR = XSS(LDIS + LOCC + 2)
table % nu_d_edist(i) % law = LAW
nuc % nu_d_edist(i) % law = LAW
! allocate space for ENDF interpolation parameters
if (NR > 0) then
allocate(table % nu_d_edist(i) % nbt(NR))
allocate(table % nu_d_edist(i) % int(NR))
allocate(nuc % nu_d_edist(i) % nbt(NR))
allocate(nuc % nu_d_edist(i) % int(NR))
end if
! read ENDF interpolation parameters
XSS_index = LDIS + LOCC + 3
table % nu_d_edist(i) % nbt = get_real(NR)
table % nu_d_edist(i) % int = get_real(NR)
nuc % nu_d_edist(i) % nbt = get_real(NR)
nuc % nu_d_edist(i) % int = get_real(NR)
! allocate space for law validity data
NE = XSS(LDIS + LOCC + 3 + 2*NR)
allocate(table % nu_d_edist(i) % energy(NE))
allocate(table % nu_d_edist(i) % pvalid(NE))
allocate(nuc % nu_d_edist(i) % energy(NE))
allocate(nuc % nu_d_edist(i) % pvalid(NE))
length_interp_data = 5 + 2*(NR + NE)
! read law validity data
XSS_index = LDIS + LOCC + 4 + 2*NR
table % nu_d_edist(i) % energy = get_real(NE)
table % nu_d_edist(i) % pvalid = get_real(NE)
nuc % nu_d_edist(i) % energy = get_real(NE)
nuc % nu_d_edist(i) % pvalid = get_real(NE)
! Set index to beginning of IDAT array
loc = LDIS + IDAT - 2
@ -465,11 +465,11 @@ contains
length = length_energy_dist(loc, LAW, LOCC, length_interp_data)
! allocate secondary energy distribution array
allocate(table % nu_d_edist(i) % data(length))
allocate(nuc % nu_d_edist(i) % data(length))
! read secondary energy distribution
XSS_index = loc + 1
table % nu_d_edist(i) % data = get_real(length)
nuc % nu_d_edist(i) % data = get_real(length)
end do
! =======================================================================
@ -485,14 +485,14 @@ contains
end do
! allocate space for precusor data
allocate(table % nu_d_precursor_data(length))
allocate(nuc % nu_d_precursor_data(length))
! read delayed neutron precursor data
XSS_index = loc
table % nu_d_precursor_data = get_real(length)
nuc % nu_d_precursor_data = get_real(length)
else
table % nu_d_type = NU_NONE
nuc % nu_d_type = NU_NONE
end if
end subroutine read_nu_data
@ -503,9 +503,9 @@ contains
! multiplicities, and cross-sections.
!===============================================================================
subroutine read_reactions(table)
subroutine read_reactions(nuc)
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
integer :: i ! loop indices
integer :: LMT ! index of MT list in XSS
@ -527,20 +527,20 @@ contains
! allocate array of reactions. Add one since we need to include an elastic
! scattering channel
table%n_reaction = NMT + 1
allocate(table%reactions(NMT+1))
nuc%n_reaction = NMT + 1
allocate(nuc%reactions(NMT+1))
! Store elastic scattering cross-section on reaction one
rxn => table%reactions(1)
rxn => nuc%reactions(1)
rxn%MT = 2
rxn%Q_value = ZERO
rxn%TY = 1
rxn%IE = 1
allocate(rxn%sigma(table%n_grid))
rxn%sigma = table%sigma_el
allocate(rxn%sigma(nuc%n_grid))
rxn%sigma = nuc%sigma_el
do i = 1, NMT
rxn => table%reactions(i+1)
rxn => nuc%reactions(i+1)
! read MT number, Q-value, and neutrons produced
rxn % MT = XSS(LMT+i-1)
@ -567,9 +567,9 @@ contains
! secondary neutrons
!===============================================================================
subroutine read_angular_dist(table)
subroutine read_angular_dist(nuc)
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
integer :: JXS8 ! location of angular distribution locators
integer :: JXS9 ! location of angular distributions
@ -588,7 +588,7 @@ contains
! loop over all reactions with secondary neutrons -- NXS(5) does not include
! elastic scattering
do i = 1, NXS(5) + 1
rxn => table%reactions(i)
rxn => nuc%reactions(i)
! find location of angular distribution
LOCB = XSS(JXS8 + i - 1)
@ -610,7 +610,7 @@ contains
allocate(rxn % adist % type(NE))
allocate(rxn % adist % location(NE))
! read incoming energy grid and location of tables
! read incoming energy grid and location of nucs
XSS_index = JXS9 + LOCB
rxn % adist % energy = get_real(NE)
rxn % adist % location = get_int(NE)
@ -657,9 +657,9 @@ contains
! with seconary neutrons (except elastic scattering)
!===============================================================================
subroutine read_energy_dist(table)
subroutine read_energy_dist(nuc)
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
integer :: LED ! location of energy distribution locators
integer :: LDIS ! location of all energy distributions
@ -680,7 +680,7 @@ contains
! Loop over all reactions
do i = 1, NXS(5)
rxn => table % reactions(i+1) ! skip over elastic scattering
rxn => nuc % reactions(i+1) ! skip over elastic scattering
rxn % has_energy_dist = .true.
! find location of energy distribution data
@ -877,9 +877,9 @@ contains
! READ_UNR_RES reads in unresolved resonance probability tables if present.
!===============================================================================
subroutine read_unr_res(table)
subroutine read_unr_res(nuc)
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
integer :: JXS23 ! location of URR data
integer :: loc ! locator
@ -894,36 +894,36 @@ contains
! check if URR data is present
if (JXS23 /= 0) then
table % urr_present = .true.
allocate(table % urr_params(6))
nuc % urr_present = .true.
allocate(nuc % urr_params(6))
loc = JXS23
else
table % urr_present = .false.
nuc % urr_present = .false.
return
end if
! read parameters
table % urr_params(1) = XSS(loc) ! # of incident energies
table % urr_params(2) = XSS(loc + 1) ! # of probabilities
table % urr_params(3) = XSS(loc + 2) ! interpolation parameter
table % urr_params(4) = XSS(loc + 3) ! inelastic competition flag
table % urr_params(5) = XSS(loc + 4) ! other absorption flag
table % urr_params(6) = XSS(loc + 5) ! factors flag
nuc % urr_params(1) = XSS(loc) ! # of incident energies
nuc % urr_params(2) = XSS(loc + 1) ! # of probabilities
nuc % urr_params(3) = XSS(loc + 2) ! interpolation parameter
nuc % urr_params(4) = XSS(loc + 3) ! inelastic competition flag
nuc % urr_params(5) = XSS(loc + 4) ! other absorption flag
nuc % urr_params(6) = XSS(loc + 5) ! factors flag
! allocate incident energies and probability tables
N = table % urr_params(1)
M = table % urr_params(2)
allocate(table % urr_energy(N))
allocate(table % urr_prob(N,6,M))
N = nuc % urr_params(1)
M = nuc % urr_params(2)
allocate(nuc % urr_energy(N))
allocate(nuc % urr_prob(N,6,M))
! read incident energies
XSS_index = loc + 6
table % urr_energy = get_real(N)
nuc % urr_energy = get_real(N)
! read probability tables
do i = 1, N
do j = 1, 6
table % urr_prob(i,j,1:M) = get_real(M)
nuc % urr_prob(i,j,1:M) = get_real(M)
end do
end do
@ -1171,14 +1171,14 @@ contains
real(8) :: xs
integer :: i, j
integer :: n_isotopes
integer :: n_nuclides
integer :: IE
real(8) :: density
real(8) :: density_i
real(8) :: sigma_i
real(8) :: f
type(Nuclide), pointer :: table => null()
type(Reaction), pointer :: rxn => null()
type(Nuclide), pointer :: nuc => null()
type(Reaction), pointer :: rxn => null()
! initialize xs
xs = ZERO
@ -1186,28 +1186,28 @@ contains
! find material atom density
density = cMaterial%atom_density
! loop over all isotopes in material
n_isotopes = cMaterial % n_isotopes
do i = 1, n_isotopes
table => nuclides(cMaterial % table(i))
! loop over all nuclides in material
n_nuclides = cMaterial % n_nuclides
do i = 1, n_nuclides
nuc => nuclides(cMaterial % nuclide(i))
! determine nuclide atom density
density_i = cMaterial%atom_percent(i) * density
! search nuclide energy grid
IE = table%grid_index(p % IE)
f = (p%E - table%energy(IE))/(table%energy(IE+1) - table%energy(IE))
IE = nuc%grid_index(p % IE)
f = (p%E - nuc%energy(IE))/(nuc%energy(IE+1) - nuc%energy(IE))
! handle special case of total cross section
if (MT == 1) then
xs = xs + density * (ONE-f) * table%sigma_t(IE) + &
& f * (table%sigma_t(IE+1))
xs = xs + density * (ONE-f) * nuc%sigma_t(IE) + &
& f * (nuc%sigma_t(IE+1))
cycle
end if
! loop over reactions in isotope
do j = 1, table % n_reaction
rxn => table % reactions(i)
do j = 1, nuc % n_reaction
rxn => nuc % reactions(i)
! check for matching MT
if (MT /= rxn % MT) cycle
@ -1343,7 +1343,7 @@ contains
subroutine material_total_xs()
type(Material), pointer :: mat => null()
type(Nuclide), pointer :: table => null()
type(Nuclide), pointer :: nuc => null()
integer :: i, j, k
integer :: index
integer :: IE
@ -1373,9 +1373,9 @@ contains
mat % total_xs = ZERO
! loop over each isotope
do j = 1, mat%n_isotopes
index = mat % table(j)
table => nuclides(index)
do j = 1, mat % n_nuclides
index = mat % nuclide(j)
nuc => nuclides(index)
! find atom density of isotope
density_i = density * mat % atom_percent(j)
@ -1383,16 +1383,16 @@ contains
! loop over points in union energy grid
do k = 1, n_grid
! find nuclide grid index
IE = table % grid_index(k)
IE = nuc % grid_index(k)
! get nuclide grid energy and cross-section value
E_i = table % energy(IE)
sigma_i = table % sigma_t(IE)
E_i = nuc % energy(IE)
sigma_i = nuc % sigma_t(IE)
! interpolate as necessary
if (E_i < e_grid(k)) then
E_i1 = table % energy(IE + 1)
sigma_i1 = table % sigma_t(IE + 1)
E_i1 = nuc % energy(IE + 1)
sigma_i1 = nuc % sigma_t(IE + 1)
r = (e_grid(k) - E_i)/(E_i1 - E_i)
val = (ONE - r)*sigma_i + r*sigma_i1

View file

@ -21,7 +21,7 @@ contains
type(ListReal), pointer :: list => null()
type(ListReal), pointer :: current => null()
type(Material), pointer :: mat => null()
type(Nuclide), pointer :: table => null()
type(Nuclide), pointer :: nuc => null()
type(Reaction), pointer :: rxn => null()
integer :: i, j
integer :: n
@ -35,12 +35,12 @@ contains
mat => materials(i)
! loop over all isotopes
do j = 1, mat%n_isotopes
table => nuclides(mat%table(j))
do j = 1, mat % n_nuclides
nuc => nuclides(mat % nuclide(j))
! loop over energy points
n = size(table%energy)
call add_grid_points(list, table%energy, n)
n = size(nuc % energy)
call add_grid_points(list, nuc % energy, n)
end do
end do
@ -49,8 +49,8 @@ contains
allocate(e_grid(n_grid))
current => list
do i = 1,n_grid
e_grid(i) = current%data
current => current%next
e_grid(i) = current % data
current => current % next
end do
! delete linked list
@ -84,13 +84,13 @@ contains
allocate(list)
current => list
do index = 1, n_energy
current%data = energy(index)
current % data = energy(index)
if (index == n_energy) then
current%next => null()
current % next => null()
return
end if
allocate(current%next)
current => current%next
allocate(current % next)
current => current % next
end do
end if
@ -105,9 +105,9 @@ contains
if (.not. associated(current)) then
! finish remaining energies
do while (index <= n_energy)
allocate(previous%next)
current => previous%next
current%data = energy(index)
allocate(previous % next)
current => previous % next
current % data = energy(index)
previous => current
index = index + 1
end do
@ -115,13 +115,13 @@ contains
exit
end if
if (E < current%data) then
if (E < current % data) then
! create new element and insert it in energy grid list
allocate(tmp)
tmp%data = E
tmp%next => current
tmp % data = E
tmp % next => current
if (associated(previous)) then
previous%next => tmp
previous % next => tmp
previous => tmp
else
previous => tmp
@ -133,14 +133,14 @@ contains
index = index + 1
E = energy(index)
elseif (E == current%data) then
elseif (E == current % data) then
! found the exact same energy, no need to store duplicates so just
! skip and move to next index
index = index + 1
E = energy(index)
else
previous => current
current => current%next
current => current % next
end if
end do
@ -165,20 +165,20 @@ contains
real(8) :: energy
do i = 1, n_nuclides
do i = 1, n_nuclides_total
nuc => nuclides(i)
allocate(nuc%grid_index(n_grid))
allocate(nuc % grid_index(n_grid))
index = 1
energy = nuc%energy(index)
energy = nuc % energy(index)
do j = 1, n_grid
union_energy = e_grid(j)
if (union_energy >= energy) then
index = index + 1
energy = nuc%energy(index)
energy = nuc % energy(index)
end if
nuc%grid_index(j) = index-1
nuc % grid_index(j) = index-1
end do
end do

View file

@ -979,9 +979,9 @@ contains
character(*), intent(in) :: words(n_words) ! words on material entry
integer, intent(in) :: n_words ! number of words
integer :: i ! index over isotopes
integer :: i ! index over nuclides
integer :: ioError ! error status for file access
integer :: n_isotopes ! number of isotopes in material
integer :: n_nuclides ! number of nuclides in material
character(MAX_LINE_LEN) :: msg ! output/error message
type(Material), pointer :: mat => null()
@ -991,10 +991,10 @@ contains
call fatal_error(msg)
end if
! Determine and set number of isotopes
n_isotopes = (n_words-3)/2
! Determine and set number of nuclides
n_nuclides = (n_words-3)/2
mat => materials(index)
mat % n_isotopes = n_isotopes
mat % n_nuclides = n_nuclides
! Read surface identifier
read(words(2), FMT='(I8)', IOSTAT=ioError) mat % uid
@ -1008,13 +1008,13 @@ contains
mat % atom_density = str_to_real(words(3))
! allocate isotope and density list
allocate(mat%names(n_isotopes))
allocate(mat%isotopes(n_isotopes))
allocate(mat%table(n_isotopes))
allocate(mat%atom_percent(n_isotopes))
allocate(mat%names(n_nuclides))
allocate(mat%xsdata(n_nuclides))
allocate(mat%nuclide(n_nuclides))
allocate(mat%atom_percent(n_nuclides))
! read names and percentage
do i = 1, n_isotopes
do i = 1, n_nuclides
mat % names(i) = words(2*i+2)
mat % atom_percent(i) = str_to_real(words(2*i+3))
end do
@ -1029,14 +1029,14 @@ contains
integer :: index ! index used for several purposes
integer :: i ! index in materials array
integer :: j ! index over isotopes in material
integer :: j ! index over nuclides in material
real(8) :: sum_percent !
real(8) :: awr ! atomic weight ratio
real(8) :: w ! weight percent
real(8) :: x ! atom percent
logical :: percent_in_atom ! isotopes specified in atom percent?
logical :: percent_in_atom ! nuclides specified in atom percent?
logical :: density_in_atom ! density specified in atom/b-cm?
character(10) :: key ! name of isotopes, e.g. 92235.03c
character(10) :: key ! name of nuclide, e.g. 92235.03c
character(MAX_LINE_LEN) :: msg ! output/error message
type(xsData), pointer :: iso => null()
type(Material), pointer :: mat => null()
@ -1058,8 +1058,8 @@ contains
density_in_atom = (mat%atom_density > ZERO)
sum_percent = ZERO
do j = 1, mat % n_isotopes
! Set indices for isotopes
do j = 1, mat % n_nuclides
! Set indices for nuclides
key = mat % names(j)
index = dict_get_key(xsdata_dict, key)
if (index == DICT_NULL) then
@ -1067,7 +1067,7 @@ contains
&xsdata file."
call fatal_error(msg)
end if
mat % isotopes(j) = index
mat % xsdata(j) = index
! determine atomic weight ratio
awr = xsdatas(index) % awr
@ -1090,8 +1090,8 @@ contains
! percent, the sum needs to be re-evaluated as 1/sum(x*awr)
if (.not. density_in_atom) then
sum_percent = ZERO
do j = 1, mat % n_isotopes
index = mat % isotopes(j)
do j = 1, mat % n_nuclides
index = mat % xsdata(j)
awr = xsdatas(index) % awr
x = mat % atom_percent(j)
sum_percent = sum_percent + x*awr

View file

@ -47,7 +47,7 @@ module global
! Cross section arrays
type(Nuclide), allocatable, target :: nuclides(:)
type(SAB_Table), allocatable, target :: sab_tables(:)
integer :: n_nuclides
integer :: n_nuclides_total
integer :: n_sab_tables
! Current cell, surface, material

View file

@ -8,10 +8,10 @@ module material_header
type Material
integer :: uid
integer :: n_isotopes
integer :: n_nuclides
character(10), allocatable :: names(:) ! isotope names
integer, allocatable :: isotopes(:) ! index in xsdata list
integer, allocatable :: table(:) ! index in nuclides array
integer, allocatable :: xsdata(:) ! index in xsdata list
integer, allocatable :: nuclide(:) ! index in nuclides array
real(8) :: atom_density ! total atom density in atom/b-cm
real(8), allocatable :: atom_percent(:)
real(8), allocatable :: total_xs(:) ! macroscopic cross-section

View file

@ -481,21 +481,21 @@ contains
integer :: n_lines
real(8) :: density
character(max_line_len) :: string
type(Nuclide), pointer :: table => null()
type(Nuclide), pointer :: nuc => null()
write(ou,*) 'Material ' // int_to_str(mat % uid)
write(ou,*) ' Atom Density = ' // trim(real_to_str(mat % atom_density)) &
& // ' atom/b-cm'
do i = 1, mat % n_isotopes
table => nuclides(mat % table(i))
do i = 1, mat % n_nuclides
nuc => nuclides(mat % nuclide(i))
density = mat % atom_density * mat % atom_percent(i)
string = ' ' // trim(table % name) // ' = ' // &
string = ' ' // trim(nuc % name) // ' = ' // &
& trim(real_to_str(density)) // ' atom/b-cm'
write(ou,*) trim(string)
end do
write(ou,*)
nullify(table)
nullify(nuc)
end subroutine print_material

View file

@ -150,8 +150,8 @@ contains
type(Particle), pointer :: p
integer :: i ! index over isotopes in a material
integer :: n_isotopes ! number of isotopes in material
integer :: i ! index over nuclides in a material
integer :: n_nuclides ! number of nuclides in material
integer :: IE ! index on nuclide energy grid
real(8) :: f ! interpolation factor
real(8) :: sigma ! microscopic total xs for nuclide
@ -163,45 +163,45 @@ contains
real(8) :: flux ! collision estimator of flux
real(8), allocatable :: Sigma_rxn(:) ! macroscopic xs for each nuclide
character(MAX_LINE_LEN) :: msg ! output/error message
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
type(Reaction), pointer :: rxn
density = cMaterial%atom_density
! calculate total cross-section for each nuclide at current energy in order
! to create discrete pdf for sampling nuclide
n_isotopes = cMaterial%n_isotopes
allocate(Sigma_rxn(n_isotopes))
do i = 1, n_isotopes
table => nuclides(cMaterial%table(i))
n_nuclides = cMaterial % n_nuclides
allocate(Sigma_rxn(n_nuclides))
do i = 1, n_nuclides
nuc => nuclides(cMaterial % nuclide(i))
! determine nuclide atom density
density_i = cMaterial%atom_percent(i) * density
density_i = cMaterial % atom_percent(i) * density
! search nuclide energy grid
IE = table%grid_index(p % IE)
f = (p%E - table%energy(IE))/(table%energy(IE+1) - table%energy(IE))
IE = nuc%grid_index(p % IE)
f = (p%E - nuc%energy(IE))/(nuc%energy(IE+1) - nuc%energy(IE))
! calculate nuclide macroscopic cross-section
Sigma_rxn(i) = density_i*((ONE-f)*table%Sigma_t(IE) + &
& f*(table%Sigma_t(IE+1)))
Sigma_rxn(i) = density_i*((ONE-f)*nuc%Sigma_t(IE) + &
& f*(nuc%Sigma_t(IE+1)))
end do
total = sum(Sigma_rxn)
! sample nuclide
r1 = rang()
prob = ZERO
do i = 1, n_isotopes
do i = 1, n_nuclides
prob = prob + Sigma_rxn(i) / total
if (r1 < prob) exit
end do
! Get table, total xs, interpolation factor
density_i = cMaterial%atom_percent(i)*density
table => nuclides(cMaterial%table(i))
nuc => nuclides(cMaterial%nuclide(i))
sigma = Sigma_rxn(i) / density_i
IE = table%grid_index(p % IE)
f = (p%E - table%energy(IE))/(table%energy(IE+1) - table%energy(IE))
IE = nuc%grid_index(p % IE)
f = (p%E - nuc%energy(IE))/(nuc%energy(IE+1) - nuc%energy(IE))
! free memory
deallocate(Sigma_rxn)
@ -209,8 +209,8 @@ contains
! sample reaction channel
r1 = rang()*sigma
prob = ZERO
do i = 1, table%n_reaction
rxn => table%reactions(i)
do i = 1, nuc%n_reaction
rxn => nuc%reactions(i)
! some materials have gas production cross sections with MT > 200 that
! are duplicates. Also MT=4 is total level inelastic scattering which
@ -228,21 +228,21 @@ contains
if (verbosity >= 10) then
msg = " " // trim(reaction_name(rxn%MT)) // " with nuclide " // &
& trim(table%name)
& trim(nuc%name)
call message(msg)
end if
! call appropriate subroutine
select case (rxn % MT)
case (ELASTIC)
call elastic_scatter(p, table, rxn)
call elastic_scatter(p, nuc, rxn)
case (N_2ND, N_2N, N_3N, N_2NA, N_3NA, N_2N2A, N_4N, N_2NP, N_3NP)
call n_xn(p, table, rxn)
call n_xn(p, nuc, rxn)
case (N_NA, N_N3A, N_NP, N_N2A, N_ND, N_NT, N_N3HE, &
& N_NT2A, N_N2P, N_NPA, N_N1 : N_NC)
call inelastic_scatter(p, table, rxn)
call inelastic_scatter(p, nuc, rxn)
case (FISSION, N_F, N_NF, N_2NF, N_3NF)
call n_fission(p, table, rxn)
call n_fission(p, nuc, rxn)
case (N_GAMMA : N_DA)
call n_absorption(p)
case default
@ -274,10 +274,10 @@ contains
! need to be fixed
!===============================================================================
subroutine elastic_scatter(p, table, rxn)
subroutine elastic_scatter(p, nuc, rxn)
type(Particle), pointer :: p
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
type(Reaction), pointer :: rxn
real(8) :: awr ! atomic weight ratio of target
@ -296,7 +296,7 @@ contains
real(8) :: E ! energy
vel = sqrt(p % E)
awr = table % awr
awr = nuc % awr
vx = vel * p%uvw(1)
vy = vel * p%uvw(2)
@ -349,10 +349,10 @@ contains
! with implicit absorption, namely sampling of the number of neutrons!
!===============================================================================
subroutine n_fission(p, table, rxn)
subroutine n_fission(p, nuc, rxn)
type(Particle), pointer :: p
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
type(Reaction), pointer :: rxn
integer :: i ! loop index
@ -384,24 +384,24 @@ contains
! ==========================================================================
! DETERMINE TOTAL NU
if (table % nu_t_type == NU_NONE) then
msg = "No neutron emission data for table: " // table % name
if (nuc % nu_t_type == NU_NONE) then
msg = "No neutron emission data for table: " // nuc % name
call fatal_error(msg)
elseif (table % nu_t_type == NU_POLYNOMIAL) then
elseif (nuc % nu_t_type == NU_POLYNOMIAL) then
! determine number of coefficients
NC = int(table % nu_t_data(1))
NC = int(nuc % nu_t_data(1))
! sum up polynomial in energy
nu_total = ZERO
do i = 0, NC - 1
c = table % nu_t_data(i+2)
c = nuc % nu_t_data(i+2)
nu_total = nu_total + c * E**i
end do
elseif (table % nu_t_type == NU_TABULAR) then
elseif (nuc % nu_t_type == NU_TABULAR) then
! determine number of interpolation regions -- as far as I can tell, no
! nu data has multiple interpolation regions. Furthermore, it seems all
! are lin-lin.
NR = int(table % nu_t_data(1))
NR = int(nuc % nu_t_data(1))
if (NR /= 0) then
msg = "Multiple interpolation regions not supported while &
&attempting to determine total nu."
@ -410,47 +410,47 @@ contains
! determine number of energies
loc = 2 + 2*NR
NE = int(table % nu_t_data(loc))
NE = int(nuc % nu_t_data(loc))
! do binary search over tabuled energies to determine appropriate index
! and interpolation factor
if (E < table % nu_t_data(loc+1)) then
if (E < nuc % nu_t_data(loc+1)) then
j = 1
f = ZERO
elseif (E > table % nu_t_data(loc+NE)) then
elseif (E > nuc % nu_t_data(loc+NE)) then
j = NE - 1
f = ONE
else
j = binary_search(table % nu_t_data(loc+1), NE, E)
f = (E - table % nu_t_data(loc+j)) / &
& (table % nu_t_data(loc+j+1) - table % nu_t_data(loc+j))
j = binary_search(nuc % nu_t_data(loc+1), NE, E)
f = (E - nuc % nu_t_data(loc+j)) / &
& (nuc % nu_t_data(loc+j+1) - nuc % nu_t_data(loc+j))
end if
! determine nu total
loc = loc + NE
nu_total = table % nu_t_data(loc+j) + f * &
& (table % nu_t_data(loc+j+1) - table % nu_t_data(loc+j))
nu_total = nuc % nu_t_data(loc+j) + f * &
& (nuc % nu_t_data(loc+j+1) - nuc % nu_t_data(loc+j))
end if
! ==========================================================================
! DETERMINE PROMPT NU
if (table % nu_p_type == NU_NONE) then
if (nuc % nu_p_type == NU_NONE) then
! since no prompt or delayed data is present, this means all neutron
! emission is prompt
nu_prompt = nu_total
elseif (table % nu_p_type == NU_POLYNOMIAL) then
elseif (nuc % nu_p_type == NU_POLYNOMIAL) then
! determine number of coefficients
NC = int(table % nu_p_data(1))
NC = int(nuc % nu_p_data(1))
! sum up polynomial in energy
nu_prompt = ZERO
do i = 0, NC - 1
c = table % nu_p_data(i+2)
c = nuc % nu_p_data(i+2)
nu_prompt = nu_prompt + c * E**i
end do
elseif (table % nu_p_type == NU_TABULAR) then
elseif (nuc % nu_p_type == NU_TABULAR) then
! determine number of interpolation regions
NR = int(table % nu_p_data(1))
NR = int(nuc % nu_p_data(1))
if (NR /= 0) then
msg = "Multiple interpolation regions not supported while &
&attempting to determine prompt nu."
@ -459,35 +459,35 @@ contains
! determine number of energies
loc = 2 + 2*NR
NE = int(table % nu_p_data(loc))
NE = int(nuc % nu_p_data(loc))
! do binary search over tabuled energies to determine appropriate index
! and interpolation factor
if (E < table % nu_p_data(loc+1)) then
if (E < nuc % nu_p_data(loc+1)) then
j = 1
f = ZERO
elseif (E > table % nu_p_data(loc+NE)) then
elseif (E > nuc % nu_p_data(loc+NE)) then
j = NE - 1
f = ONE
else
j = binary_search(table % nu_p_data(loc+1), NE, E)
f = (E - table % nu_p_data(loc+j)) / &
& (table % nu_p_data(loc+j+1) - table % nu_p_data(loc+j))
j = binary_search(nuc % nu_p_data(loc+1), NE, E)
f = (E - nuc % nu_p_data(loc+j)) / &
& (nuc % nu_p_data(loc+j+1) - nuc % nu_p_data(loc+j))
end if
! determine nu total
loc = loc + NE
nu_prompt = table % nu_p_data(loc+j) + f * &
& (table % nu_p_data(loc+j+1) - table % nu_p_data(loc+j))
nu_prompt = nuc % nu_p_data(loc+j) + f * &
& (nuc % nu_p_data(loc+j+1) - nuc % nu_p_data(loc+j))
end if
! ==========================================================================
! DETERMINE DELAYED NU
if (table % nu_d_type == NU_NONE) then
if (nuc % nu_d_type == NU_NONE) then
nu_delay = ZERO
elseif (table % nu_d_type == NU_TABULAR) then
elseif (nuc % nu_d_type == NU_TABULAR) then
! determine number of interpolation regions
NR = int(table % nu_d_data(1))
NR = int(nuc % nu_d_data(1))
if (NR /= 0) then
msg = "Multiple interpolation regions not supported while &
&attempting to determine delayed nu."
@ -496,26 +496,26 @@ contains
! determine number of energies
loc = 2 + 2*NR
NE = int(table % nu_d_data(loc))
NE = int(nuc % nu_d_data(loc))
! do binary search over tabuled energies to determine appropriate index
! and interpolation factor
if (E < table % nu_d_data(loc+1)) then
if (E < nuc % nu_d_data(loc+1)) then
j = 1
f = ZERO
elseif (E > table % nu_d_data(loc+NE)) then
elseif (E > nuc % nu_d_data(loc+NE)) then
j = NE - 1
f = ONE
else
j = binary_search(table % nu_d_data(loc+1), NE, E)
f = (E - table % nu_d_data(loc+j)) / &
& (table % nu_d_data(loc+j+1) - table % nu_d_data(loc+j))
j = binary_search(nuc % nu_d_data(loc+1), NE, E)
f = (E - nuc % nu_d_data(loc+j)) / &
& (nuc % nu_d_data(loc+j+1) - nuc % nu_d_data(loc+j))
end if
! determine nu total
loc = loc + NE
nu_delay = table % nu_d_data(loc+j) + f * &
& (table % nu_d_data(loc+j+1) - table % nu_d_data(loc+j))
nu_delay = nuc % nu_d_data(loc+j) + f * &
& (nuc % nu_d_data(loc+j+1) - nuc % nu_d_data(loc+j))
end if
beta = nu_delay / nu_total
@ -549,10 +549,10 @@ contains
xi = rang()
loc = 1
prob = ZERO
do j = 1, table % n_precursor
do j = 1, nuc % n_precursor
! determine number of interpolation regions and energies
NR = table % nu_d_precursor_data(loc + 1)
NE = table % nu_d_precursor_data(loc + 2 + 2*NR)
NR = nuc % nu_d_precursor_data(loc + 1)
NE = nuc % nu_d_precursor_data(loc + 2 + 2*NR)
if (NR > 0) then
msg = "Multiple interpolation regions not supported while &
&sampling delayed neutron precursor yield."
@ -561,24 +561,24 @@ contains
! interpolate on energy grid
loc = loc + 2 + 2*NR
if (E < table%nu_d_precursor_data(loc+1)) then
if (E < nuc%nu_d_precursor_data(loc+1)) then
k = 1
f = ZERO
elseif (E > table%nu_d_precursor_data(loc+NE)) then
elseif (E > nuc%nu_d_precursor_data(loc+NE)) then
k = NE - 1
f = ONE
else
k = binary_search(table%nu_d_precursor_data(loc+1), NE, E)
f = (E - table%nu_d_precursor_data(loc+k)) / &
& (table%nu_d_precursor_data(loc+k+1) - &
& table%nu_d_precursor_data(loc+k))
k = binary_search(nuc%nu_d_precursor_data(loc+1), NE, E)
f = (E - nuc%nu_d_precursor_data(loc+k)) / &
& (nuc%nu_d_precursor_data(loc+k+1) - &
& nuc%nu_d_precursor_data(loc+k))
end if
! determine delayed neutron precursor yield for group j
loc = loc + NE
yield = table%nu_d_precursor_data(loc+k) + f * &
(table%nu_d_precursor_data(loc+k+1) - &
& table%nu_d_precursor_data(loc+k))
yield = nuc%nu_d_precursor_data(loc+k) + f * &
(nuc%nu_d_precursor_data(loc+k+1) - &
& nuc%nu_d_precursor_data(loc+k))
prob = prob + yield
if (xi < prob) exit
@ -587,12 +587,12 @@ contains
end do
! sample from energy distribution for group j
law = table % nu_d_edist(j) % law
law = nuc % nu_d_edist(j) % law
do
if (law == 44 .or. law == 61) then
call sample_energy(table%nu_d_edist(j), E, E_out, mu)
call sample_energy(nuc%nu_d_edist(j), E, E_out, mu)
else
call sample_energy(table%nu_d_edist(j), E, E_out)
call sample_energy(nuc%nu_d_edist(j), E, E_out)
end if
! resample if energy is >= 20 MeV
if (E_out < 20) exit
@ -639,10 +639,10 @@ contains
! than fission), i.e. level scattering, (n,np), (n,na), etc.
!===============================================================================
subroutine inelastic_scatter(p, table, rxn)
subroutine inelastic_scatter(p, nuc, rxn)
type(Particle), pointer :: p
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
type(Reaction), pointer :: rxn
integer :: law
@ -659,7 +659,7 @@ contains
E_in = p % E
! determine A and Q
A = table % awr
A = nuc % awr
Q = rxn % Q_value
! determine secondary energy distribution law
@ -708,10 +708,10 @@ contains
! N_XN
!===============================================================================
subroutine n_xn(p, table, rxn)
subroutine n_xn(p, nuc, rxn)
type(Particle), pointer :: p
type(Nuclide), pointer :: table
type(Nuclide), pointer :: nuc
type(Reaction), pointer :: rxn
integer :: i ! loop index
@ -730,7 +730,7 @@ contains
E_in = p % E
! determine A and Q
A = table % awr
A = nuc % awr
Q = rxn % Q_value
! determine secondary energy distribution law