Skip to content

Commit

Permalink
Slurm set up, jupyter now outdated, merge to main
Browse files Browse the repository at this point in the history
  • Loading branch information
@chinahg
chinahg committed Oct 20, 2022
1 parent 21b4d80 commit 80b63f6
Show file tree
Hide file tree
Showing 9 changed files with 2,634 additions and 459 deletions.
1 change: 1 addition & 0 deletions .gitignore
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,3 +22,4 @@ docs/site/
# committed for packages, but should be committed for applications that require a static
# environment.
Manifest.toml
.out
17 changes: 17 additions & 0 deletions Project.toml
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,17 @@
[deps]
CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab"
HDF5 = "f67ccb44-e63f-5c2f-98bd-6dc0ccc4ba2f"
Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59"
NBInclude = "0db19996-df87-5ea3-a455-e3a50d440464"
OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
PlotlyJS = "f0f68f2c-4968-5e81-91da-67840de0976a"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
PyCall = "438e738f-606a-5dbb-bf0a-cddfbfd45ab0"
PyPlot = "d330b81b-6aea-500a-939a-2ce795aea3ee"
Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
Revise = "295af30f-e4ad-537b-8983-00126c2a3abe"
Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
StructArrays = "09ab397b-f2b6-538f-b94a-2f83cf4a842a"
Trapz = "592b5752-818d-11e9-1e9a-2b8ca4a44cd1"
YAML = "ddb6d928-2868-570f-bddf-ab3f9cf99eb6"
238 changes: 177 additions & 61 deletions plumeSSMEJP.ipynb → UpdatedplumeSSMEJP.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -7,18 +7,18 @@
"metadata": {},
"outputs": [],
"source": [
"#import Pkg\n",
"#\n",
"#Pkg.add(\"Interpolations\")\n",
"#Pkg.add(\"PyPlot\")\n",
"#Pkg.add(\"PyCall\")\n",
"#Pkg.add(\"OrdinaryDiffEq\")\n",
"#Pkg.add(\"YAML\")\n",
"#Pkg.add(\"DelimitedFiles\")\n",
"#Pkg.add(\"CSV\")\n",
"#Pkg.add(\"HDF5\")\n",
"#Pkg.add(\"StructArrays\")\n",
"#Pkg.add(\"Random\")"
"import Pkg\n",
"\n",
"Pkg.add(\"Interpolations\")\n",
"Pkg.add(\"PyPlot\")\n",
"Pkg.add(\"PyCall\")\n",
"Pkg.add(\"OrdinaryDiffEq\")\n",
"Pkg.add(\"YAML\")\n",
"Pkg.add(\"DelimitedFiles\")\n",
"Pkg.add(\"CSV\")\n",
"Pkg.add(\"HDF5\")\n",
"Pkg.add(\"StructArrays\")\n",
"Pkg.add(\"Random\")\n"
]
},
{
Expand All @@ -41,7 +41,7 @@
{
"cell_type": "code",
"execution_count": null,
"id": "0b92c61f",
"id": "aba7c4a1",
"metadata": {},
"outputs": [],
"source": [
Expand Down Expand Up @@ -173,7 +173,8 @@
"end\n",
"\n",
"\"\"\"\n",
" compute_y(u, T, Δψ, R, p)\n",
" compute_y\n",
" (u, T, Δψ, R, p)\n",
"Computes y locations corresponding to particular states, based on the \n",
"grid spacing in ϕ-ψ space. \n",
"# Inputs \n",
Expand Down Expand Up @@ -345,7 +346,7 @@
" end\n",
"\n",
" heatmap((x[2:x_lim]), y, clim=(0, Inf),\n",
" var[1:y_len, 2:x_lim], colorbar_title=clabel, size=(700, 500), dpi=300, c=colormap)\n",
" var[1:y_len, 2:x_lim], colorbar_title=clabel, size=(700, 250), dpi=300, c=colormap)\n",
" xlabel!(xlabel)\n",
" ylabel!(ylabel)\n",
" #xticks!([-2, -1, 0, 1, 2, 3], [\"0.1\", \"1\", \"10\", \"100\",])\n",
Expand Down Expand Up @@ -506,8 +507,8 @@
"h = [16000] #Int.(LinRange(16000, 20000, space))\n",
"g = 1\n",
"\n",
"n = 200 #90 #n steps in x dir\n",
"s = 900 #80 #s steps in y direction\n",
"n = 1000 #90 #n steps in x dir\n",
"s = 2000 #80 #s steps in y direction\n",
"T_save = zeros(length(h), s, n)\n",
"u_save = zeros(length(h), s, n)\n",
"χ_save = zeros(length(h), s, n, n_species)\n",
Expand Down Expand Up @@ -675,8 +676,8 @@
"# Define altitude range to plot over\n",
"h = [16000] #Int.(LinRange(16000, 20000, space)) #make an array and put in loop for all alts\n",
"h_string_arr = string(h)\n",
"\n",
"for m = 1:length(h) #For all specified altitudes\n",
"m = 1\n",
"#for m = 1:length(h) #For all specified altitudes\n",
"\n",
" # Initialize arrays to save results\n",
" T = T_save[m, :, :] #temperature\n",
Expand All @@ -690,19 +691,19 @@
"\n",
" #REGRID SOLUTION\n",
" xx, yy, u_g, T_g, χ_gO2 = regrid_solution(x, y, u, T, χ[:,:,4], 0.01)\n",
" #xx, yy, u_g, T_g, χ_gN2 = regrid_solution(x, y, u, T, χ[:,:,48], 0.01)\n",
" xx, yy, u_g, T_g, χ_gN2 = regrid_solution(x, y, u, T, χ[:,:,48], 0.01)\n",
" #xx, yy, u_g, T_g, χ_gNO2 = regrid_solution(x, y, u, T, χ[:, :, 37], 0.01)\n",
" #xx, yy, u_g, T_g, χ_gN2O = regrid_solution(x, y, u, T, χ[:, :, 38], 0.01)\n",
" #xx, yy, u_g, T_g, χ_gNO = regrid_solution(x, y, u, T, χ[:, :, 36], 0.01)\n",
" #xx, yy, u_g, T_g, χ_gAr = regrid_solution(x, y, u, T, χ[:, :, 49], 0.01)\n",
" xx, yy, u_g, T_g, χ_gAr = regrid_solution(x, y, u, T, χ[:, :, 49], 0.01)\n",
"\n",
" #PLOT 2D MAPS AND SAVE\n",
" x_max = 20 #[m]\n",
"\n",
" #CALCULATE NO EI FOR ALTITUDE\n",
" # Sum of all NOx species as a function of x\n",
" n = 200\n",
" s = 900\n",
" n = 1000\n",
" s = 2000\n",
" sumNO = zeros(n)\n",
" sumNO2 = zeros(n)\n",
" sumN2O = zeros(n)\n",
Expand Down Expand Up @@ -736,6 +737,22 @@
" mdot_tot = zeros(n)\n",
" ppm_tot = zeros(n)\n",
"\n",
" #IN XY\n",
" r = length(yy)\n",
" k = length(xx)\n",
" mdot_Ar_sum_g = zeros(k)\n",
" mdot_N2_sum_g = zeros(k)\n",
" mdot_O2_sum_g = zeros(k)\n",
" mdot_Ar_g = zeros(r,k)\n",
" mdot_N2_g = zeros(r,k)\n",
" mdot_O2_g = zeros(r,k)\n",
" Xarea_g = zeros(r,k)\n",
" rho_tot_g = zeros(r,k)\n",
" MFar_g = zeros(r,k)\n",
" MFn2_g = zeros(r,k)\n",
" MFo2_g = zeros(r,k)\n",
" mdot_tot_g = zeros(k)\n",
"\n",
" Δψ = 0.001 * ones(s) #s vertical grid points\n",
" for i = 2:s\n",
" Δψ[i] = 1.01 * Δψ[i-1] #enlarge with each step by 1.1\n",
Expand Down Expand Up @@ -782,14 +799,41 @@
" mdot_N2[s,:] = mdot_N2[s-1,:]\n",
" mdot_O2[s,:] = mdot_O2[s-1,:]\n",
"\n",
"\n",
" # IN XY COORDINATES\n",
" for i = 1:length(xx), j = 1:(length(yy)-1)\n",
" rho_tot_g[j,i] = P_atm/(R*T_g[j,i]) #[kg/m^3] #need to convert to xx yy coordinates\n",
" MFar_g[j,i] = (χ_gAr[j,i]*(1/28.97)*(39.9))/1000000 #[kg/kg tot]\n",
" MFn2_g[j,i] = (χ_gN2[j,i]*(1/28.97)*(28))/1000000 #[kg/kg tot]\n",
" MFo2_g[j,i] = (χ_gO2[j,i]*(1/28.97)*(32))/1000000 #[kg/kg tot]\n",
" Xarea_g[j,i] = pi*(yy[j+1]^2 - yy[j]^2)\n",
"\n",
" mdot_Ar_g[j,i] = Xarea_g[j,i]*MFar_g[j,i]*rho_tot_g[j,i]*u_g[j,i]\n",
" mdot_N2_g[j,i] = Xarea_g[j,i]*MFn2_g[j,i]*rho_tot_g[j,i]*u_g[j,i]\n",
" mdot_O2_g[j,i] = Xarea_g[j,i]*MFo2_g[j,i]*rho_tot_g[j,i]*u_g[j,i]\n",
"\n",
" #TODO: populate the last entry of sum arrays (stops at s-1)\n",
" mdot_O2_sum_g[i] += mdot_O2_g[j,i]\n",
" mdot_N2_sum_g[i] += mdot_N2_g[j,i]\n",
" mdot_Ar_sum_g[i] += mdot_Ar_g[j,i]\n",
" \n",
" mdot_tot_g[i] += mdot_Ar_g[j,i] + mdot_N2_g[j,i] + mdot_O2_g[j,i]\n",
" end\n",
"\n",
" mdot_Ar_g[length(yy),:] = mdot_Ar_g[length(yy)-1,:]\n",
" mdot_N2_g[length(yy),:] = mdot_N2_g[length(yy)-1,:]\n",
" mdot_O2_g[length(yy),:] = mdot_O2_g[length(yy)-1,:]\n",
"\n",
"\n",
"#PLOTTING\n",
" fig,axc = plt.subplots()\n",
" axc.plot(ϕ, u[s,:], label=\"300/300\")\n",
" axc.plot(ϕ, u[s-1,:], label=\"299/300\")\n",
" axc.plot(ϕ, u[s-50,:], label=\"250/300\")\n",
" axc.plot(ϕ, u[s-50,:], label=\"200/300\")\n",
" axc.plot(ϕ, u[s-100,:], label=\"200/300\")\n",
" axc.plot(ϕ, u[s-150,:], label=\"150/300\")\n",
" axc.plot(ϕ, u[s-200,:], label=\"100/300\")\n",
" axc.plot(ϕ, u[s-250,:], label=\"50/300\")\n",
" axc.plot(ϕ, u[s-200,:], label=\"50/300\")\n",
" axc.plot(ϕ, u[s-299,:], label=\"1/300\")\n",
"#\n",
" axc.set_xlabel(\"ϕ\")\n",
Expand Down Expand Up @@ -819,69 +863,69 @@
" savefig(\"/home/chinahg/GCresearch/rocketemissions/rockettests/\" * h_string * \"m/ring1.png\")\n",
"\n",
" fig,axT = plt.subplots()\n",
" im = axT.imshow(T[:,:,m], cmap=\"summer\")\n",
" im = axT.imshow(T[:,1:50,m], cmap=\"summer\")\n",
" axT.set_ylabel(\"ψ\")\n",
" axT.set_xlabel(\"ϕ\")\n",
" axT.set_title(\"Temperature\")\n",
" cbar = fig.colorbar(im,ax=axT)\n",
" savefig(\"/home/chinahg/GCresearch/rocketemissions/rockettests/\" * h_string * \"m/T.png\")\n",
"\n",
" fig,axu = plt.subplots()\n",
" im = axu.imshow(u[:,:,m], cmap=\"viridis\")\n",
" im = axu.imshow(u[:,1:50,m], cmap=\"viridis\")\n",
" axu.set_ylabel(\"ψ\")\n",
" axu.set_xlabel(\"ϕ\")\n",
" cbar = fig.colorbar(im,ax=axu)\n",
" axu.set_title(\"Velocity\")\n",
" savefig(\"/home/chinahg/GCresearch/rocketemissions/rockettests/\" * h_string * \"m/u.png\")\n",
"\n",
"\n",
" #fig,axs = plt.subplots(2,2)\n",
" #axs[1,1].plot(ϕ, mdot_tot)\n",
" #axs[1,1].set_xlabel(\"ϕ\")\n",
" #axs[1,1].set_ylabel(\"Total mass flow rate [kg/s]\")\n",
" #axs[1,1].set_title(\"Total mdot\")\n",
" #axs[1,1].set_xscale(\"log\")\n",
"#\n",
" #axs[1,2].plot(ϕ, mdot_N2_sumTrunc)\n",
" #axs[1,2].set_xlabel(\"ϕ\")\n",
" #axs[1,2].set_ylabel(\"N2 mass flow rate [kg/s]\")\n",
" #axs[1,2].set_title(\"N2\")\n",
" #axs[1,2].set_xscale(\"log\")\n",
"#\n",
" #axs[2,1].plot(ϕ, mdot_O2_sumTrunc)\n",
" #axs[2,1].set_xlabel(\"ϕ\")\n",
" #axs[2,1].set_ylabel(\"O2 mass flow rate [kg/s]\")\n",
" #axs[2,1].set_xscale(\"log\")\n",
"#\n",
" #axs[2,2].plot(ϕ, mdot_Ar_sumTrunc)\n",
" #axs[2,2].set_xlabel(\"ϕ\")\n",
" #axs[2,2].set_ylabel(\"Ar mass flow rate [kg/s]\")\n",
" #axs[2,2].set_title(\"Ar\")\n",
" #axs[2,2].set_xscale(\"log\")\n",
" #\n",
" #fig.suptitle(\"Truncated Mass Flow Rates, s = 250/300\")\n",
" #fig.tight_layout()\n",
" #savefig(\"/home/chinahg/GCresearch/rocketemissions/rockettests/\" * h_string * \"m/mdot_trunc.png\")\n",
" fig,axs = plt.subplots(2,2)\n",
" axs[1,1].plot(ϕ, mdot_tot)\n",
" axs[1,1].set_xlabel(\"ϕ\")\n",
" axs[1,1].set_ylabel(\"Total mass flow rate [kg/s]\")\n",
" axs[1,1].set_title(\"Total mdot\")\n",
" axs[1,1].set_xscale(\"log\")\n",
"\n",
" axs[1,2].plot(ϕ, mdot_N2_sumTrunc)\n",
" axs[1,2].set_xlabel(\"ϕ\")\n",
" axs[1,2].set_ylabel(\"N2 mass flow rate [kg/s]\")\n",
" axs[1,2].set_title(\"N2\")\n",
" axs[1,2].set_xscale(\"log\")\n",
"\n",
" axs[2,1].plot(ϕ, mdot_O2_sumTrunc)\n",
" axs[2,1].set_xlabel(\"ϕ\")\n",
" axs[2,1].set_ylabel(\"O2 mass flow rate [kg/s]\")\n",
" axs[2,1].set_xscale(\"log\")\n",
"\n",
" axs[2,2].plot(ϕ, mdot_Ar_sumTrunc)\n",
" axs[2,2].set_xlabel(\"ϕ\")\n",
" axs[2,2].set_ylabel(\"Ar mass flow rate [kg/s]\")\n",
" axs[2,2].set_title(\"Ar\")\n",
" axs[2,2].set_xscale(\"log\")\n",
" \n",
" fig.suptitle(\"Truncated Mass Flow Rates, s = 200/300\")\n",
" fig.tight_layout()\n",
" savefig(\"/home/chinahg/GCresearch/rocketemissions/rockettests/\" * h_string * \"m/mdot_trunc.png\")\n",
"\n",
" fig,axX = plt.subplots(2,2)\n",
" axX[1,1].plot(ϕ, mdot_tot)\n",
" axX[1,1].plot(ϕ[1:50], mdot_tot[1:50])\n",
" axX[1,1].set_xlabel(\"ϕ\")\n",
" axX[1,1].set_ylabel(\"Total mass flow rate [kg/s]\")\n",
" axX[1,1].set_title(\"Total mdot\")\n",
" axX[1,1].set_xscale(\"log\")\n",
"\n",
" axX[1,2].plot(ϕ, mdot_N2_sum)\n",
" axX[1,2].plot(ϕ[1:50], mdot_N2_sum[1:50])\n",
" axX[1,2].set_xlabel(\"ϕ\")\n",
" axX[1,2].set_ylabel(\"N2 mass flow rate [kg/s]\")\n",
" axX[1,2].set_title(\"N2\")\n",
" axX[1,2].set_xscale(\"log\")\n",
"\n",
" axX[2,1].plot(ϕ, mdot_O2_sum)\n",
" axX[2,1].plot(ϕ[1:50], mdot_O2_sum[1:50])\n",
" axX[2,1].set_xlabel(\"ϕ\")\n",
" axX[2,1].set_ylabel(\"O2 mass flow rate [kg/s]\")\n",
" axX[2,1].set_xscale(\"log\")\n",
"\n",
" axX[2,2].plot(ϕ, mdot_Ar_sum)\n",
" axX[2,2].plot(ϕ[1:50], mdot_Ar_sum[1:50])\n",
" axX[2,2].set_xlabel(\"ϕ\")\n",
" axX[2,2].set_ylabel(\"Ar mass flow rate [kg/s]\")\n",
" axX[2,2].set_title(\"Ar\")\n",
Expand Down Expand Up @@ -919,19 +963,91 @@
" fig.tight_layout()\n",
" savefig(\"/home/chinahg/GCresearch/rocketemissions/rockettests/\" * h_string * \"m/ppm_tot.png\")\n",
"\n",
" fig,axXg = plt.subplots(2,2)\n",
" custom_xlim = (0,250)\n",
" plt.setp(axXg, xlim=custom_xlim)\n",
" axXg[1,1].plot(xx, mdot_tot_g)\n",
" axXg[1,1].set_xlabel(\"x\")\n",
" axXg[1,1].set_ylabel(\"Total mass flow rate [kg/s]\")\n",
" axXg[1,1].set_title(\"Total mdot\")\n",
" \n",
"\n",
" axXg[1,2].plot(xx, mdot_N2_sum_g)\n",
" axXg[1,2].set_xlabel(\"x\")\n",
" axXg[1,2].set_ylabel(\"N2 mass flow rate [kg/s]\")\n",
" axXg[1,2].set_title(\"N2\")\n",
"\n",
"\n",
" axXg[2,1].plot(xx, mdot_O2_sum_g)\n",
" axXg[2,1].set_xlabel(\"x\")\n",
" axXg[2,1].set_ylabel(\"O2 mass flow rate [kg/s]\")\n",
"\n",
"\n",
" axXg[2,2].plot(xx, mdot_Ar_sum_g)\n",
" axXg[2,2].set_xlabel(\"x\")\n",
" axXg[2,2].set_ylabel(\"Ar mass flow rate [kg/s]\")\n",
" axXg[2,2].set_title(\"Ar\")\n",
"\n",
"\n",
" fig.suptitle(\"Total Mass Flow (all rings)\")\n",
" fig.tight_layout()\n",
" savefig(\"/home/chinahg/GCresearch/rocketemissions/rockettests/\" * h_string * \"m/mdot_tot_g.png\")\n",
"\n",
"\n",
" println(\"Δ mdot Ar total is \", mdot_Ar_sum[n]-mdot_Ar_sum[1])\n",
" println(\"Δ mdot Ar/Δ mdot total is \", (mdot_Ar_sum[n]-mdot_Ar_sum[1])/(mdot_tot[n]-mdot_tot[1]))\n",
" println(mdot_Ar_sum[n], mdot_tot[n])\n",
"\n",
" println(\"plotted \", h_string[m], \" altitude!\\n\")\n",
" #println(mdot_Ar_sum[1])\n",
" #println(mdot_Ar[:,1])\n",
" println(yy[473])\n",
"end\n",
" println(mdot_Ar[:,1])\n",
" println(xx[250])\n",
" #println(yy[167772])\n",
"#end\n",
"\n",
"println(\"done plotting!\")\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "4daf0a88",
"metadata": {},
"outputs": [],
"source": [
"\n",
"close(\"all\")\n",
"fig,axXg = plt.subplots(2,2)\n",
"custom_xlim = (0,250)\n",
"plt.setp(axXg, xlim=custom_xlim)\n",
"axXg[1,1].plot(xx[1:250], mdot_tot_g[1:250])\n",
"axXg[1,1].set_xlabel(\"x\")\n",
"axXg[1,1].set_ylabel(\"Total mass flow rate [kg/s]\")\n",
"axXg[1,1].set_title(\"Total mdot\")\n",
"\n",
"\n",
"axXg[1,2].plot(xx[1:250], mdot_N2_sum_g[1:250])\n",
"axXg[1,2].set_xlabel(\"x\")\n",
"axXg[1,2].set_ylabel(\"N2 mass flow rate [kg/s]\")\n",
"axXg[1,2].set_title(\"N2\")\n",
"\n",
"\n",
"axXg[2,1].plot(xx[1:250], mdot_O2_sum_g[1:250])\n",
"axXg[2,1].set_xlabel(\"x\")\n",
"axXg[2,1].set_ylabel(\"O2 mass flow rate [kg/s]\")\n",
"\n",
"\n",
"axXg[2,2].plot(xx[1:250], mdot_Ar_sum_g[1:250])\n",
"axXg[2,2].set_xlabel(\"x\")\n",
"axXg[2,2].set_ylabel(\"Ar mass flow rate [kg/s]\")\n",
"axXg[2,2].set_title(\"Ar\")\n",
"\n",
"\n",
"fig.suptitle(\"Total Mass Flow (all rings)\")\n",
"fig.tight_layout()\n",
"savefig(\"/home/chinahg/GCresearch/rocketemissions/rockettests/\" * h_string * \"m/mdot_tot_g.png\")"
]
},
{
"cell_type": "code",
"execution_count": null,
Expand Down
Loading

0 comments on commit 80b63f6

Please sign in to comment.