all fixed, scaling plots

SSME.py

@@ -97,6 +97,7 @@
     h = altitudes[g] #altitude [m]
     #######################################################################
     #P_atm = 101325*(1 - 2.25577*(10**-5) * h)**5.25588
+
     if h >= 25000:
         T_atm = (-131.21 + 0.00299*h) + 273.14 #[K]
         P_atm = (2.488 * ((T_atm)/ 216.6)**-11.388)*1000 #[Pa]
@@ -112,6 +113,17 @@
     else:
         print("ERROR: OUT OF ALTITUDE RANGE")
 
+    #calculate exit velocity
+    #Nozzle Geometry
+    Area_ratio = 19.8
+    A_throat = 0.0599 #[m]
+    A_exit = A_throat*Area_ratio #[m]
+
+    F = 2188080 #[N] thrust at 104.5% RPL
+    P_e = 13798.5 #[Pa] exit pressure of nozzle
+    u_e = (F-((P_e-ambient_P)*A_exit))/(mdot_f+mdot_ox)
+    print(u_e)
+
     #Calculate composition of O2, and N2
     X_N2 = 0.78084
     X_O2 = 0.2095
@@ -177,11 +189,8 @@
     comp_Noz1 = state.X[n]
     mdot_Noz = mdot_f+mdot_ox
 
-    #Nozzle Geometry
-    L_Noz = 0.01 #[m]
-    A_throat = 0.0599 #[m]
-    A_exit = A_throat*35 #[m]
-
+    L_Noz = (math.sqrt(1/3.1415 *math.sqrt(A_throat*Area_ratio)) - (1/3.1415 *math.sqrt(A_throat)))**2 #[m]
+    print(L_Noz)
     #Call nozzle function
     Noz_states = nozzle(T_Noz1, P_Noz1, comp_Noz1, A_throat, A_exit, L_Noz, mdot_ox, mdot_f)
     n = len(Noz_states.T)-1
@@ -190,7 +199,7 @@
     results_T[1] = Noz_states.T[n]
     results_P[1] = Noz_states.P[n]
     results_X[1,:] = Noz_states.X[n]
-    results_u[1] = (mdot_f+mdot_ox)/(Noz_states.density[n]*A_exit)
+    results_u[1] = u_e #[m/s] from F = m_dot*u_e (rocket equation)
 
     #SAVE NEW STATES TO YAML
     #Load YAML file to append new data
@@ -212,11 +221,12 @@
 
     #Area function
     A_throat = 0.0599
+    r_in = math.sqrt(A_throat)/3.1415
     dAdx = np.zeros(len(Noz_states.x)-1)
     A = np.zeros(len(Noz_states.x)-1)
     t = 0
     for t in range(n):
-        A[t] = 3.1415*(np.sqrt(Noz_states.x[t])+A_throat)**2
+        A[t] = 3.1415*(np.sqrt(t)+r_in)**2
 
     #############################################################################
     ### SHOCKS/EXPANSION ###
@@ -228,7 +238,7 @@
     u = results_u[1] #velocity at exit of nozzle (m/s)
     gamma = 1.1
     print(u)
-    P1 = 13789.5 #Noz_states.P[n] PLACEHOLDER, need design exit pressure or design alt
+    P1 = P_e #Noz_states.P[n] PLACEHOLDER, need design exit pressure or design alt
     T1 = Noz_states.T[n]
     M1 = u/math.sqrt(gamma*P1/Noz_states.density[n])
     P2 = P_atm #Pa

UpstreamPlotting.ipynb

@@ -68,6 +68,7 @@
     }
    ],
    "source": [
+    "#COMBUSTION CHAMBER\n",
     "#PLOT TEMPERATURE\n",
     "T_lit = np.full((len(state.T),1),3588.7)\n",
     "OxInj_T = np.full((len(state.T),1),145.9)\n",
@@ -226,7 +227,18 @@
    "id": "9969d84b",
    "metadata": {},
    "outputs": [],
-   "source": []
+   "source": [
+    "#NOZZLE\n",
+    "#velocity vs isentropic\n",
+    "\n",
+    "#species vs isentropic\n",
+    "\n",
+    "#pressure vs isentropic\n",
+    "\n",
+    "#temp vs isentropic\n",
+    "\n",
+    "#"
+   ]
   }
  ],
  "metadata": {

nozzle.py

@@ -30,12 +30,13 @@ def __call__(self, t, y):
         nsp = 53 #number of species in mechanism
 
         if t == 0:
-            dAdx = 0.1888181/math.sqrt(t+0.001) + 3.1415
+            dAdx = (3.1415*math.sqrt(t+0.001) + 0.244745)/math.sqrt(t+0.001)
         else:
-            dAdx = 0.1888181/math.sqrt(t) + 3.1415
+            dAdx = (3.1415*math.sqrt(t) + 0.244745)/math.sqrt(t)
 
         mdot = 67.35 + 404.79
         A_in = 0.0599
+        r_in = math.sqrt(A_in)/3.1415
 
         # State vector is [T, Y_1, Y_2, ... Y_K]
         self.gas.TDY = y[1], y[0], y[2:nsp+2]
@@ -46,7 +47,7 @@ def __call__(self, t, y):
 
         #converging
         #create new function to find dAdx and A etc
-        A = 3.1415*(np.sqrt(t)+A_in)**2
+        A = 3.1415*(np.sqrt(t)+r_in)**2
 
         MW_mix = self.gas.mean_molecular_weight
         Ru = ct.gas_constant
@@ -97,7 +98,7 @@ def nozzle(T_Noz1, P_Noz1, comp_Noz1, A_throat, A_exit, L_Noz, mdot_ox, mdot_f):
     P_throat = P_t*(2*gamma-1)**(-gamma/(gamma-1))
     T_throat = T_t*(1/(2*gamma-1))
 
-    #print("\nMACH NUMBER AT NOZZLE INLET: ", M_CC)
+    print("\nMACH NUMBER AT NOZZLE INLET: ", M_CC)
     #print("\nTHROAT TEMPERATURE: ", T_throat)
     #print("\nTHROAT PRESSURE: ", P_throat)