Fix missing JEFF-4.0 directory-structure TSL PuO2, ThO2, UO2, Zy4

generate_jeff.py

@@ -15,6 +15,87 @@
 import openmc.data
 from utils import download, extract, process_neutron, process_thermal
 
+# Directory-based JEFF-4.0 thermal evaluations.
+DIRECTORY_TSL = {
+    'PuO2': ('O16', 'Pu238', 'Pu239', 'Pu240', 'Pu241', 'Pu242'),
+    'ThO2': ('O16', 'Th232'),
+    'UO2':  ('O16', 'U238'),
+    'Zy4':  ('Sn112', 'Sn114', 'Sn115', 'Sn116', 'Sn117', 'Sn118', 'Sn119', 'Sn120', 'Sn122', 'Sn124',
+             'Zr90', 'Zr91', 'Zr92', 'Zr94', 'Zr96')}
+
+
+def directory_tsl_args(neutron_dir, thermal_dir, output_dir, libver):
+    """Yield arguments for directory-based JEFF-4.0 TSL evaluations."""
+    for material, nuclides in DIRECTORY_TSL.items():
+        material_dir = thermal_dir / f'tsl_{material}'
+
+        burnup_dirs = []
+        if material == 'Zy4':
+            burnup_dirs = sorted(
+                (material_dir / 'Burnup').glob('[0-9]*GWdt'),
+                key=lambda path: int(path.name[:-4]),
+            )
+            if not burnup_dirs:
+                raise FileNotFoundError('No Zy4 burnup evaluations found')
+
+        for nuclide in nuclides:
+            symbol = nuclide.rstrip('0123456789')
+            Z, A, _ = openmc.data.zam(nuclide)
+            path_neutron = (
+                neutron_dir / f'n_{Z}-{symbol}-{A:03d}g.jeff'
+            )
+
+            evaluations = [
+                ('', sorted(
+                    material_dir.glob(
+                        f'[0-9]*K/'
+                        f'tsl_{nuclide}_{material}_[0-9]*K.jeff'
+                    ),
+                    key=lambda path: int(path.parent.name[:-1]),
+                ))
+            ]
+            evaluations.extend(
+                (f'_{path.name}', [
+                    path / f'tsl_{nuclide}_{material}.jeff'
+                ])
+                for path in burnup_dirs
+            )
+
+            if nuclide == 'O16':
+                table_name = f'o{material.lower()}'
+            elif nuclide == 'U238' and material == 'UO2':
+                table_name = 'uuo2'
+            else:
+                table_name = (
+                    f'{nuclide.lower()}{material.lower()}'[:6]
+                )
+
+            isotope_specific = (
+                material == 'Zy4' or (material == 'PuO2' and symbol == 'Pu')
+            )
+            name_part = nuclide if isotope_specific else symbol
+            name = f'c_{name_part}_in_{material}'
+
+            for suffix, paths_thermal in evaluations:
+                missing = [p for p in paths_thermal if not p.is_file()]
+                if not paths_thermal or missing:
+                    raise FileNotFoundError(
+                        f'No TSL evaluations found for '
+                        f'{nuclide} in {material}{suffix}'
+                    )
+
+                yield (
+                    path_neutron,
+                    paths_thermal,
+                    output_dir,
+                    libver,
+                    name + suffix,
+                    table_name,
+                    1000*Z + A,
+                    nuclide,
+                    material != 'Zy4',
+                )
+
 
 class CustomFormatter(argparse.ArgumentDefaultsHelpFormatter,
                       argparse.RawDescriptionHelpFormatter):
@@ -359,6 +440,14 @@ def sort_key(path):
                              args.destination / particle, args.libver)
                 r = pool.apply_async(process_thermal, func_args)
                 results.append(r)
+            # special treatment for directory organized tsl (PuO2, ThO2, UO2, Zy4)
+            for func_args in directory_tsl_args(
+                    neutron_dir,
+                    thermal_dir,
+                    args.destination / particle,
+                    args.libver):
+                r = pool.apply_async(process_thermal, func_args)
+                results.append(r)
 
             for r in results:
                 r.wait()

utils.py

@@ -32,12 +32,8 @@ def process_neutron(path, output_dir, libver, temperatures=None):
     data.export_to_hdf5(h5_file, 'w', libver=libver)
 
 
-def process_thermal(path_neutron, path_thermal, output_dir, libver):
-    """Process ENDF thermal scattering sublibrary file into HDF5 and write into a
-    specified output directory."""
-    print(f'Converting: {path_thermal}')
-
-    # Check if divide_incoherent_elastic should be set
+def _thermal_from_njoy(path_neutron, path_thermal, **kwargs):
+    """Process one thermal scattering evaluation with NJOY."""
     divide_incoherent_elastic = False
     with warnings.catch_warnings(action='error', category=UserWarning):
         try:
@@ -48,13 +44,65 @@ def process_thermal(path_neutron, path_thermal, output_dir, libver):
 
     try:
         with warnings.catch_warnings(action='ignore', category=UserWarning):
-            data = openmc.data.ThermalScattering.from_njoy(
-                path_neutron, path_thermal,
-                divide_incoherent_elastic=divide_incoherent_elastic
+            return openmc.data.ThermalScattering.from_njoy(
+                path_neutron,
+                path_thermal,
+                divide_incoherent_elastic=divide_incoherent_elastic,
+                **kwargs,
             )
     except Exception as e:
         print(path_neutron, path_thermal, e)
         raise
+
+
+def process_thermal(path_neutron, path_thermal, output_dir, libver,
+                    name=None, table_name=None, zaid=None, nuclide=None,
+                    use_endf_data=True):
+    """Process ENDF thermal scattering sublibrary file into HDF5 and write into a
+    specified output directory."""
+
+    if isinstance(path_thermal, (str, Path)):
+        paths_thermal = [path_thermal]
+    else:
+        paths_thermal = path_thermal
+
+    data = None
+
+    for thermal_path in paths_thermal:
+        print(f'Converting: {thermal_path}')
+        # Needed for Zy4 mixed-elastic handling.
+        kwargs = {'use_endf_data': use_endf_data}
+        if table_name is not None:
+            kwargs.update(table_name=table_name, zaids=[zaid], nmix=1)
+
+        new_data = _thermal_from_njoy(path_neutron, thermal_path, **kwargs)
+
+        if name is not None:
+            new_data.name = name
+            new_data.nuclides = [nuclide]
+
+        if data is None:
+            data = new_data
+            continue
+
+        overlap = set(data.temperatures) & set(new_data.temperatures)
+        if overlap:
+            raise ValueError(
+                f'Duplicate temperatures: {sorted(overlap)}'
+            )
+
+        data.kTs.extend(new_data.kTs)
+
+        for reaction_name in ('elastic', 'inelastic'):
+            reaction = getattr(data, reaction_name)
+            new_reaction = getattr(new_data, reaction_name)
+
+            if new_reaction is not None and new_reaction.xs is not None:
+                reaction.xs.update(new_reaction.xs)
+                reaction.distribution.update(new_reaction.distribution)
+
+    data.kTs.sort()
+
     h5_file = output_dir / f'{data.name}.h5'
     print(f'Writing {h5_file} ...')
     data.export_to_hdf5(h5_file, 'w', libver=libver)