1
0
mirror of https://gitlab.com/lander-team/air-prop-simulation.git synced 2025-07-23 22:51:35 +00:00

added monte carlo simulation

This commit is contained in:
2021-03-20 15:15:37 -07:00
parent bc86dede11
commit 6c5f6b8459

View File

@@ -1,37 +1,42 @@
using Unitful using Unitful
using DataFrames using DataFrames
using Plots using Plots
# using Gadfly
using UnitfulRecipes using UnitfulRecipes
using Roots using Roots
using CSV using CSV
using MonteCarloMeasurements
unsafe_comparisons(true)
# Tank https://www.amazon.com/Empire-Paintball-BASICS-Pressure-Compressed/dp/B07B6M48SR/ # Tank https://www.amazon.com/Empire-Paintball-BASICS-Pressure-Compressed/dp/B07B6M48SR/
V = 90u"inch^3" V = (85 ± 5)u"inch^3"
P0 = 4500.0u"psi" P0 = (4200.0 ± 300)u"psi"
Wtank = 2.3u"lb" Wtank = (2.3 ± 0.2)u"lb"
Pmax = 100u"psi" Pmax = (100 ± 5)u"psi"
Wsolenoid = 1.5u"kg" Wsolenoid = 1.5u"kg"
# Params # Params
d_nozzle = (1 // 16) * u"inch" d_nozzle = ((1 // 16) ± 0.001)u"inch"
a_nozzle = (pi / 4) * d_nozzle^2 a_nozzle = (pi / 4) * d_nozzle^2
Poutmax = 800u"psi" Poutmax = (800 ± 50)u"psi"
# Universal Stuff # Universal Stuff
P_amb = 1u"atm" P_amb = (1 ± 0.2)u"atm"
γ = 1.4 γ = 1.4 ± 0.05
R = 287.05u"J/(kg * K)" R = 287.05u"J/(kg * K)"
T = 300u"K" T = (300 ± 20)u"K"
# Setting up while loop # Setting up while loop
t = 1.0u"s" t = 0.0u"s"
P = P0 |> u"Pa" P = P0 |> u"Pa"
M = V * (P / (R * T)) |> u"kg" M = V * (P / (R * T)) |> u"kg"
df = DataFrame(Thrust = 0.0u"N", Pressure = P, Time = t, Mass = M) # df = DataFrame(Thrust = (0 ± 0)u"N", Pressure = P, Time = t, Mass = M)
ts = 0.001u"s" df = DataFrame()
while M > 0.01u"kg" ts = 10u"ms"
while M > 0.1u"kg"
# while t < 30u"s" # while t < 30u"s"
# Calculate what is leaving tank # Calculate what is leaving tank
P = minimum([P, Pmax]) P = minimum([P, Pmax])
@@ -44,16 +49,25 @@ while M > 0.01u"kg"
# Calculate what is still in the tank # Calculate what is still in the tank
M = M - * ts |> u"kg" M = M - * ts |> u"kg"
P = (M * R * T) / V |> u"Pa" P = (M * R * T) / V |> u"Pa"
t = t + ts # df_step = DataFrame(Thrust = Thrust, Pressure = P, Time = t, Mass = M)
df_step = DataFrame(Thrust = Thrust, Pressure = P, Time = t, Mass = M) if df == DataFrame()
append!(df, df_step) df = DataFrame(Thrust = Thrust, Pressure = P, Time = t, Mass = M)
else
append!(df, DataFrame(Thrust = Thrust, Pressure = P, Time = t, Mass = M))
end
t = t + ts
end end
impulse = sum(df.Thrust) * ts |> u"N*s"
println("---------------------------\n\n\n\n") println("---------------------------\n\n\n\n")
println("Total Impulse: ", sum(df.Thrust) * ts) println("Total Impulse: ", impulse)
println("Burn Time: ", t) println("Burn Time: ", t)
plot(impulse)
println("Mass Total: ", Wtank + Wsolenoid + maximum(df.Mass)) println("Mass Total: ", Wtank + Wsolenoid + maximum(df.Mass))
print(describe(df)) # plot(sum(df.Thrust) * ts |> u"N*s")
plot(df.Time, df.Thrust, title = "Thrust Over Time") # print(describe(df))
CSV.write("thrustdata.csv", df .|> ustrip) # df.Time = round(u"s", df.Time)
# plot(df.Time, df.Thrust, title = "Thrust Over Time")
# CSV.write("thrustdata.csv", df .|> ustrip)
# plot(df.Time .|> u"s" .|> ustrip, df.Thrust .|> u"N" .|> ustrip)> ustrip, df.Thrust .|> u"N" .|> ustrip)