mirror of
https://gitlab.com/lander-team/air-prop-simulation.git
synced 2025-07-22 22:21:36 +00:00
119 lines
2.6 KiB
Julia
119 lines
2.6 KiB
Julia
using Unitful
|
||
using DataFrames
|
||
using Plots
|
||
theme(:ggplot2);
|
||
using CSV
|
||
using Measurements
|
||
using Measurements: value, uncertainty
|
||
using Printf
|
||
|
||
# Tank https://www.amazon.com/Empire-Paintball-BASICS-Pressure-Compressed/dp/B07B6M48SR/
|
||
V = (85 ± 5)u"inch^3"
|
||
P0 = (4200.0 ± 300)u"psi"
|
||
Wtank = (2.3 ± 0.2)u"lb"
|
||
Pmax = (250 ± 50)u"psi" # Max Pressure that can come out the nozzle
|
||
|
||
Wsolenoid = 1.5u"kg"
|
||
|
||
# Params
|
||
d_nozzle = ((1 // 18) ± 0.001)u"inch"
|
||
a_nozzle = (pi / 4) * d_nozzle^2
|
||
|
||
# Universal Stuff
|
||
P_amb = (1 ± 0.2)u"atm"
|
||
γ = 1.4 ± 0.05
|
||
R = 287.05u"J/(kg * K)"
|
||
T = (300 ± 20)u"K"
|
||
|
||
t = 0.0u"s"
|
||
P = P0 |> u"Pa"
|
||
M = V * (P / (R * T)) |> u"kg"
|
||
ts = 1u"ms"
|
||
df = DataFrame(Thrust = (0 ± 0)u"N", Pressure = P, Time = t, Mass = M)
|
||
|
||
while M > 0.005u"kg"
|
||
# Calculate what is leaving tank
|
||
P = minimum([P, Pmax])
|
||
ve = sqrt((2 * γ / (γ - 1)) * R * T * (1 - P_amb / P)^((γ - 1) / γ)) |> u"m/s"
|
||
ρ = P / (R * T) |> u"kg/m^3"
|
||
ṁ = ρ * a_nozzle * ve |> u"kg/s"
|
||
|
||
Thrust = ṁ * ve + a_nozzle * (P - P_amb) |> u"N"
|
||
|
||
M = M - ṁ * ts |> u"kg"
|
||
P = (M * R * T) / V |> u"Pa"
|
||
t = t + ts
|
||
|
||
df_step = DataFrame(Thrust = Thrust, Pressure = P, Time = t, Mass = M)
|
||
append!(df, df_step)
|
||
end
|
||
final_time = t |> u"s"
|
||
impulse = sum(df.Thrust) * ts |> u"N*s"
|
||
println("---------------------------\n\n\n\n")
|
||
println("Burn Time: ", t)
|
||
|
||
println("Mass Total: ", Wtank + Wsolenoid + maximum(df.Mass))
|
||
|
||
print(describe(df))
|
||
|
||
thrust_values = df.Thrust .|> ustrip .|> value
|
||
thrust_uncertainties = df.Thrust .|> ustrip .|> uncertainty
|
||
|
||
plot(
|
||
df.Time .|> ustrip,
|
||
thrust_values,
|
||
title = "Thrust Over Time",
|
||
ribbon = (thrust_uncertainties, thrust_uncertainties),
|
||
fillalpha = 0.2,
|
||
label = "Thrust",
|
||
xlabel = "Time (s)",
|
||
ylabel = "Thrust (N)",
|
||
)
|
||
|
||
out = DataFrame(
|
||
Thrust = thrust_values,
|
||
Uncertainty = thrust_uncertainties,
|
||
Time = df.Time .|> u"s" .|> ustrip,
|
||
)
|
||
CSV.write("AirProp.csv", out)
|
||
|
||
### Save data to readme.md
|
||
|
||
savefig("ThrustCurve.png")
|
||
|
||
include("prop_compare.jl")
|
||
|
||
b = IOBuffer();
|
||
t = TextDisplay(b);
|
||
display(t, "text/html", describe(df));
|
||
table = String(take!(b)); # https://stackoverflow.com/a/60443621/8774114
|
||
|
||
t = TextDisplay(b);
|
||
display(t, "text/html", describe(df));
|
||
table = String(take!(b)); # https://stackoverflow.com/a/60443621/8774114
|
||
|
||
readme = Printf.@sprintf """
|
||
# Air Prop Thrust
|
||
|
||
`readme auto generated, run air_prop_sim.jl to recreate.`
|
||
|
||
|
||
|
||
Total Impulse: %s
|
||
|
||
Burn TIme: %s
|
||
|
||
More Data:
|
||
|
||
%s
|
||
|
||
## Comparison to Rocket Motors
|
||
|
||
|
||
|
||
https://www.thrustcurve.org/
|
||
|
||
""" impulse final_time table
|
||
|
||
write("readme.md", readme);
|