diff --git a/Project.toml b/Project.toml index fa13fe8..2485fb2 100644 --- a/Project.toml +++ b/Project.toml @@ -3,6 +3,10 @@ uuid = "77d2bf28-a3e9-4b9c-9fcf-b85f74cc8a50" authors = ["Alex Fleming and contributors"] version = "1.0.0-DEV" +[deps] +PhysicalConstants = "5ad8b20f-a522-5ce9-bfc9-ddf1d5bda6ab" +Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" + [compat] julia = "1.6" diff --git a/src/ShockwaveProperties.jl b/src/ShockwaveProperties.jl index 6a1bf27..25d4f0f 100644 --- a/src/ShockwaveProperties.jl +++ b/src/ShockwaveProperties.jl @@ -1,5 +1,6 @@ module ShockwaveProperties -# Write your package code here. +include("cpg.jl") +include("normal_shocks.jl") end diff --git a/src/billig.jl b/src/billig.jl new file mode 100644 index 0000000..5aa2edb --- /dev/null +++ b/src/billig.jl @@ -0,0 +1,3 @@ +module BilligShock + +end \ No newline at end of file diff --git a/src/cpg.jl b/src/cpg.jl new file mode 100644 index 0000000..cd78022 --- /dev/null +++ b/src/cpg.jl @@ -0,0 +1,56 @@ +using Unitful + +const _units_cvcp = u"J/kg/K" +const _units_ℳ = u"kg/mol" +const _units_ρ = u"kg/m^3" +const _dimension_ρ = Unitful.dimension(_units_ρ) + +""" +Provides the properties of a calorically perfect gas (or mixture of gases). +""" +struct CaloricallyPerfectGas{T} + c_p::Quantity{T,Unitful.dimension(_units_cvcp),typeof(_units_cvcp)} + c_v::Quantity{T,Unitful.dimension(_units_cvcp),typeof(_units_cvcp)} + ℳ::Quantity{T,Unitful.dimension(_units_ℳ),typeof(_units_ℳ)} + + γ::T + R::Quantity{T,Unitful.dimension(_units_cvcp),typeof(_units_cvcp)} +end + +function CaloricallyPerfectGas(c_p::T, c_v::T, ℳ::T) where {T} + return CaloricallyPerfectGas{T}( + Quantity(c_p, _units_cvcp), + Quantity(c_v, _units_cvcp), + Quantity(ℳ, _units_ℳ), + c_p / c_v, + Quantity(c_p - c_v, _units_cvcp) + ) +end + +const DRY_AIR = CaloricallyPerfectGas(1.0049, 0.7178, 0.0289647) + +enthalpy(gas::CaloricallyPerfectGas, T::Real) = gas.c_p * Quantity(T, u"K") +enthalpy(gas::CaloricallyPerfectGas, T::Quantity{T1,Unitful.𝚯,Units}) where {T1,Units} = gas.c_p * T + +int_energy(gas::CaloricallyPerfectGas, T::Real) = gas.c_v * Quantity(T, u"K") +int_energy(gas::CaloricallyPerfectGas, T::Quantity{T1,Unitful.𝚯,Units}) where {T1,Units} = gas.c_v * T + +function speed_of_sound(gas::CaloricallyPerfectGas, T::Real) + return sqrt(gas.γ * gas.R * Quantity(T, u"K")) +end + +function speed_of_sound(gas::CaloricallyPerfectGas, T::Quantity{T1,Unitful.𝚯,Units}) where {T1,Units} + return sqrt(gas.γ * gas.R * T) +end + +function pressure(gas::CaloricallyPerfectGas, ρ::Real, T::Real) + # calculate ρe then + # calculate p from calorically perfect gas relations + return (gas.γ - 1) * Quantity(ρ, _units_ρ) * int_energy(gas, T) +end + +function pressure(gas::CaloricallyPerfectGas, + ρ::Quanity{U,_dimension_ρ,Units1}, + T::Quantity{U,Unitful.𝚯,Units2}) where {U,Units1,Units2} + return (gas.γ-1) * ρ * T +end diff --git a/src/normal_shocks.jl b/src/normal_shocks.jl new file mode 100644 index 0000000..797c22c --- /dev/null +++ b/src/normal_shocks.jl @@ -0,0 +1,71 @@ +""" +**Equation 4.8** from Anderson&Anderson + +Computes the density change across a shock wave where free stream mach number is ``M`` +and the outward (away from body) normal is ``n̂`` +""" +function shock_density_ratio(M, n̂; gas::CaloricallyPerfectGas=DRY_AIR) + Mn2 = (M * n̂[1])^2 + return ((gas.γ + 1) * Mn2) / ((gas.γ - 1) * Mn2 + 2) +end + +""" +**Equation 4.9** from Anderson&Anderson + +Computes the pressure ratio across a shock wave where the free stream mach number is ``M`` +and the outward (away from body) normal is ``n̂``. +""" +function shock_pressure_ratio(M, n̂; gas::CaloricallyPerfectGas=DRY_AIR) + Mn2 = (M * n̂[1])^2 + return 1 + (2 * gas.γ) / (gas.γ + 1) * (Mn2 - 1) +end + +""" +**Equation 4.10** from Anderson&Anderson + +Computes the normal Mach number ratio across a shock wave where the free stream mach number is ``M`` +and the outward (away from body) normal is ``n̂``. +""" +function shock_normal_mach_ratio(M, n̂; gas::CaloricallyPerfectGas=DRY_AIR) + MnL2 = (M * n̂[1])^2 + Mn2sqr = (MnL2 + (2 / (gas.γ - 1))) / (2 * gas.γ / (gas.γ - 1) * MnL2 - 1) + Mn2 = sqrt(Mn2sqr) + return Mn2 +end + + +""" +Computes the normal velocity ratio across a shock wave where the free stream mach number is ``M`` +and the outward (away from body) normal is ``n̂``. + +Derived from speed of sound proportional to the square root of temperature. + +Useful for computing momentum ratio. +""" +function shock_normal_velocity_ratio(M, n̂; gas::CaloricallyPerfectGas=DRY_AIR) + return sqrt(1 / shock_temperature_ratio(M, n̂; gas=gas)) * shock_normal_mach_ratio(M, n̂; gas=gas) +end + +""" +Computes the tangential mach number ratio across a shock wave, where the free stream mach number is ``M`` +and the outward (away from body) normal is ``n̂``. + +Derived from speed of sound proportional to the square root of temperature. +""" +function shock_tangent_mach_ratio(M, n̂; gas::CaloricallyPerfectGas=DRY_AIR) + return sqrt(shock_temperature_ratio(M, n̂; gas=gas)) +end + +""" +**Equation 4.11** from Anderson&Anderson + +Computes the temperature across a shock wave where the free stream mach number is ``M`` +and the outward (away from body) normal is ``n̂``. +""" +function shock_temperature_ratio(M, n̂; gas::CaloricallyPerfectGas=DRY_AIR) + return shock_pressure_ratio(M, n̂; gas=gas) / shock_density_ratio(M, n̂; gas=gas) +end + +""" +Computes the +""" \ No newline at end of file