mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-27 21:55:41 -04:00
Improved performance of cross section reading. Only one pass is made through each ACE library now.
This commit is contained in:
Paul Romano
parent
cbb23472f3
commit
c63ae4755b
2 changed files with 198 additions and 251 deletions
|
|
@ -3,8 +3,8 @@ module cross_section
|
|||
use constants
|
||||
use cross_section_header, only: Nuclide, Reaction, SAB_Table, XsListing
|
||||
use datatypes, only: dict_create, dict_add_key, dict_get_key, &
|
||||
dict_has_key, dict_delete
|
||||
use datatypes_header, only: DictionaryCI
|
||||
dict_has_key, dict_delete, dict_keys
|
||||
use datatypes_header, only: DictionaryCI, ListKeyValueCI
|
||||
use endf, only: reaction_name
|
||||
use error, only: fatal_error
|
||||
use fileio, only: read_line, skip_lines
|
||||
|
|
@ -34,16 +34,24 @@ contains
|
|||
|
||||
subroutine read_xs()
|
||||
|
||||
integer :: i ! index in materials array
|
||||
integer :: j ! index over nuclides in material
|
||||
integer :: index ! index in xs_listings array
|
||||
integer :: index_nuclides ! index in nuclides
|
||||
integer :: index_sab ! index in sab_tables
|
||||
character(10) :: key ! name of isotope, e.g. 92235.03c
|
||||
type(Material), pointer :: mat => null()
|
||||
type(Nuclide), pointer :: nuc => null()
|
||||
type(SAB_Table), pointer :: sab => null()
|
||||
type(DictionaryCI), pointer :: temp_dict => null()
|
||||
integer :: i ! index in materials array
|
||||
integer :: j ! index over nuclides in material
|
||||
integer :: index ! index in xs_listings array
|
||||
integer :: index_nuclides ! index in nuclides
|
||||
integer :: index_sab ! index in sab_tables
|
||||
character(10) :: name ! name of isotope, e.g. 92235.03c
|
||||
character(10) :: alias ! alias of nuclide, e.g. U-235.03c
|
||||
character(MAX_FILE_LEN) :: path ! path to ACE cross section library
|
||||
type(Material), pointer :: mat => null()
|
||||
type(Nuclide), pointer :: nuc => null()
|
||||
type(SAB_Table), pointer :: sab => null()
|
||||
type(DictionaryCI), pointer :: library_dict => null()
|
||||
type(ListKeyValueCI), pointer :: list => null()
|
||||
|
||||
call dict_create(library_dict)
|
||||
|
||||
! ==========================================================================
|
||||
! COUNT NUMBER OF TABLES AND CREATE DICTIONARIES
|
||||
|
||||
! First we need to determine how many continuous-energy tables and how many
|
||||
! S(a,b) thermal scattering tables there are -- this loop doesn't actually
|
||||
|
|
@ -58,89 +66,109 @@ contains
|
|||
! materials -- if not, then increment the count for the number of
|
||||
! nuclides and add to dictionary
|
||||
do j = 1, mat % n_nuclides
|
||||
key = mat % names(j)
|
||||
call lower_case(key)
|
||||
if (.not. dict_has_key(nuclide_dict, key)) then
|
||||
index_nuclides = index_nuclides + 1
|
||||
call dict_add_key(nuclide_dict, key, index_nuclides)
|
||||
name = mat % names(j)
|
||||
|
||||
! First check that this nuclide is listed in the cross_sections.xml
|
||||
! file
|
||||
if (.not. dict_has_key(xs_listing_dict, name)) then
|
||||
message = "Could not find nuclide " // trim(name) // &
|
||||
" in cross_sections.xml file!"
|
||||
call fatal_error()
|
||||
end if
|
||||
|
||||
! Find index in xs_listing and set the name and alias according to the
|
||||
! listing
|
||||
index = dict_get_key(xs_listing_dict, name)
|
||||
name = xs_listings(index) % name
|
||||
alias = xs_listings(index) % alias
|
||||
|
||||
! If this nuclide hasn't been encountered yet, we need to add its name
|
||||
! and alias to the nuclide_dict
|
||||
if (.not. dict_has_key(nuclide_dict, name)) then
|
||||
index_nuclides = index_nuclides + 1
|
||||
mat % nuclide(j) = index_nuclides
|
||||
|
||||
call dict_add_key(nuclide_dict, name, index_nuclides)
|
||||
call dict_add_key(nuclide_dict, alias, index_nuclides)
|
||||
else
|
||||
mat % nuclide(j) = dict_get_key(nuclide_dict, key)
|
||||
mat % nuclide(j) = dict_get_key(nuclide_dict, name)
|
||||
end if
|
||||
|
||||
! Determine the set of ACE files which contain cross-sections used
|
||||
! in this problem
|
||||
path = xs_listings(index) % path
|
||||
if (.not. dict_has_key(library_dict, path)) then
|
||||
call dict_add_key(library_dict, path, 0)
|
||||
end if
|
||||
end do
|
||||
|
||||
! Check if S(a,b) exists on other materials
|
||||
if (mat % has_sab_table) then
|
||||
key = mat % sab_name
|
||||
call lower_case(key)
|
||||
if (.not. dict_has_key(sab_dict, key)) then
|
||||
index_sab = index_sab + 1
|
||||
call dict_add_key(sab_dict, key, index_sab)
|
||||
mat % sab_table = index_sab
|
||||
else
|
||||
mat % sab_table = dict_get_key(sab_dict, key)
|
||||
name = mat % sab_name
|
||||
|
||||
! First check that this nuclide is listed in the cross_sections.xml
|
||||
! file
|
||||
if (.not. dict_has_key(xs_listing_dict, name)) then
|
||||
message = "Could not find S(a,b) table " // trim(name) // &
|
||||
" in cross_sections.xml file!"
|
||||
call fatal_error()
|
||||
end if
|
||||
|
||||
! Find index in xs_listing and set the name and alias according to the
|
||||
! listing
|
||||
index = dict_get_key(xs_listing_dict, name)
|
||||
name = xs_listings(index) % name
|
||||
alias = xs_listings(index) % alias
|
||||
|
||||
! If this sab hasn't been encountered yet, we need to add its name
|
||||
! and alias to the sab_dict
|
||||
if (.not. dict_has_key(sab_dict, name)) then
|
||||
index_sab = index_sab + 1
|
||||
mat % sab_table = index_sab
|
||||
|
||||
call dict_add_key(sab_dict, name, index_sab)
|
||||
call dict_add_key(sab_dict, alias, index_sab)
|
||||
else
|
||||
mat % sab_table = dict_get_key(sab_dict, name)
|
||||
end if
|
||||
|
||||
! Determine the set of ACE files which contain cross-sections used
|
||||
! in this problem
|
||||
path = xs_listings(index) % path
|
||||
if (.not. dict_has_key(library_dict, path)) then
|
||||
call dict_add_key(library_dict, path, 0)
|
||||
end if
|
||||
else
|
||||
mat % sab_table = 0
|
||||
end if
|
||||
end do
|
||||
|
||||
n_nuclides_total = index_nuclides
|
||||
n_sab_tables = index_sab
|
||||
|
||||
! allocate arrays for ACE table storage
|
||||
! allocate arrays for ACE table storage and cross section cache
|
||||
allocate(nuclides(n_nuclides_total))
|
||||
allocate(sab_tables(n_sab_tables))
|
||||
|
||||
! allocate array for microscopic cross section cache
|
||||
allocate(micro_xs(n_nuclides_total))
|
||||
|
||||
call dict_create(temp_dict)
|
||||
! ==========================================================================
|
||||
! READ ALL ACE CROSS SECTION TABLES
|
||||
|
||||
! Get list to list of nuclides
|
||||
list => dict_keys(library_dict)
|
||||
|
||||
! Loop over all files
|
||||
do while (associated(list))
|
||||
! Read all cross sections from this file
|
||||
path = list % data % key
|
||||
call read_ace_tables(path)
|
||||
|
||||
! Move pointer to next item
|
||||
list => list % next
|
||||
end do
|
||||
|
||||
! Now that the nuclides and sab_tables arrays have been allocated, we can
|
||||
! repeat the same loop through each table in each material and read the
|
||||
! cross sections
|
||||
index_nuclides = 0
|
||||
index_sab = 0
|
||||
do i = 1, n_materials
|
||||
mat => materials(i)
|
||||
do j = 1, mat % n_nuclides
|
||||
! Get index in xs_listings array for this nuclide
|
||||
index = mat % xs_listing(j)
|
||||
|
||||
! Get name of nuclide
|
||||
key = mat % names(j)
|
||||
call lower_case(key)
|
||||
|
||||
if (.not. dict_has_key(temp_dict, key)) then
|
||||
index_nuclides = index_nuclides + 1
|
||||
call read_ACE_continuous(index_nuclides, index)
|
||||
call dict_add_key(temp_dict, key, index_nuclides)
|
||||
end if
|
||||
end do
|
||||
|
||||
! Read S(a,b) table if this material has one and it hasn't been read
|
||||
! already
|
||||
if (mat % has_sab_table) then
|
||||
! Get name of S(a,b) table
|
||||
key = mat % sab_name
|
||||
|
||||
if (.not. dict_has_key(temp_dict, key)) then
|
||||
! Find the entry in xs_listings for this table
|
||||
if (dict_has_key(xs_listing_dict, key)) then
|
||||
index = dict_get_key(xs_listing_dict, key)
|
||||
else
|
||||
message = "Cannot find cross-section " // trim(key) // &
|
||||
" in specified cross_sections.xml file."
|
||||
call fatal_error()
|
||||
end if
|
||||
|
||||
! Read the table and add entry to dictionary
|
||||
index_sab = index_sab + 1
|
||||
call read_ACE_thermal(index_sab, index)
|
||||
call dict_add_key(temp_dict, key, index_sab)
|
||||
end if
|
||||
|
||||
! In order to know which nuclide the S(a,b) table applies to, we need
|
||||
! to search through the list of nuclides for one which has a matching
|
||||
! zaid
|
||||
|
|
@ -162,51 +190,35 @@ contains
|
|||
end if
|
||||
end do
|
||||
|
||||
! delete dictionary
|
||||
call dict_delete(temp_dict)
|
||||
call dict_delete(library_dict)
|
||||
|
||||
end subroutine read_xs
|
||||
|
||||
!===============================================================================
|
||||
! READ_ACE_CONTINUOUS reads in a single ACE continuous-energy neutron cross
|
||||
! section table. This routine reads the header data and then calls appropriate
|
||||
! subroutines to parse the actual data.
|
||||
! READ_ACE_TABLES reads in all nuclides and S(a,b) tables that are present in
|
||||
! the problem from a single ACE library. This routine reads the header data for
|
||||
! each table and then calls appropriate subroutines to parse the actual data.
|
||||
!===============================================================================
|
||||
|
||||
subroutine read_ACE_continuous(index_table, index)
|
||||
subroutine read_ace_tables(filename)
|
||||
|
||||
integer, intent(in) :: index_table ! index in nuclides array
|
||||
integer, intent(in) :: index ! index in xs_listings array
|
||||
character(*) :: filename ! name of ACE library file
|
||||
|
||||
integer :: in = 7 ! unit to read from
|
||||
integer :: ioError ! error status for file access
|
||||
integer :: lines ! number of lines (data
|
||||
integer :: n ! number of data values
|
||||
real(8) :: kT ! ACE table temperature
|
||||
logical :: file_exists ! does ACE library exist?
|
||||
logical :: found_xs ! did we find table in library?
|
||||
integer :: index_nuclides ! index in nuclides array
|
||||
integer :: index_sab ! index in sab_tables array
|
||||
integer :: in = 7 ! unit to read from
|
||||
integer :: ioError ! error status for file access
|
||||
integer :: lines ! number of lines (data
|
||||
integer :: n ! number of data values
|
||||
real(8) :: kT ! ACE table temperature
|
||||
logical :: file_exists ! does ACE library exist?
|
||||
logical :: found_nuclide ! did we find nuclide in library?
|
||||
logical :: found_sab ! did we find S(a,b) table in library?
|
||||
character(7) :: readable ! is ACE library readable?
|
||||
character(MAX_LINE_LEN) :: line ! single line to read
|
||||
character(MAX_WORD_LEN) :: words(MAX_WORDS) ! words on a line
|
||||
character(MAX_FILE_LEN) :: filename ! name of ACE library file
|
||||
character(10) :: tablename ! name of cross section table
|
||||
type(Nuclide), pointer :: nuc => null()
|
||||
|
||||
! Check to make sure index in nuclides array and xs_listings arrays are
|
||||
! valid
|
||||
if (index_table > size(nuclides)) then
|
||||
message = "Index of table to read is greater than length of nuclides."
|
||||
call fatal_error()
|
||||
elseif (index > size(xs_listings)) then
|
||||
message = "Index of xs_listing entry is greater than length of " // &
|
||||
"xs_listings."
|
||||
call fatal_error()
|
||||
end if
|
||||
|
||||
filename = xs_listings(index) % path
|
||||
tablename = xs_listings(index) % name
|
||||
|
||||
nuc => nuclides(index_table)
|
||||
type(Nuclide), pointer :: nuc => null()
|
||||
type(SAB_Table), pointer :: sab => null()
|
||||
|
||||
! Check if input file exists and is readable
|
||||
inquire(FILE=filename, EXIST=file_exists, READ=readable)
|
||||
|
|
@ -219,37 +231,73 @@ contains
|
|||
call fatal_error()
|
||||
end if
|
||||
|
||||
! display message
|
||||
message = "Loading ACE cross section table: " // tablename
|
||||
call write_message(6)
|
||||
|
||||
! open file
|
||||
open(file=filename, unit=in, status='old', &
|
||||
& action='read', iostat=ioError)
|
||||
open(FILE=filename, UNIT=in, STATUS='old', &
|
||||
& ACTION='read', IOSTAT=ioError)
|
||||
if (ioError /= 0) then
|
||||
message = "Error while opening file: " // filename
|
||||
call fatal_error()
|
||||
end if
|
||||
|
||||
found_xs = .false.
|
||||
do while (.not. found_xs)
|
||||
do
|
||||
! Reset flags for finding table to read
|
||||
found_nuclide = .false.
|
||||
found_sab = .false.
|
||||
|
||||
! Read next line and tokenize it into strings. If we have reached the end
|
||||
! of the file, then the loop is exited
|
||||
call read_line(in, line, ioError)
|
||||
if (ioError < 0) then
|
||||
message = "Could not find ACE table " // tablename // "."
|
||||
call fatal_error()
|
||||
end if
|
||||
if (ioError < 0) exit
|
||||
call split_string(line, words, n)
|
||||
if (trim(words(1)) == trim(tablename)) then
|
||||
found_xs = .true.
|
||||
|
||||
if (dict_has_key(nuclide_dict, words(1))) then
|
||||
found_nuclide = .true.
|
||||
index_nuclides = dict_get_key(nuclide_dict, words(1))
|
||||
nuc => nuclides(index_nuclides)
|
||||
|
||||
! Copy header data
|
||||
nuc % name = words(1)
|
||||
nuc % awr = str_to_real(words(2))
|
||||
kT = str_to_real(words(3))
|
||||
kT = str_to_real(words(3))
|
||||
nuc % temp = kT / K_BOLTZMANN
|
||||
|
||||
! display message
|
||||
message = "Loading ACE cross section table: " // nuc % name
|
||||
call write_message(6)
|
||||
|
||||
! Skip 5 lines
|
||||
call skip_lines(in, 5, ioError)
|
||||
|
||||
elseif (dict_has_key(sab_dict, words(1))) then
|
||||
found_sab = .true.
|
||||
index_sab = dict_get_key(sab_dict, words(1))
|
||||
sab => sab_tables(index_sab)
|
||||
|
||||
! Copy header data
|
||||
sab % name = words(1)
|
||||
sab % awr = str_to_real(words(2))
|
||||
kT = str_to_real(words(3))
|
||||
sab % temp = kT / K_BOLTZMANN
|
||||
|
||||
! display message
|
||||
message = "Loading ACE cross section table: " // sab % name
|
||||
call write_message(6)
|
||||
|
||||
! Skip 1 line
|
||||
call skip_lines(in, 1, ioError)
|
||||
|
||||
! Read corresponding ZAID
|
||||
call read_line(in, line, ioError)
|
||||
call split_string(line, words, n)
|
||||
sab % zaid = str_to_int(words(1))
|
||||
|
||||
! Skip remaining 3 lines
|
||||
call skip_lines(in, 3, ioError)
|
||||
|
||||
else
|
||||
call skip_lines(in, 5, ioError)
|
||||
end if
|
||||
|
||||
! Skip 5 lines
|
||||
call skip_lines(in, 5, ioError)
|
||||
|
||||
! Read NXS and JXS data
|
||||
read(UNIT=in, FMT=*) NXS
|
||||
read(UNIT=in, FMT=*) JXS
|
||||
|
|
@ -258,7 +306,7 @@ contains
|
|||
n = NXS(1)
|
||||
lines = (n + 3)/4
|
||||
|
||||
if (found_xs) then
|
||||
if (found_nuclide) then
|
||||
! Set ZAID of nuclide
|
||||
nuc % zaid = NXS(2)
|
||||
|
||||
|
|
@ -267,26 +315,38 @@ contains
|
|||
|
||||
! Read XSS
|
||||
read(UNIT=in, FMT=*) XSS
|
||||
|
||||
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
|
||||
deallocate(XSS)
|
||||
elseif (found_sab) then
|
||||
! allocate storage for XSS array
|
||||
allocate(XSS(n))
|
||||
|
||||
! Read XSS
|
||||
read(UNIT=in, FMT=*) XSS
|
||||
|
||||
call read_thermal_data(sab)
|
||||
|
||||
! Free memory from XSS array
|
||||
deallocate(XSS)
|
||||
else
|
||||
call skip_lines(in, lines, ioError)
|
||||
end if
|
||||
|
||||
end do
|
||||
|
||||
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(nuc)) nullify(nuc)
|
||||
if(associated(sab)) nullify(sab)
|
||||
|
||||
close(unit=in)
|
||||
close(UNIT=in)
|
||||
|
||||
end subroutine read_ACE_continuous
|
||||
end subroutine read_ace_tables
|
||||
|
||||
!===============================================================================
|
||||
! READ_ESZ - reads through the ESZ block. This block contains the energy grid,
|
||||
|
|
@ -1057,119 +1117,6 @@ contains
|
|||
|
||||
end subroutine read_unr_res
|
||||
|
||||
!===============================================================================
|
||||
! READ_ACE_THERMAL reads in a single ACE thermal scattering S(a,b) table. This
|
||||
! routine in particular reads the header data and then calls read_thermal_data
|
||||
! to parse the actual data blocks.
|
||||
!===============================================================================
|
||||
|
||||
subroutine read_ACE_thermal(index_table, index)
|
||||
|
||||
integer, intent(in) :: index_table ! index in sab_tables array
|
||||
integer, intent(in) :: index ! index in xs_listings array
|
||||
|
||||
integer :: in = 7 ! unit to read from
|
||||
integer :: ioError ! error status for file access
|
||||
integer :: lines ! number of lines (data
|
||||
integer :: n ! number of data values
|
||||
real(8) :: kT ! ACE table temperature
|
||||
logical :: file_exists ! does ACE library exist?
|
||||
logical :: found_xs ! did we find table in library?
|
||||
character(7) :: readable ! is ACE library readable?
|
||||
character(MAX_LINE_LEN) :: line ! single line to read
|
||||
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(SAB_Table), pointer :: table => null()
|
||||
|
||||
filename = xs_listings(index) % path
|
||||
tablename = xs_listings(index) % name
|
||||
|
||||
table => sab_tables(index_table)
|
||||
|
||||
! Check if input file exists and is readable
|
||||
inquire(FILE=filename, EXIST=file_exists, READ=readable)
|
||||
if (.not. file_exists) then
|
||||
message = "ACE library '" // trim(filename) // "' does not exist!"
|
||||
call fatal_error()
|
||||
elseif (readable(1:3) == 'NO') then
|
||||
message = "ACE library '" // trim(filename) // "' is not readable! &
|
||||
&Change file permissions with chmod command."
|
||||
call fatal_error()
|
||||
end if
|
||||
|
||||
! display message
|
||||
message = "Loading ACE cross section table: " // tablename
|
||||
call write_message(6)
|
||||
|
||||
! open file
|
||||
open(file=filename, unit=in, status='old', &
|
||||
& action='read', iostat=ioError)
|
||||
if (ioError /= 0) then
|
||||
message = "Error while opening file: " // filename
|
||||
call fatal_error()
|
||||
end if
|
||||
|
||||
found_xs = .false.
|
||||
do while (.not. found_xs)
|
||||
call read_line(in, line, ioError)
|
||||
if (ioError < 0) then
|
||||
message = "Could not find ACE table " // tablename // "."
|
||||
call fatal_error()
|
||||
end if
|
||||
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))
|
||||
kT = str_to_real(words(3))
|
||||
table%temp = kT / K_BOLTZMANN
|
||||
|
||||
! Skip 1 line
|
||||
call skip_lines(in, 1, ioError)
|
||||
|
||||
! Read corresponding ZAID
|
||||
call read_line(in, line, ioError)
|
||||
call split_string(line, words, n)
|
||||
table % zaid = str_to_int(words(1))
|
||||
|
||||
! Skip remaining 3 lines
|
||||
call skip_lines(in, 3, ioError)
|
||||
else
|
||||
! Skip 5 lines
|
||||
call skip_lines(in, 5, ioError)
|
||||
end if
|
||||
|
||||
! Read NXS and JXS data
|
||||
read(UNIT=in, FMT=*) NXS
|
||||
read(UNIT=in, FMT=*) JXS
|
||||
|
||||
! Calculate how many data points and lines in the XSS array
|
||||
n = NXS(1)
|
||||
lines = (n + 3)/4
|
||||
|
||||
if (found_xs) then
|
||||
! allocate storage for XSS array
|
||||
allocate(XSS(n))
|
||||
|
||||
! Read XSS
|
||||
read(UNIT=in, FMT=*) XSS
|
||||
else
|
||||
call skip_lines(in, lines, ioError)
|
||||
end if
|
||||
|
||||
end do
|
||||
|
||||
call read_thermal_data(table)
|
||||
|
||||
! Free memory from XSS array
|
||||
if(allocated(XSS)) deallocate(XSS)
|
||||
if(associated(table)) nullify(table)
|
||||
|
||||
close(unit=in)
|
||||
|
||||
end subroutine read_ACE_thermal
|
||||
|
||||
!===============================================================================
|
||||
! READ_THERMAL_DATA reads elastic and inelastic cross sections and corresponding
|
||||
! secondary energy/angle distributions derived from experimental S(a,b)
|
||||
|
|
|
|||
|
|
@ -17,7 +17,7 @@ module datatypes_header
|
|||
!===============================================================================
|
||||
|
||||
! Key length for dictionary
|
||||
integer, parameter :: DICT_KEY_LENGTH = 20
|
||||
integer, parameter :: DICT_KEY_LENGTH = 255
|
||||
|
||||
type KeyValueCI
|
||||
character(len=DICT_KEY_LENGTH) :: key
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue