crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-28 06:05:58 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Slowly getting the tests to work

Browse source
This commit is contained in:
Adam Nelson 2016-09-05 21:04:50 -04:00
parent 1305e99ce3
commit 1625927b44
10 changed files with 106 additions and 54 deletions
Download patch file Download diff file Expand all files Collapse all files

12
scripts/openmc-update-mgxs.py
View file

@ -42,18 +42,20 @@ def parse_args():
# Create argument parser
parser = argparse.ArgumentParser(description=description,
formatter_class=argparse.RawTextHelpFormatter)
parser.add_argument('-c', '--compression', type=int,
help='input XML file')
parser.add_argument('-i', '--input', type=argparse.FileType('r'),
help='input XML file')
parser.add_argument('-o', '--output', nargs='?', default='',
help='output file, in HDF5 format')
parser.add_argument('-c', '--compression', type=int,
help='HDF5 Compression Level')
args = vars(parser.parse_args())
if args['output'] == '':
fname = args['input'].name.split('.')
fname[-1] = '.h5'
args['output'] = ''.join(fname)
filename = args['input'].name
extension = filename[filename.rfind('.'):]
if extension == '.xml':
filename = filename[:filename.rfind('.')] + '.h5'
args['output'] = filename
# Parse and return commandline arguments.
return args

20
src/input_xml.F90
View file

@ -49,10 +49,18 @@ contains
call read_settings_xml()
call read_geometry_xml()
call read_tallies_xml()
call read_materials()
call read_tallies_xml()
if (cmfd_run) call configure_cmfd()
if (.not. run_CE) then
! Create material macroscopic data for MGXS
call time_read_xs % start()
call read_mgxs()
call create_macro_xs()
call time_read_xs % stop()
end if
end subroutine read_input_xml
!===============================================================================
@ -2092,15 +2100,9 @@ contains
call get_temperatures(nuc_temps, sab_temps)
! Read continuous-energy cross sections
if (run_mode /= MODE_PLOTTING) then
if (run_CE .and. run_mode /= MODE_PLOTTING) then
call time_read_xs % start()
if (run_CE) then
call read_ce_cross_sections(libraries, library_dict, nuc_temps, sab_temps)
else
! Create material macroscopic data for MGXS
call read_mgxs(nuc_temps)
call create_macro_xs()
end if
call read_ce_cross_sections(libraries, library_dict, nuc_temps, sab_temps)
call time_read_xs % stop()
end if

54
src/mgxs_data.F90
View file

@ -20,9 +20,7 @@ contains
! nuclides and sab_tables arrays
!===============================================================================
subroutine read_mgxs(temps)
type(VectorReal), allocatable :: temps(:)
subroutine read_mgxs()
integer :: i ! index in materials array
integer :: j ! index over nuclides in material
integer :: i_xsdata ! index in <xsdata> list
@ -38,6 +36,7 @@ contains
logical :: get_kfiss, get_fiss
integer :: l
type(DictCharInt) :: xsdata_dict
type(VectorReal), allocatable :: temps(:)
! Check if cross_sections.xml exists
inquire(FILE=path_cross_sections, EXIST=file_exists)
@ -49,6 +48,9 @@ contains
call write_message("Loading Cross Section Data...", 5)
! Get temperatures
call get_temperatures(temps)
! Open file for reading
file_id = file_open(path_cross_sections, 'r', parallel=.true.)
@ -234,5 +236,51 @@ contains
end subroutine get_mat_kTs
!===============================================================================
! GET_TEMPERATURES returns a list of temperatures that each MGXS table
! appears at in the model. Later, this list is used to determine the actual
! temperatures to read (which may be different if interpolation is used)
!===============================================================================
subroutine get_temperatures(temps)
type(VectorReal), allocatable, intent(out) :: temps(:)
integer :: i, j, k
integer :: i_nuclide ! index in nuclides array
integer :: i_material
real(8) :: temperature ! temperature in Kelvin
allocate(temps(n_nuclides_total))
do i = 1, size(cells)
do j = 1, size(cells(i) % material)
! Skip any non-material cells and void materials
if (cells(i) % material(j) == NONE .or. &
cells(i) % material(j) == MATERIAL_VOID) cycle
! Get temperature of cell (rounding to nearest integer)
if (size(cells(i) % sqrtkT) > 1) then
temperature = cells(i) % sqrtkT(j)**2 / K_BOLTZMANN
else
temperature = cells(i) % sqrtkT(1)**2 / K_BOLTZMANN
end if
i_material = material_dict % get_key(cells(i) % material(j))
associate (mat => materials(i_material))
NUC_NAMES_LOOP: do k = 1, size(mat % names)
! Get index in temps array
i_nuclide = nuclide_dict % get_key(to_lower(mat % names(k)))
! Add temperature if it hasn't already been added
if (find(temps(i_nuclide), temperature) == -1) then
call temps(i_nuclide) % push_back(temperature)
end if
end do NUC_NAMES_LOOP
end associate
end do
end do
end subroutine get_temperatures
end module mgxs_data

