diff --git a/examples/xml/basic/materials.xml b/examples/xml/basic/materials.xml index 2f88731ff..573b96cb2 100644 --- a/examples/xml/basic/materials.xml +++ b/examples/xml/basic/materials.xml @@ -6,13 +6,15 @@ + 294K - + + 294K diff --git a/openmc/data/neutron.py b/openmc/data/neutron.py index 041978742..87415d427 100644 --- a/openmc/data/neutron.py +++ b/openmc/data/neutron.py @@ -399,7 +399,9 @@ class IncidentNeutron(EqualityMixin): g.attrs['A'] = self.mass_number g.attrs['metastable'] = self.metastable g.attrs['atomic_weight_ratio'] = self.atomic_weight_ratio - g.attrs['kTs'] = self.kTs + ktg = g.create_group('kTs') + for i, temperature in enumerate(self.temperatures): + ktg.create_dataset(temperature, data=self.kTs[i]) # Write energy grid eg = g.create_group('energy') @@ -457,7 +459,10 @@ class IncidentNeutron(EqualityMixin): mass_number = group.attrs['A'] metastable = group.attrs['metastable'] atomic_weight_ratio = group.attrs['atomic_weight_ratio'] - kTs = group.attrs['kTs'].tolist() + kTg = group['kTs'] + kTs = [] + for temp in kTg: + kTs.append(temp.value) temperatures = [str(int(round(kT_to_K(kT)))) + "K" for kT in kTs] data = cls(name, atomic_number, mass_number, metastable, diff --git a/openmc/data/thermal.py b/openmc/data/thermal.py index a28dcb504..6f4aadd74 100644 --- a/openmc/data/thermal.py +++ b/openmc/data/thermal.py @@ -225,9 +225,11 @@ class ThermalScattering(EqualityMixin): # Write basic data g = f.create_group(self.name) g.attrs['atomic_weight_ratio'] = self.atomic_weight_ratio - g.attrs['kTs'] = self.kTs g.attrs['zaids'] = self.zaids g.attrs['secondary_mode'] = np.string_(self.secondary_mode) + ktg = g.create_group('kTs') + for i, temperature in enumerate(self.temperatures): + ktg.create_dataset(temperature, data=self.kTs[i]) for T in self.temperatures: Tg = g.create_group(T) @@ -436,7 +438,10 @@ class ThermalScattering(EqualityMixin): name = group.name[1:] atomic_weight_ratio = group.attrs['atomic_weight_ratio'] - kTs = group.attrs['kTs'].tolist() + kTg = group['kTs'] + kTs = [] + for temp in kTg: + kTs.append(temp.value) temperatures = [str(int(round(kT_to_K(kT)))) + "K" for kT in kTs] table = cls(name, atomic_weight_ratio, kTs) diff --git a/src/global.F90 b/src/global.F90 index c22fec25c..0fc1bc0c5 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -77,9 +77,6 @@ module global ! Dictionaries to look up cross sections and listings type(DictCharInt) :: nuclide_dict - ! Default xs identifier (e.g. 70c or 300K) - character(5):: default_xs - ! ============================================================================ ! CONTINUOUS-ENERGY CROSS SECTION RELATED VARIABLES diff --git a/src/input_xml.F90 b/src/input_xml.F90 index a60d8196e..b20927743 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -2159,16 +2159,9 @@ contains &exist!") end if - ! Initialize default cross section variable - default_xs = "" - ! Parse materials.xml file call open_xmldoc(doc, filename) - ! Copy default cross section if present - if (check_for_node(doc, "default_xs")) & - call get_node_value(doc, "default_xs", default_xs) - ! Get pointer to list of XML call get_node_list(doc, "material", node_mat_list) @@ -2204,6 +2197,11 @@ contains call get_node_value(node_mat, "name", mat % name) end if + ! Copy material temperature + if (check_for_node(node_mat, "temperature")) then + call get_node_value(node_mat, "temperature", mat % temperature) + end if + ! ======================================================================= ! READ AND PARSE TAG @@ -2301,22 +2299,9 @@ contains // trim(to_str(mat % id))) end if - ! Check for cross section - 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))) - else - name = to_lower(trim(default_xs)) - end if - end if - - ! store full name - call get_node_value(node_nuc, "name", temp_str) - if (check_for_node(node_nuc, "xs")) & - call get_node_value(node_nuc, "xs", name) - name = trim(temp_str) // "." // trim(name) - name = to_lower(name) + ! store nuclide name + call get_node_value(node_nuc, "name", name) + name = trim(name) ! save name and density to list call names % push_back(name) @@ -2345,16 +2330,6 @@ contains // trim(to_str(mat % id))) end if - ! Check for cross section - 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))) - else - name = to_lower(trim(default_xs)) - end if - end if - ! Check enforced isotropic lab scattering if (run_CE) then if (check_for_node(node_nuc, "scattering")) then @@ -2372,11 +2347,9 @@ contains end if end if - ! store full name - call get_node_value(node_nuc, "name", temp_str) - if (check_for_node(node_nuc, "xs")) & - call get_node_value(node_nuc, "xs", name) - name = trim(temp_str) // "." // trim(name) + ! store nuclide name + call get_node_value(node_nuc, "name", name) + name = trim(name) ! save name and density to list call names % push_back(name) @@ -2424,18 +2397,6 @@ contains end if call get_node_value(node_ele, "name", name) - ! Check for cross section - if (check_for_node(node_ele, "xs")) then - call get_node_value(node_ele, "xs", temp_str) - else - if (default_xs == '') then - call fatal_error("No cross section specified for nuclide in & - &material " // trim(to_str(mat % id))) - else - temp_str = to_lower(trim(default_xs)) - end if - end if - ! 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. & @@ -2581,14 +2542,11 @@ 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 - call fatal_error("Need to specify and for S(a,b) & - &table.") + if (.not. check_for_node(node_sab, "name")) then + call fatal_error("Need to specify for S(a,b) table.") end if call get_node_value(node_sab, "name", name) - call get_node_value(node_sab, "xs", temp_str) - name = trim(name) // "." // trim(temp_str) + name = trim(name) mat % sab_names(j) = name ! Check that this nuclide is listed in the cross_sections.xml file @@ -5856,7 +5814,8 @@ contains ! Read nuclide data from HDF5 file_id = file_open(libraries(i_library) % path, 'r') group_id = open_group(file_id, name) - call nuclides(i_nuclide) % from_hdf5(group_id) + call nuclides(i_nuclide) % from_hdf5(group_id, & + materials(i) % temperature) call close_group(group_id) call file_close(file_id) @@ -5988,7 +5947,7 @@ contains ! Read nuclide data from HDF5 file_id = file_open(libraries(i_library) % path, 'r') group_id = open_group(file_id, name) - call resonant_nuc % from_hdf5(group_id) + call resonant_nuc % from_hdf5(group_id, '0K') call close_group(group_id) call file_close(file_id) diff --git a/src/material_header.F90 b/src/material_header.F90 index be4c860e2..27a0e61fe 100644 --- a/src/material_header.F90 +++ b/src/material_header.F90 @@ -8,11 +8,14 @@ module material_header type Material integer :: id ! unique identifier - character(len=104) :: name = "" ! User-defined name + character(len=104) :: name = "" ! User-defined name integer :: n_nuclides ! number of nuclides integer, allocatable :: nuclide(:) ! index in nuclides array real(8) :: density ! total atom density in atom/b-cm real(8), allocatable :: atom_density(:) ! nuclide atom density in atom/b-cm + character(6) :: temperature ! Temperature of the material + ! as presented in the HDF5 library; + ! e.g., "300K" ! Energy grid information integer :: n_grid ! # of union material grid points diff --git a/src/nuclide_header.F90 b/src/nuclide_header.F90 index f32a93030..fd043d789 100644 --- a/src/nuclide_header.F90 +++ b/src/nuclide_header.F90 @@ -180,9 +180,10 @@ module nuclide_header end subroutine nuclide_clear - subroutine nuclide_from_hdf5(this, group_id) + subroutine nuclide_from_hdf5(this, group_id, temperature) class(Nuclide), intent(inout) :: this - integer(HID_T), intent(in) :: group_id + integer(HID_T), intent(in) :: group_id + character(6), intent(in) :: temperature integer :: i integer :: Z @@ -192,12 +193,14 @@ module nuclide_header integer :: n_links integer :: hdf5_err integer(HID_T) :: urr_group, nu_group - integer(HID_T) :: energy_dset + integer(HID_T) :: energy_group, energy_dset + integer(HID_T) :: kT_group, kT_dset integer(HID_T) :: rxs_group integer(HID_T) :: rx_group integer(HID_T) :: total_nu integer(HID_T) :: fer_group ! fission_energy_release group integer(HID_T) :: fer_dset + integer(HID_T) :: temp_dset integer(SIZE_T) :: name_len, name_file_len integer(HSIZE_T) :: j integer(HSIZE_T) :: dims(1) @@ -217,15 +220,21 @@ module nuclide_header call read_attribute(this % metastable, group_id, 'metastable') this % zaid = 1000*Z + A + 400*this % metastable call read_attribute(this % awr, group_id, 'atomic_weight_ratio') - call read_attribute(this % kT, group_id, 'temperature') + kT_group = open_group(group_id, 'kTs') + kT_dset = open_dataset(kT_group, temperature) + call read_dataset(this % kT, kT_dset) + call close_dataset(kT_dset) + call close_group(kT_group) ! Read energy grid - energy_dset = open_dataset(group_id, 'energy') + energy_group = open_group(group_id, 'energy') + energy_dset = open_dataset(energy_group, temperature) call get_shape(energy_dset, dims) this % n_grid = int(dims(1), 4) allocate(this % energy(this % n_grid)) call read_dataset(this % energy, energy_dset) call close_dataset(energy_dset) + call close_group(energy_group) ! Get MT values based on group names rxs_group = open_group(group_id, 'reactions') @@ -243,7 +252,7 @@ module nuclide_header do i = 1, size(this % reactions) rx_group = open_group(rxs_group, 'reaction_' // trim(& zero_padded(MTs % data(i), 3))) - call this % reactions(i) % from_hdf5(rx_group) + call this % reactions(i) % from_hdf5(rx_group, temperature) call close_group(rx_group) end do call close_group(rxs_group) diff --git a/src/reaction_header.F90 b/src/reaction_header.F90 index 055dcb392..4289e2160 100644 --- a/src/reaction_header.F90 +++ b/src/reaction_header.F90 @@ -29,9 +29,10 @@ module reaction_header contains - subroutine reaction_from_hdf5(this, group_id) + subroutine reaction_from_hdf5(this, group_id, temperature) class(Reaction), intent(inout) :: this integer(HID_T), intent(in) :: group_id + character(6), intent(in) :: temperature integer :: i integer :: cm @@ -41,7 +42,7 @@ contains integer :: n_links integer :: hdf5_err integer(HID_T) :: pgroup - integer(HID_T) :: xs + integer(HID_T) :: xs, xs_group integer(SIZE_T) :: name_len integer(HSIZE_T) :: dims(1) integer(HSIZE_T) :: j @@ -49,16 +50,18 @@ contains call read_attribute(this % Q_value, group_id, 'Q_value') call read_attribute(this % MT, group_id, 'mt') - call read_attribute(this % threshold, group_id, 'threshold_idx') call read_attribute(cm, group_id, 'center_of_mass') this % scatter_in_cm = (cm == 1) - ! Read cross section - xs = open_dataset(group_id, 'xs') + ! Read cross section and threshold_idx data + xs_group = open_group(group_id, temperature) + call read_attribute(this % threshold, xs_group, 'threshold_idx') + xs = open_dataset(xs_group, 'xs') call get_shape(xs, dims) allocate(this % sigma(dims(1))) call read_dataset(this % sigma, xs) call close_dataset(xs) + call close_group(xs_group) ! Determine number of products call h5gget_info_f(group_id, storage_type, n_links, max_corder, hdf5_err)