sets up test config and get analytical concentrations

src/box_model.py

@@ -285,32 +285,31 @@ def solve(self, path_to_output = None):
             next_conditions = self.evolving_conditions.conditions[0]
             next_conditions_time = self.evolving_conditions.times[0]
 
-        #initializes file headers if output file is present
+        #initalizes output headers
         output_array = []
-        if(path_to_output != None):
-            headers = []
-            headers.append("time")
-            headers.append("ENV.temperature")
-            headers.append("ENV.pressure")
-            for spec in self.species_list.species:
-                headers.append("CONC." + spec.name)
-            
-            output_array.append(headers)
+
+        headers = []
+        headers.append("time")
+        headers.append("ENV.temperature")
+        headers.append("ENV.pressure")
+        for spec in self.species_list.species:
+            headers.append("CONC." + spec.name)
+
+        output_array.append(headers)
 
         #runs the simulation at each timestep
+
         curr_time = 0
         while(curr_time <= self.box_model_options.simulation_length):
 
-            #appends row to output if file is present
-            if(path_to_output != None):
-                row = []
-                row.append(curr_time)
-                row.append(curr_conditions.temperature)
-                row.append(curr_conditions.pressure)
-                for conc in curr_concentrations:
-                    row.append(conc)
-                output_array.append(row)
-
+            row = []
+            row.append(curr_time)
+            row.append(curr_conditions.temperature)
+            row.append(curr_conditions.pressure)
+            for conc in curr_concentrations:
+                row.append(conc)
+            output_array.append(row)
+
             #iterates evolvings conditons if enough time has elapsed
             if(next_conditions != None and next_conditions_time <= curr_time):
                 curr_conditions = next_conditions
@@ -330,9 +329,7 @@ def solve(self, path_to_output = None):
             #solves and updates concentration values in concentration array
             musica.micm_solve(self.solver, self.box_model_options.chem_step_time, curr_conditions.temperature, curr_conditions.pressure, curr_concentrations)
 
-
-
-
+
             #increments time
             curr_time += self.box_model_options.chem_step_time  
 
@@ -342,6 +339,9 @@ def solve(self, path_to_output = None):
                 writer = csv.writer(output)
                 writer.writerows(output_array)
 
+        #returns output_array 
+        return output_array
+
     def readFromUIJson(self, path_to_json):
         """
         TODO: Read the box model configuration from json and sets config

src/configs/test_config_1/camp_data/reactions.json

@@ -1 +1 @@
-{"camp-data": [{"type": "MECHANISM", "name": "music box interactive configuration", "reactions": [{"type": "ARRHENIUS", "A": 0.00012, "Ea": 1.0354868e-21, "B": 7, "D": 50, "E": 0.5, "reactants": {"B": {"qty": 1}}, "products": {"C": {"yield": 1}, "irr__089f1f45-4cd8-4278-83d5-d638e98e4315": {"yield": 1}}}, {"type": "ARRHENIUS", "A": 0.004, "Ea": -6.903245e-22, "B": 0, "D": 300, "E": 0, "reactants": {"A": {"qty": 1}}, "products": {"B": {"yield": 1}, "irr__2a109b21-bb24-41ae-8f06-7485fd36f1a7": {"yield": 1}}}]}]}
+{"camp-data": [{"type": "MECHANISM", "name": "music box interactive configuration", "reactions": [{"type": "ARRHENIUS", "A": 0.00012, "B": 7, "C" : 75, "D": 50, "E": 0.5, "reactants": {"B": {"qty": 1}}, "products": {"C": {"yield": 1}, "irr__089f1f45-4cd8-4278-83d5-d638e98e4315": {"yield": 1}}}, {"type": "ARRHENIUS", "A": 0.004, "C" : 50, "reactants": {"A": {"qty": 1}}, "products": {"B": {"yield": 1}, "irr__2a109b21-bb24-41ae-8f06-7485fd36f1a7": {"yield": 1}}}]}]}

src/python_test_1.py

@@ -1,6 +1,8 @@
 
 from box_model import BoxModel
 
+import math
+
 
 def __main__():
     box_model = BoxModel()
@@ -12,9 +14,52 @@ def __main__():
     box_model.readConditionsFromJson(conditions_path)
     box_model.create_solver(camp_path)
 
-    #solves and outputs to file
-    output_path = "test_config_1_out.csv"
-    box_model.solve(path_to_output=output_path)
+    #solves and saves output
+    output = box_model.solve()
+
+    conc_a_index = output[0].index('CONC.A')
+    conc_b_index = output[0].index('CONC.B')
+    conc_c_index = output[0].index('CONC.C')
+
+    #extracts model concentratins from data output
+    model_concentrations = [[row[conc_a_index], row[conc_b_index], row[conc_c_index]] for row in output[1:]]
+
+
+
+
+    #initalizes concentrations 
+    analytical_concentrations = []
+    analytical_concentrations.append([1,0,0])
+
+    time_step = box_model.box_model_options.chem_step_time
+    sim_length = box_model.box_model_options.simulation_length
+
+    temperature = box_model.initial_conditions.temperature
+    pressure = box_model.initial_conditions.pressure
+
+    k1 = 4.0e-3 * math.exp(50 / temperature)
+    k2 = 1.2e-4 * math.exp(75 / temperature) * (temperature / 50)**7* (1.0 + 0.5 * pressure)
+
+
+    #gets analytical concentrations 
+    curr_time = time_step
+
+    idx_A = 0
+
+    while curr_time < sim_length:
+
+        row = []
+        initial_A = analytical_concentrations[0][idx_A]
+        A_conc = initial_A * math.exp(-(k1)* curr_time)
+        B_conc = initial_A * (k1 / (k2 - k1)) * (math.exp(-k1 * curr_time) - math.exp(-k2 * curr_time))
+        C_conc =  initial_A * (1.0 + (k1 * math.exp(-k2 * curr_time) - k2 * math.exp(-k1 * curr_time)) / (k2 - k1))
+
+        analytical_concentrations.append([A_conc, B_conc, C_conc])
+        curr_time += time_step
+
+
+    #TODO: Assert concentrations
+
 
 if __name__ == "__main__":
     __main__()