58
src/mgxs_header.F90
View file

@ -691,7 +691,7 @@ module mgxs_header
if (check_dataset(xsdata_grp, "nu-fission")) then
allocate(temp_arr(groups * groups * this % n_azi * this % n_pol))
call read_dataset(temp_arr, xsdata_grp, "nu-fission")
allocate(temp_4d(groups, groups, this % n_azi,this % n_pol))
allocate(temp_4d(groups, groups, this % n_azi, this % n_pol))
temp_4d(:, :, :, :) = reshape(temp_arr, (/groups, groups, &
this % n_azi, this % n_pol/))
deallocate(temp_arr)
@ -1241,25 +1241,25 @@ module mgxs_header
! Create the Xs Data for each temperature
TEMP_LOOP: do t = 1, temps % size()
! Allocate and initialize the data needed for macro_xs(i_mat) object
allocate(this % xs(t) % total(n_azi, n_pol, groups))
allocate(this % xs(t) % total(groups, n_azi, n_pol))
this % xs(t) % total = ZERO
allocate(this % xs(t) % absorption(n_azi, n_pol, groups))
allocate(this % xs(t) % absorption(groups, n_azi, n_pol))
this % xs(t) % absorption = ZERO
allocate(this % xs(t) % fission(n_azi, n_pol, groups))
allocate(this % xs(t) % fission(groups, n_azi, n_pol))
this % xs(t) % fission = ZERO
allocate(this % xs(t) % k_fission(n_azi, n_pol, groups))
allocate(this % xs(t) % k_fission(groups, n_azi, n_pol))
this % xs(t) % k_fission = ZERO
allocate(this % xs(t) % nu_fission(n_azi, n_pol, groups))
allocate(this % xs(t) % nu_fission(groups, n_azi, n_pol))
this % xs(t) % nu_fission = ZERO
allocate(this % xs(t) % chi(n_azi, n_pol, groups, groups))
allocate(this % xs(t) % chi(groups, groups, n_azi, n_pol))
this % xs(t) % chi = ZERO
allocate(temp_mult(n_azi, n_pol, groups, groups))
allocate(temp_mult(groups, groups, n_azi, n_pol))
temp_mult = ZERO
allocate(mult_num(n_azi, n_pol, groups, groups))
allocate(mult_num(groups, groups, n_azi, n_pol))
mult_num = ZERO
allocate(mult_denom(n_azi, n_pol, groups, groups))
allocate(mult_denom(groups, groups, n_azi, n_pol))
mult_denom = ZERO
allocate(scatt_coeffs(n_azi, n_pol, order_dim, groups, groups))
allocate(scatt_coeffs(order_dim, groups, groups, n_azi, n_pol))
scatt_coeffs = ZERO
allocate(this % xs(t) % scatter(n_azi, n_pol))
@ -1342,18 +1342,18 @@ module mgxs_header
nuc % xs(nuc_t) % scatter(iazi, ipol) % obj % gmax(gin)
nuscatt = nuc % xs(nuc_t) % scatter(iazi, ipol) % obj % scattxs(gin) * &
nuc % xs(nuc_t) % scatter(iazi, ipol) % obj % energy(gin) % data(gout)
mult_num(iazi, ipol, gout, gin) = &
mult_num(iazi, ipol, gout, gin) + atom_density * &
mult_num(gout, gin, iazi, ipol) = &
mult_num(gout, gin, iazi, ipol) + atom_density * &
nuscatt
if (nuc % xs(nuc_t) % scatter(iazi, ipol) % obj % mult(gin) % data(gout) > ZERO) then
mult_denom(iazi, ipol, gout, gin) = &
mult_denom(iazi, ipol, gout, gin) + atom_density * &
mult_denom(gout, gin, iazi, ipol) = &
mult_denom(gout, gin, iazi, ipol) + atom_density * &
nuscatt / &
nuc % xs(nuc_t) % scatter(iazi, ipol) % obj % mult(gin) % data(gout)
else
! Avoid division by zero
mult_denom(iazi, ipol, gout, gin) = &
mult_denom(iazi, ipol, gout, gin) + atom_density
mult_denom(gout, gin, iazi, ipol) = &
mult_denom(gout, gin, iazi, ipol) + atom_density
end if
end do
end do
@ -1361,8 +1361,8 @@ module mgxs_header
! Get the complete scattering matrix
nuc_order_dim = &
size(nuc % xs(nuc_t) % scatter(iazi, ipol) % obj % dist(1) % data, dim=1)
scatt_coeffs(iazi, ipol, 1:min(nuc_order_dim, order_dim), :, :) = &
scatt_coeffs(iazi, ipol, 1:min(nuc_order_dim, order_dim), :, :) + &
scatt_coeffs(1:min(nuc_order_dim, order_dim), :, :, iazi, ipol) = &
scatt_coeffs(1:min(nuc_order_dim, order_dim), :, :, iazi, ipol) + &
atom_density * &
nuc % xs(nuc_t) % scatter(iazi, ipol) % obj % get_matrix(min(nuc_order_dim, order_dim))
end do
@ -1376,12 +1376,12 @@ module mgxs_header
do iazi = 1, n_azi
do gin = 1, groups
do gout = 1, groups
if (mult_denom(iazi, ipol, gout, gin) > ZERO) then
temp_mult(iazi, ipol, gout, gin) = &
mult_num(iazi, ipol, gout, gin) / &
mult_denom(iazi, ipol, gout, gin)
if (mult_denom(gout, gin, iazi, ipol) > ZERO) then
temp_mult(gout, gin, iazi, ipol) = &
mult_num(gout, gin, iazi, ipol) / &
mult_denom(gout, gin, iazi, ipol)
else
temp_mult(iazi, ipol, gout, gin) = ONE
temp_mult(gout, gin, iazi, ipol) = ONE
end if
end do
end do
@ -1392,8 +1392,8 @@ module mgxs_header
do ipol = 1, n_pol
do iazi = 1, n_azi
call this % xs(t) % scatter(iazi, ipol) % obj % init(&
temp_mult(iazi, ipol, :, :), &
scatt_coeffs(iazi, ipol, :, :, :))
temp_mult(:, :, iazi, ipol), &
scatt_coeffs(:, :, :, iazi, ipol))
end do
end do
@ -1402,10 +1402,10 @@ module mgxs_header
do ipol = 1, n_pol
do iazi = 1, n_azi
do gin = 1, groups
norm = sum(this % xs(t) % chi(iazi, ipol, :, gin))
norm = sum(this % xs(t) % chi(:, gin, iazi, ipol))
if (norm > ZERO) then
this % xs(t) % chi(iazi, ipol, :, gin) = &
this % xs(t) % chi(iazi, ipol, :, gin) / norm
this % xs(t) % chi(:, gin, iazi, ipol) = &
this % xs(t) % chi(:, gin, iazi, ipol) / norm
end if
end do
end do

