1from collections.abc import Iterable 2from numbers import Real 3 4import numpy as np 5 6import openmc.checkvalue as cv 7from openmc.mixin import EqualityMixin 8from .angle_energy import AngleEnergy 9from .function import Tabulated1D, Polynomial, Function1D 10 11 12class Product(EqualityMixin): 13 """Secondary particle emitted in a nuclear reaction 14 15 Parameters 16 ---------- 17 particle : str, optional 18 The particle type of the reaction product. Defaults to 'neutron'. 19 20 Attributes 21 ---------- 22 applicability : Iterable of openmc.data.Tabulated1D 23 Probability of sampling a given distribution for this product. 24 decay_rate : float 25 Decay rate in inverse seconds 26 distribution : Iterable of openmc.data.AngleEnergy 27 Distributions of energy and angle of product. 28 emission_mode : {'prompt', 'delayed', 'total'} 29 Indicate whether the particle is emitted immediately or whether it 30 results from the decay of reaction product (e.g., neutron emitted from a 31 delayed neutron precursor). A special value of 'total' is used when the 32 yield represents particles from prompt and delayed sources. 33 particle : str 34 The particle type of the reaction product 35 yield_ : openmc.data.Function1D 36 Yield of secondary particle in the reaction. 37 38 """ 39 40 def __init__(self, particle='neutron'): 41 self.applicability = [] 42 self.decay_rate = 0.0 43 self.distribution = [] 44 self.emission_mode = 'prompt' 45 self.particle = particle 46 self.yield_ = Polynomial((1,)) # 0-order polynomial, i.e., a constant 47 48 def __repr__(self): 49 if isinstance(self.yield_, Tabulated1D): 50 if np.all(self.yield_.y == self.yield_.y[0]): 51 return "<Product: {}, emission={}, yield={}>".format( 52 self.particle, self.emission_mode, self.yield_.y[0]) 53 else: 54 return "<Product: {}, emission={}, yield=tabulated>".format( 55 self.particle, self.emission_mode) 56 else: 57 return "<Product: {}, emission={}, yield=polynomial>".format( 58 self.particle, self.emission_mode) 59 60 @property 61 def applicability(self): 62 return self._applicability 63 64 @property 65 def decay_rate(self): 66 return self._decay_rate 67 68 @property 69 def distribution(self): 70 return self._distribution 71 72 @property 73 def emission_mode(self): 74 return self._emission_mode 75 76 @property 77 def particle(self): 78 return self._particle 79 80 @property 81 def yield_(self): 82 return self._yield 83 84 @applicability.setter 85 def applicability(self, applicability): 86 cv.check_type('product distribution applicability', applicability, 87 Iterable, Tabulated1D) 88 self._applicability = applicability 89 90 @decay_rate.setter 91 def decay_rate(self, decay_rate): 92 cv.check_type('product decay rate', decay_rate, Real) 93 cv.check_greater_than('product decay rate', decay_rate, 0.0, True) 94 self._decay_rate = decay_rate 95 96 @distribution.setter 97 def distribution(self, distribution): 98 cv.check_type('product angle-energy distribution', distribution, 99 Iterable, AngleEnergy) 100 self._distribution = distribution 101 102 @emission_mode.setter 103 def emission_mode(self, emission_mode): 104 cv.check_value('product emission mode', emission_mode, 105 ('prompt', 'delayed', 'total')) 106 self._emission_mode = emission_mode 107 108 @particle.setter 109 def particle(self, particle): 110 cv.check_type('product particle type', particle, str) 111 self._particle = particle 112 113 @yield_.setter 114 def yield_(self, yield_): 115 cv.check_type('product yield', yield_, Function1D) 116 self._yield = yield_ 117 118 def to_hdf5(self, group): 119 """Write product to an HDF5 group 120 121 Parameters 122 ---------- 123 group : h5py.Group 124 HDF5 group to write to 125 126 """ 127 group.attrs['particle'] = np.string_(self.particle) 128 group.attrs['emission_mode'] = np.string_(self.emission_mode) 129 if self.decay_rate > 0.0: 130 group.attrs['decay_rate'] = self.decay_rate 131 132 # Write yield 133 self.yield_.to_hdf5(group, 'yield') 134 135 # Write applicability/distribution 136 group.attrs['n_distribution'] = len(self.distribution) 137 for i, d in enumerate(self.distribution): 138 dgroup = group.create_group('distribution_{}'.format(i)) 139 if self.applicability: 140 self.applicability[i].to_hdf5(dgroup, 'applicability') 141 d.to_hdf5(dgroup) 142 143 @classmethod 144 def from_hdf5(cls, group): 145 """Generate reaction product from HDF5 data 146 147 Parameters 148 ---------- 149 group : h5py.Group 150 HDF5 group to read from 151 152 Returns 153 ------- 154 openmc.data.Product 155 Reaction product 156 157 """ 158 particle = group.attrs['particle'].decode() 159 p = cls(particle) 160 161 p.emission_mode = group.attrs['emission_mode'].decode() 162 if 'decay_rate' in group.attrs: 163 p.decay_rate = group.attrs['decay_rate'] 164 165 # Read yield 166 p.yield_ = Function1D.from_hdf5(group['yield']) 167 168 # Read applicability/distribution 169 n_distribution = group.attrs['n_distribution'] 170 distribution = [] 171 applicability = [] 172 for i in range(n_distribution): 173 dgroup = group['distribution_{}'.format(i)] 174 if 'applicability' in dgroup: 175 applicability.append(Tabulated1D.from_hdf5( 176 dgroup['applicability'])) 177 distribution.append(AngleEnergy.from_hdf5(dgroup)) 178 179 p.distribution = distribution 180 p.applicability = applicability 181 182 return p 183