1
0
mirror of https://gitlab.com/lander-team/air-prop-simulation.git synced 2025-07-23 14:41:38 +00:00
Files
Air-Prop-Simulation/thrust.jl
2021-03-22 11:39:26 -07:00

106 lines
2.3 KiB
Julia
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

using Unitful
using DataFrames
using Plots
theme(:ggplot2);
using UnitfulRecipes
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 // 16) ± 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"
# Setting up while loop
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"
# while t < 30u"s"
# 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"
# Calculate what is still in the tank
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("Total Impulse: ", impulse)
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=.2,label="Thrust",
xlabel="Time (s)",
ylabel="Thrust (N)",
palette=:leonardo,
)
### Save data to readme.md
savefig("ThrustCurve.png")
b = IOBuffer();
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
Total Impulse: %s
Burn TIme: %s
![Thrust Over Time](ThrustCurve.png)
More Data:
%s
`readme auto generated`
""" impulse final_time table
println(readme)
write("readme.md", readme)