2
tests/input_set.py
View file

@ -893,7 +893,7 @@ class MGInputSet(InputSet):
self.settings.source = Source(space=Box([0.0, 0.0, 0.0],
[10.0, 10.0, 2.0]))
self.settings.energy_mode = "multi-group"
self.settings.cross_sections = "../1d_mgxs.xml"
self.settings.cross_sections = "../1d_mgxs.h5"
def build_defualt_plots(self):
plot = openmc.Plot()

2
tests/test_mg_basic/inputs_true.dat
View file

@ -1 +1 @@
d90692c6bf8db3672d70103ded682326d7f14d416e79e18ec055bdeafe1dbc1dbe642e8410cad95bf247f27c907b10b1567127bd2868d466ce2586aacff0411f
ed4d8a6c1a00d4dad10bc5522bf42d0c7a13ef5f86fc93a3f249e77681cda4424ce8879919015b05272c276ffea41e3aecef882592e8dbd97085273a356fb630

2
tests/test_mg_max_order/inputs_true.dat
View file

@ -1 +1 @@
7d508b1f3a2661566b8e8cb76fee61aecb96e8b60d633b06f13c4600bd854ea3366cdebd033c0a71ffb8adb90a9aeb64fe5ac0ef3235260921f5689c93e54305
8752d3167989d876c6f29b4b1ca7eed4b91793d8579c92de85633fd90235c67d63ee76d49081006047e0078f6f86eb6e3420c1f8fe44564e7a1fb62e828587d0

