1
0
mirror of https://gitlab.com/Anson-Projects/projects.git synced 2025-06-16 06:56:46 +00:00
Projects/_posts/2021-04-01-air-propulsion-simulation/air-propulsion-simulation.Rmd
2021-04-01 22:53:44 -07:00

119 lines
2.9 KiB
Plaintext
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.

---
title: "Air Propulsion Simulation"
description: |
Simulating the performace of an air propulsion system as an alternative to solid rocket motors.
author:
- name: Anson Biggs
url: https://ansonbiggs.com
repository_url: https://gitlab.com/lander-team/air-prop-simulation
date: 04-01-2021
fig_width: 6
fig_align: "center"
output:
distill::distill_article:
self_contained: false
categories:
- Julia
- Capstone
---
Boilerplate intro about why all of this was done
```{r setup, include=FALSE}
library(ggplot2)
knitr::opts_chunk$set(echo = TRUE, results = 'hide')
library(JuliaCall)
julia_setup(JULIA_HOME = "/opt/julia-1.6.0/bin/")
```
```{julia, code_folding=TRUE}
using Plots
plotly()
theme(:ggplot2); # In true R spirit
using Unitful
using DataFrames
using Measurements
using Measurements: value, uncertainty
```
This code is just the setup the setup, using values scraped from various parts of the world wide web.
```{julia}
# 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"
```
This is the actual simulation. Maybe throw some references in and explain some equations.
The rocket equation is pretty sick:
$$T = \dot{m} \cdot v_\text{Exit} + A_\text{Nozzle} \cdot (P - P_\text{Ambient}) $$
And thats about all you need to get to the moon.
```{julia}
let
t = 0.0u"s"
P = P0 |> u"Pa"
M = V * (P / (R * T)) |> u"kg"
ts = 1u"ms"
global df = DataFrame(Thrust=(0 ± 0)u"N", Pressure=P0, Time=0.0u"s", 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
end
```
Heres the results plotted. Notice the massive error once the tank starts running low.
```{julia, echo=FALSE, results='show'}
thrust_values = df.Thrust .|> ustrip .|> value;
thrust_uncertainties = df.Thrust .|> ustrip .|> uncertainty;
out = DataFrame(Thrust=thrust_values, Uncertainty=thrust_uncertainties, Time=df.Time .|> u"s" .|> ustrip);
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)",
)
```
Big conclusion about things.