Update gas file output to use alpha/eta SST values

PyBoltz/Boltz.pxd

@@ -11,6 +11,8 @@ cdef class Boltz:
     cpdef Start(self)
 
     cdef public:
+        bint Crossed_SST
+        '''This is flag to mark if the steady state threshold has been crossed.'''
         double EFieldOverBField
         '''This is a constant that is equal to the electric field / magentic field * 1e-9.'''
         double AngularSpeedOfRotation
@@ -414,4 +416,4 @@ cdef class Boltz:
         double IonCollisionFreq[4000]
         '''Array that adds the ionisation collision frequencies at energy step I. Used in the Friedland estimation of Alpha.'''
         double AttCollisionFreq[4000]
-        '''Array that adds the attachment collision frequencies at energy step I. Used in the Friedland estimation of Alpha.'''
+        '''Array that adds the attachment collision frequencies at energy step I. Used in the Friedland estimation of Alpha.'''

PyBoltz/Boltz.pyx

@@ -195,6 +195,7 @@ cdef class Boltz:
         self.ReducedIonizationErr = 0.0
         self.ReducedAttachmentErr = 0.0
         self.Steady_State_Threshold = 40.0
+        self.Crossed_SST = False
         self.MixObject = Gasmix()
 
     def reset(self):
@@ -377,10 +378,13 @@ cdef class Boltz:
         self.end()
         # Steady state
         if abs(self.ReducedIonization - self.ReducedAttachment) >= self.Steady_State_Threshold:
-            if self.ReducedIonization ==0:
+            if self.ReducedIonization == 0:
                 print("Steady State Threshold has been crossed. Will not run the SST simulation as the ionisation rate is zero.")
                 return
-            if self.Console_Output_Flag: print("\n**Crossed the set Steady state simulation threshold = {}\n".format(self.Steady_State_Threshold))
+            if self.Console_Output_Flag:
+                print("\n**Crossed the set Steady state simulation threshold = {}\n".format(self.Steady_State_Threshold))
+
             TownsendFunc.run(self)
+            self.Crossed_SST = True
 
         return

PyBoltz/OdieRun.py

@@ -4,7 +4,6 @@
 from PyBoltz.Boltz import Boltz
 from PyBoltz.PyBoltzRun import PBRes
 
-
 class OdieRun:
     '''Class to run PyBoltz and provide output for the Garfield++ package'''
 
@@ -29,9 +28,10 @@ class OdieRun:
     }
 
     Gases = [
-        np.nan, 'CF4', 'ARGON', 'HELIUM4', 'HELIUM3', 'NEON', 'KRYPTON', 'XENON', 'CH4', 'ETHANE',
-        'PROPANE', 'ISOBUTANE', 'CO2', np.nan, 'H2O', 'OXYGEN', 'NITROGEN', np.nan, np.nan, np.nan,
-        np.nan, 'HYDROGEN', 'DEUTERIUM', np.nan, np.nan, 'DME'
+        [], ['CF4'], ['ARGON', 'AR'], ['HELIUM4', 'HE4'], ['HELIUM3', 'HE3'], ['NEON', 'NE'],
+        ['KRYPTON', 'KR'], ['XENON', 'XE'], ['METHANE', 'CH4'], ['ETHANE', 'C2H6'], ['PROPANE', 'C3H8'],
+        ['ISOBUTANE', 'C4H10'], ['CO2'], [], ['WATER', 'H2O'], ['OXYGEN', 'O2'], ['NITROGEN', 'N2'],
+        [], [], [], [], ['HYDROGEN', 'H2'], ['DEUTERIUM', 'D2'], [], [], ['DME']
     ]
 
     GridSettings = {
@@ -51,11 +51,13 @@ class OdieRun:
 
     def ListGases(self):
         for idx, gas in enumerate(self.Gases):
-            if type(gas) is str:
+            if gas:
                 print("{} {}".format(idx, gas))
 
     def GasCode(self, GasName):
-        return self.Gases.index(GasName)
+        for idx, names in enumerate(self.Gases):
+            if GasName.upper() in names:
+                return idx
 
     def GasName(self, Code):
         return Gases[Code]
@@ -109,7 +111,8 @@ def ProcessInputs(self, MBObject, Inputs, PrintSettings=False):
         MBObject.Num_Samples = Inputs['NumSamples']
 
         if PrintSettings:
-            print(Inputs)
+            print("Simulation settings...")
+            print(json.dumps(Inputs, indent=4))
 
         return True
 
@@ -138,8 +141,12 @@ def ProcessOutputs(self, MBObject):
             [MBObject.ErrorDiffusionXZ, MBObject.ErrorDiffusionYZ, MBObject.ErrorDiffusionZ]
         ]
         Outputs['DTensor'] = PBRes(np.array(DTensor), np.array(DTensorErr))
+
         Outputs['AttachmentRate'] = PBRes(MBObject.AttachmentRate, MBObject.AttachmentRateError)
         Outputs['IonisationRate'] = PBRes(MBObject.IonisationRate, MBObject.IonisationRateError)