2
tests/test_mg_nuclide/inputs_true.dat
View file

@ -1 +1 @@
be296da93031694b2915e1a10e3e6fd663d612cdd29a84745e3ccb0065c088b4bcda45f6919962f6ee784efe52ebe178bf7d3ce019859637ae57c2da1e240a04
10b64fa97cc4feb1b1b8ba401b25b68ab682aba55ad4303a76a9a293ae79cdfff2f65d74310eb917da0aac243c24713d91eafae7ec5f0c4ae535d49fc928a869

2
tests/test_mg_tallies/inputs_true.dat
View file

@ -1 +1 @@
707285b77a091904a69720e78bf87a6779f2223c0fa2a3725353d00ca30f2624c7d53af016d78c7ad30b70a24bb53eda7dea2f592499b098f541591825b67143
34af3b6ee8f65c6257d85a7c9594bc9bbabdabcb055e42bf6dd3ebfc0e41503367674223093363798950524e04f4377eae80f6d535301e55d268ff0c74e4261f

6
tests/test_mgxs_library_ce_to_mg/test_mgxs_library_ce_to_mg.py
View file

@ -61,14 +61,14 @@ class MGXSTestHarness(PyAPITestHarness):
self._input_set.geometry = self.mgxs_lib.create_mg_mode()
# Modify settings so we can run in MG mode
self._input_set.settings.cross_sections = './mgxs.xml'
self._input_set.settings.cross_sections = './mgxs.h5'
self._input_set.settings.energy_mode = 'multi-group'
# Write modified input files
self._input_set.settings.export_to_xml()
self._input_set.geometry.export_to_xml()
self._input_set.materials.export_to_xml()
self._input_set.mgxs_file.export_to_xml()
self._input_set.mgxs_file.export_to_hdf5()
# Dont need tallies.xml, so remove the file
if os.path.exists('./tallies.xml'):
os.remove('./tallies.xml')
@ -85,7 +85,7 @@ class MGXSTestHarness(PyAPITestHarness):
def _cleanup(self):
return
super(MGXSTestHarness, self)._cleanup()
f = os.path.join(os.getcwd(), 'mgxs.xml')
f = os.path.join(os.getcwd(), 'mgxs.h5')
if os.path.exists(f):
os.remove(f)