updated materials

bright/xsgen/openmc_origen.py

@@ -1,4 +1,45 @@
 from __future__ import print_function
+import os
+import sys
+import subprocess
+
+from pyne import nucname
+
+# templates are from openmc/examples/lattice/simple
+
+SETTINGS_TEMPLATE = """<?xml version="1.0"?>
+<settings>
+  <eigenvalue>
+    <batches>{k_cycles}</batches>
+    <inactive>{k_cycles_skip}</inactive>
+    <particles>{k_particles}</particles>
+  </eigenvalue>
+  <source>
+    <space type="box">
+      <parameters>0 0 0 {unit_cell_height} {unit_cell_height} {unit_cell_height}</parameters>
+    </space>
+  </source>
+</settings>
+"""
+
+MATERIALS_TEMPLATE = """<?xml version="1.0"?>
+<materials>
+  <default_xs>71c</default_xs>
+  <material id="1">
+    <density value="{fuel_density}" units="g/cc" />
+    {_fuel_nucs}
+  </material>
+  <material id="2">
+    <density value="{cool_density}" units="g/cc" />
+    {_cool_nucs}
+    <sab name="HH2O" xs="71t" />
+  </material>
+  <material id="3">
+    <density value="{clad_density}" units="g/cc" />
+    {_clad_nucs}
+  </material>
+</materials>
+"""
 
 class OpenMCOrigen(object):
     """An that combines OpenMC for k-code calculations and ORIGEN for 
@@ -7,6 +48,56 @@ class OpenMCOrigen(object):
 
     def __init__(self, rc):
         self.rc = rc
+        self.statelibs = {}
+        self.builddir = 'build-' + rc.reactor
+        if not os.path.isdir(self.builddir):
+            os.makedirs(self.builddir)
+
+    def pwd(self, state):
+        return os.path.join(self.builddir, str(hash(state)), 'omc')
+
+    def context(self, state):
+        rc = self.rc
+        ctx = dict(rc._dict)
+        ctx.update(zip(rc.perturbation_params, state))
+        return ctx
 
     def generate(self, state):
-        pass
+        """Generates a library for a given state."""
+        if state in self.statelibs:
+            return self.statelibs[state]
+        if state.burn_times != 0.0:
+            raise ValueError("Burn must start at t=0.")
+        self.openmc(state)
+        return self.statelibs[state]
+
+    def openmc(self, state):
+        pwd = self.pwd(state)
+        if not os.path.isdir(pwd):
+            os.makedirs(pwd)
+        self._make_omc_input(state)
+        # run openmc
+        # parse results
+
+    def _make_omc_input(self, state):
+        pwd = self.pwd(state)
+        ctx = self.context(state)
+        # settings
+        settings = SETTINGS_TEMPLATE.format(**ctx)
+        with open(os.path.join(pwd, 'settings.xml'), 'w') as f:
+            f.write(settings)
+        # materials
+        ctx['_fuel_nucs'] = _mat_to_nucs(rc.fuel_material)
+        ctx['_clad_nucs'] = _mat_to_nucs(rc.clad_material)
+        ctx['_cool_nucs'] = _mat_to_nucs(rc.cool_material)
+        materials = MATERIALS_TEMPLATE.format(**ctx)
+        with open(os.path.join(pwd, 'materials.xml'), 'w') as f:
+            f.write(materials)
+
+def _mat_to_nucs(mat):
+    nucs = []
+    template = '<nuclide name="{nuc}" wo="{mass}" />'
+    for nuc, mass in mat.comp:
+        nucs.append(template.format(nuc=nucname.serpent(nuc), mass=mass))
+    nucs = "\n    ".join(nucs)
+    return nucs

bright/xsgen/pre.py

@@ -9,7 +9,7 @@
 
 import numpy as np
 from pyne import nucname
-from pyne.material import Material
+from pyne.material import Material, from_atom_frac
 from pyne.utils import failure
 
 import bright.xsgen.utils as utils
@@ -30,6 +30,8 @@
     'brightlite': BrightliteWriter,
     }
 
+ensure_mat = lambda m: m is isinstance(m, Material) else Material(m)
+
 def run_ui():
     """Runs the cross section user interface."""
     # Test to see if ui library is installed
@@ -87,11 +89,11 @@ def ensure_rc(self, rc):
         self._ensure_bt(rc)
         self._ensure_nl(rc)
         self._ensure_temp(rc)
-        self._ensure_fs(rc)
         self._ensure_smf(rc)
         self._ensure_av(rc)
         self._ensure_inp(rc)
         self._ensure_pp(rc)
+        self._ensure_mats(rc)
 
     def _ensure_bt(self, rc):
         # Make Time Steps
@@ -120,10 +122,6 @@ def _ensure_temp(self, rc):
         # Make temperature
         rc.temperature = rc.get('temperature', 600)
 
-    def _ensure_fs(self, rc):
-        # Make fuel stream
-        rc.ihm_mat = Material(rc.initial_heavy_metal)
-
     def _ensure_smf(self, rc):
         # make sensitivity mass fractions
         if 'sensitivity_mass_fractions' in rc:
@@ -175,6 +173,70 @@ def _ensure_pp(self, rc):
         # burn_times needs to be the last element
         rc.perturbation_params.append('burn_times')
 
+    def _ensure_mats(self, rc):
+        if 'fuel_material'in rc:
+            rc.fuel_material = ensure_mat(rc.fuel_material)
+        elif 'fuel_chemical_form' in rc and 'initial_heavy_metal' in rc:
+            ihm_mat = Material(rc.initial_heavy_metal)
+            atom_frac = {k: v for k, v in rc.fuel_chemical_form.items() if k != "IHM"}
+            atom_frac[ihm_mat] = rc.fuel_chemical_form.get("IHM", 0.0)
+            rc.fuel_material = from_atom_frac(atom_frac)
+        else:
+            raise ValueError("Please specify a fuel.")
+
+        if 'clad_material' in rc: 
+            rc.clad_material = ensure_mat(rc.clad_material)
+        else: 
+            rc.clad_material = Material({
+                # Natural Zirconium
+                400900: 0.98135 * 0.5145,
+                400910: 0.98135 * 0.1122,
+                400920: 0.98135 * 0.1715,
+                400940: 0.98135 * 0.1738,
+                400960: 0.98135 * 0.0280,
+                # The plastic is all melted and the natural Chromium too..
+                240500: 0.00100 * 0.04345,
+                240520: 0.00100 * 0.83789,
+                240530: 0.00100 * 0.09501,
+                240540: 0.00100 * 0.02365,
+                # Natural Iron
+                260540: 0.00135 * 0.05845,
+                260560: 0.00135 * 0.91754,
+                260570: 0.00135 * 0.02119,
+                260580: 0.00135 * 0.00282,
+                # Natural Nickel
+                280580: 0.00055 * 0.68077,
+                280600: 0.00055 * 0.26223,
+                280610: 0.00055 * 0.01140,
+                280620: 0.00055 * 0.03634,
+                280640: 0.00055 * 0.00926,
+                # Natural Tin
+                501120: 0.01450 * 0.0097,
+                501140: 0.01450 * 0.0065,
+                501150: 0.01450 * 0.0034,
+                501160: 0.01450 * 0.1454,
+                501170: 0.01450 * 0.0768,
+                501180: 0.01450 * 0.2422,
+                501190: 0.01450 * 0.0858,
+                501200: 0.01450 * 0.3259,
+                501220: 0.01450 * 0.0463,
+                501240: 0.01450 * 0.0579,
+                # We Need Oxygen!
+                80160:  0.00125,
+                })
+
+        if 'cool_material' in rc: 
+            rc.cool_material = ensure_mat(rc.cool_material)
+        else: 
+            MW = (2 * 1.0) + (1 * 16.0) + (0.199 * 550 * 10.0**-6 * 10.0) + \
+                                          (0.801 * 550 * 10.0**-6 * 11.0)
+            rc.cool_material = Material({{
+                10010: (2 * 1.0) / MW,
+                80160: (1 * 16.0) / MW,
+                50100: (0.199 * 550 * 10.0**-6 * 10.0) / MW,
+                50110: (0.801 * 550 * 10.0**-6 * 11.0) / MW,
+                })
+
     def make_states(self, rc):
         """Makes the reactor state table."""
         State = rc.State = namedtuple('State', rc.perturbation_params)