+        Outputs['AttachmentSST'] = PBRes(MBObject.AttachmentSST, MBObject.AttachmentSSTErr)
+        Outputs['IonisationSST'] = PBRes(MBObject.AlphaSST, MBObject.AlphaSSTErr)
+        Outputs['Crossed_SST'] = MBObject.Crossed_SST
 
         lor_angle, lor_error = self.CalcLorentzAngle(MBObject)
         Outputs['LorentzAngle'] = PBRes(lor_angle, lor_error)
@@ -427,9 +434,16 @@ def WriteGasFile(self, FileName, GridOutput):
                     dt = Output['DT1'].val * SQRTP * 1E-4
 
                     #Yes, alpha and alpha0 are supposed to be the same value.
-                    alpha = Output['IonisationRate'].val
+                    alpha  = Output['IonisationRate'].val
                     alpha0 = Output['IonisationRate'].val
-                    eta = Output['AttachmentRate'].val
+                    eta    = Output['AttachmentRate'].val
+
+                    #If the steady state threshold was crossed, use the SST calculation of the ionisation
+                    #and attachment coefficients as the default values.
+                    if Output['Crossed_SST'] is True:
+                        alpha  = Output['IonisationSST'].val
+                        alpha0 = Output['IonisationSST'].val
+                        eta    = Output['AttachmentSST'].val
 
                     #If the coefficients are zero, set to -30 as a sufficiently small power; log(-30) is basically zero.
                     #Otherwise store the logarithm of the reduced coefficients.

PyBoltz/PyBoltzRun.py

@@ -22,7 +22,7 @@ class PyBoltzRun:
     PBSettings   ={'Gases'                 :['NEON','CO2'],
                    'Fractions'             :[90,10],
                    'Max_collisions'        :4e7,
-                   'EField_Vcm'            :100, 
+                   'EField_Vcm'            :100,
                    'Max_electron_energy'   :0,
                    'Temperature_C'         :23,
                    'Pressure_Torr'         :750.062,
@@ -37,29 +37,35 @@ class PyBoltzRun:
                    'Decor_Step'            :0,
                    'NumSamples'            :10}
     '''Dictionary used to store the inputs/settings for the PyBoltz simulation.'''
+
     # Available Gases
-    Gases = [np.nan, 'CF4', 'ARGON', 'HELIUM4', 'HELIUM3', 'NEON', 'KRYPTON', 'XENON', 'CH4', 'ETHANE', 'PROPANE'
-         , 'ISOBUTANE', 'CO2', np.nan, 'H2O', 'OXYGEN', 'NITROGEN', np.nan, np.nan, np.nan, np.nan
-         , 'HYDROGEN', 'DEUTERIUM', np.nan, np.nan, 'DME']
-    '''Array of gases in PyBoltz.'''
+    Gases = [
+        [], ['CF4'], ['ARGON', 'AR'], ['HELIUM4', 'HE4'], ['HELIUM3', 'HE3'], ['NEON', 'NE'],
+        ['KRYPTON', 'KR'], ['XENON', 'XE'], ['METHANE', 'CH4'], ['ETHANE', 'C2H6'], ['PROPANE', 'C3H8'],
+        ['ISOBUTANE', 'C4H10'], ['CO2'], [], ['WATER', 'H2O'], ['OXYGEN', 'O2'], ['NITROGEN', 'N2'],
+        [], [], [], [], ['HYDROGEN', 'H2'], ['DEUTERIUM', 'D2'], [], [], ['DME']
+    ]
+    '''Array of available gases in PyBoltz.'''
 
     # Print list of available gases
     def ListGases(self):
-        '''Function used to print all the gases names in PyBoltz.'''
-        for g in self.Gases:
-            if(type(g)==str):
-                print(g,self.GasCode(g))
+        '''Function used to print all the gas names in PyBoltz.'''
+        for idx, gas in enumerate(self.Gases):
+            if gas:
+                print("{} {}".format(idx, gas))
 
-    # Convert GasName into MagBoltz GasCode            
-    def GasCode(self,GasName):
-        '''Function used to get the ID of the gas. The ID is simply the index of that gas in that array.'''
-        return self.Gases.index(GasName)
+    # Convert GasName into MagBoltz GasCode
+    def GasCode(self, GasName):
+        '''Function used to get the ID of the gas. The ID is simply the index of that gas in the array.'''
+        for idx, names in enumerate(self.Gases):
+            if GasName.upper() in names:
+                return idx
 
     # Convert MagBoltz GasCode in GasName
-    def GasName(self,Code):
-        '''Function used to return the name of the Gas ID given.'''
+    def GasName(self, Code):
+        '''Function used to return the name(s) of the given Gas ID.'''
         return Gases[Code]
-    
+
     # Load Input Dictionary into MagBoltz object
     def ProcessInputs(self,MBObject, Inputs):
         '''Function used to setup the PyBoltz Object with the given inputs in the PBSettings dictionary.'''