diff --git a/docs/src/assets/cryo_tank.svg b/docs/src/assets/cryo_tank.svg index 9c9f5d1c..aa87594f 100644 --- a/docs/src/assets/cryo_tank.svg +++ b/docs/src/assets/cryo_tank.svg @@ -2,9 +2,9 @@ Fuselage skin + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="168.94751" + y="98.494949">Fuselage skin Foam layers + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="187.54964" + y="126.66481">Foam layers Vacuum layer + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="186.02109" + y="139.97253">Vacuum layer Inner vessel + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="115.2411" + y="161.45566">Inner vessel Cryogenic fuel Outer vessel + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="187.24771" + y="156.99146">Outer vessel + Vapor barrier + + + diff --git a/docs/src/assets/doublewalled_tank.svg b/docs/src/assets/doublewalled_tank.svg index addca703..65ba8ebe 100644 --- a/docs/src/assets/doublewalled_tank.svg +++ b/docs/src/assets/doublewalled_tank.svg @@ -2,9 +2,9 @@ Cryogenic fuel + style="font-size:8.46667px;stroke-width:0.5" + x="93.597252" + y="128.4803">Cryogenic fuel Outer vessel + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="171.90758" + y="157.97603">Outer vessel Foam layers + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="121.68092" + y="108.42661">Foam layers Inner vessel + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="108.40814" + y="151.27167">Inner vessel Main stiffening rings + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="69.711891" + y="143.54659">Main stiffening rings Additional stiffening rings + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="93.1464" + y="178.75729">Additional stiffening rings Vacuum layer + style="font-size:5.64444px;fill:#000000;fill-opacity:1;stroke:none;stroke-width:0.5" + x="121.73133" + y="138.56248">Vacuum layer 0`` is a factor to account for the fact that the tank must contain some empty volume to account for ullage. The internal volume of a hemiellipsoidal cap is given by + where ``\rho_{mix} = f_{ull}\rho_g + (1-f_{ull})\rho_l`` is the density of the saturated mixture inside the tank, ``\rho_g`` and ``\rho_l`` are the densities of the fuel in saturated gas and liquid phases, and ``f_{ull}>0`` is a factor to account for the fact that the tank must contain some gas volume for ullage. The internal volume of a hemiellipsoidal cap is given by ```math V_{cap} = \frac{2πR_{t,i}^3}{3AR}. ``` @@ -143,7 +143,7 @@ However, alternate fuels such as cryogenic liquid hydrogen require additional st Once this length is known, the masses of the tank and insulation layers can be found from their respective volumes and densities. To support the tank, stiffener rings that go around the tank circumference are needed. It is assumed that the inner vessel contains two of these stiffeners, which are sized following the process below. #### Outer vessel - If the insulation layer contains a vacuum layer, an additional tank wall is needed to contain the vacuum. Unlike the inner vessel, this outer wall does not fail due to excessive stress, instead, it fails by buckling or collapse[^3]. To prevent buckling, stiffener rings can be used to reduce the effective cylindrical lenght that can collapse. The unsupported cylinder length is + If the insulation layer contains a vacuum layer, an additional tank wall is needed to contain the vacuum. Unlike the inner vessel, this outer wall does not fail due to excessive stress, instead, it fails by buckling or collapse[^3]. To prevent buckling, stiffener rings can be used to reduce the effective cylindrical length that can collapse. The unsupported cylinder length is ```math L = \frac{l_{cyl}}{N_{stiff} -1}, ``` diff --git a/src/IO/read_input.jl b/src/IO/read_input.jl index 8b25b5a5..330da82d 100644 --- a/src/IO/read_input.jl +++ b/src/IO/read_input.jl @@ -393,12 +393,15 @@ if pari[iifwing] == 0 #If fuel is stored in fuselage end #Calculate fuel temperature and density as a function of pressure - Tfuel, ρfuel = cryo_fuel_properties(uppercase(fueltype), fuse_tank.ptank) + Tfuel, ρfuel, ρgas, hvap = cryo_fuel_properties(uppercase(fueltype), fuse_tank.ptank) pare[ieTft, :, :] .= Tfuel #Temperature of fuel in fuel tank #TODO remove this and replace with the one in struct pare[ieTfuel, :, :] .= Tfuel #Initialize fuel temperature as temperature in tank parg[igrhofuel] = ρfuel fuse_tank.rhofuel = ρfuel fuse_tank.Tfuel = Tfuel + fuse_tank.hvap = hvap + parg[igrhofuelgas] = ρgas + fuse_tank.rhofuelgas = ρgas end # --------------------------------- # Wing diff --git a/src/fuel/fuel_properties.jl b/src/fuel/fuel_properties.jl index d7cb10f3..42ffddfb 100644 --- a/src/fuel/fuel_properties.jl +++ b/src/fuel/fuel_properties.jl @@ -1,25 +1,46 @@ """ cryo_fuel_properties(fuel::String, p::Float64) -Calculates the temperature and density of a cryogenic fuel inside a tank, using fits to data from NIST. +Calculates the saturation temperature and density of a cryogenic fuel inside a tank, using fits to data from NIST. !!! details "🔃 Inputs and Outputs" **Inputs:** - `fuel::String`: name of the fuel. - `p::String`: tank pressure (Pa). **Outputs:** - - `Tfuel::Float64`: fuel evaporation temperature (K) - - `ρfuel::String`: density of fuel (kg/m^3). + - `Tsat::Float64`: saturation temperature (K) + - `ρl::Float64`: liquid density (kg/m^3). + - `ρg::Float64`: gas density (kg/m^3). + - `hv::Float64`: enthalpy of vaporization (J/kg). """ function cryo_fuel_properties(fuel::String, p::Float64) x = p / 101325 #pressure in atm - #Fits to NIST data at p = [0.5,1,1.5,2,2.5,3,4,5,10] atm - if fuel == "LH2" - Tfuel = 20.3839 * x^0.182022 - ρfuel = -2.88304E-02*x^3 + 4.93319E-01*x^2 - 4.61334E+00*x + 7.51570E+01 - elseif fuel == "CH4" - Tfuel = 111.69030 * x^0.12093 - ρfuel = -1.49198E-01*x^3 + 2.85956E+00*x^2 - 2.20464E+01*x + 4.42724E+02 + #Fits to NIST data at p increments of 0.1 atm from 0.1 to 10 atm + if uppercase(fuel) == "LH2" + #Saturation temperature (K) + Tsat = -3.35772E-04*x^6 + 1.14637E-02*x^5 - 1.54871E-01*x^4 + 1.05803E+00*x^3 - 3.93770E+00*x^2 + 9.09205E+00*x + 1.42913E+01 + + #Liquid density (kg/m^3) + ρl = 2.58354E-04*x^6 - 8.90930E-03*x^5 + 1.21248E-01*x^4 - 8.37637E-01*x^3 + 3.14780E+00*x^2 - 8.42843E+00*x + 7.68629E+01 + #Liquid enthalpy (kJ/kg) + hl = -2.258517E-03*x^6 + 7.767483E-02*x^5 - 1.055739E+00*x^4 + 7.280899E+00*x^3 - 2.739412E+01*x^2 + 7.379258E+01*x - 5.277765E+01 + + #Gas density (kg/m^3) + ρg = 5.61446E-03*x^3 - 4.19141E-02*x^2 + 1.28761E+00*x + 6.91405E-02 + #Gas enthalpy (kJ/kg) + hg = -3.31730E-03*x^6 + 1.12617E-01*x^5 - 1.51499E+00*x^4 + 1.02766E+01*x^3 - 3.79360E+01*x^2 + 7.33546E+01*x + 4.04343E+02 + elseif uppercase(fuel) == "CH4" + Tsat = -9.38765E-04*x^6 + 3.25637E-02*x^5 - 4.49376E-01*x^4 + 3.16285E+00*x^3 - 1.22945E+01*x^2 + 3.00925E+01*x + 9.09645E+01 + + ρl = 1.23474E-03*x^6 - 4.28482E-02*x^5 + 5.91776E-01*x^4 - 4.17205E+00*x^3 + 1.62724E+01*x^2 - 4.14826E+01*x + 4.51398E+02 + hl = -3.14958E-03*x^6 + 1.09246E-01*x^5 - 1.50774E+00*x^4 + 1.06164E+01*x^3 - 4.13089E+01*x^2 + 1.02732E+02*x - 7.11737E+01 + + ρg = 1.66091E-03*x^3 - 2.82413E-02*x^2 + 1.69557E+00*x + 1.33041E-01 + hg = -1.88367E-03*x^6 + 6.52861E-02*x^5 - 8.99795E-01*x^4 + 6.31818E+00*x^3 - 2.44324E+01*x^2 + 5.58852E+01*x + 4.73575E+02 end - return Tfuel, ρfuel + + #hl and hg are in kJ/kg + hv = (hg - hl) * 1e3 #J/kg, Enthalpy of vaporization is enthalpy difference at saturation + + return Tsat, ρl, ρg, hv end diff --git a/src/misc/aircraft.jl b/src/misc/aircraft.jl index 832111a5..a083a667 100644 --- a/src/misc/aircraft.jl +++ b/src/misc/aircraft.jl @@ -27,6 +27,8 @@ mutable struct fuselage_tank ptank::Float64 rhofuel::Float64 Tfuel::Float64 + rhofuelgas::Float64 + hvap::Float64 boiloff_rate::Float64 ftankadd::Float64 diff --git a/src/misc/index.inc b/src/misc/index.inc index d2c02027..21faaf38 100644 --- a/src/misc/index.inc +++ b/src/misc/index.inc @@ -385,8 +385,9 @@ const igaislehalfwidth = 308 const igmdotboiloff = 309 + const igrhofuelgas = 310 - const igtotal = 309 + const igtotal = 310 # indices for turbo-electric systems - really just the electrical machines const ite_ratSM = 1 diff --git a/src/mission/mission.jl b/src/mission/mission.jl index 8ae72cba..e733c2fe 100644 --- a/src/mission/mission.jl +++ b/src/mission/mission.jl @@ -73,6 +73,8 @@ function mission!(pari, parg, parm, para, pare, Ldebug)#, iairf, initeng, ipc1) #Boiloff rate for cryogenic fuels mdot_boiloff = parg[igmdotboiloff] + ρfgas = parg[igrhofuelgas] + ρf = parg[igrhofuel] # set known operating conditions @@ -408,7 +410,10 @@ function mission!(pari, parg, parm, para, pare, Ldebug)#, iairf, initeng, ipc1) pare[ieFe, ip] = Ftotal # Store integrands for range and weight integration using a predictor-corrector scheme FoW[ip] = Ftotal / (BW * cosg) - DoL - FFC[ip] = Ftotal * TSFC / (W * V * cosg) + gee * mdot_boiloff / (W * cosg * V) #second term accounts for fuel boiloff in cryo tanks + + mfuel = Ftotal * TSFC / gee + mdot_vent = max(mdot_boiloff * (1 - ρfgas/ρf) - mfuel * ρfgas/ρf, 0) #Vent boiloff gas if excessive + FFC[ip] = (mfuel + mdot_vent) * gee / (W * V * cosg) #second term accounts for fuel venting in cryo tanks Vgi[ip] = 1.0 / (V * cosg) @@ -530,8 +535,9 @@ function mission!(pari, parg, parm, para, pare, Ldebug)#, iairf, initeng, ipc1) cosg = cos(gamVcr1) FoW[ip] = Ftotal / (BW * cosg) - DoL - FFC[ip] = Ftotal * TSFC / (W * V * cosg) + gee * mdot_boiloff / (W * cosg * V) #second term accounts for fuel boiloff in cryo tanks - + mfuel = Ftotal * TSFC / gee + mdot_vent = max(mdot_boiloff * (1 - ρfgas/ρf) - mfuel * ρfgas/ρf, 0) #Vent boiloff gas if excessive + FFC[ip] = (mfuel + mdot_vent) * gee / (W * V * cosg) #second term accounts for fuel venting in cryo tanks Vgi[ip] = 1.0 / (V * cosg) #---- set end-of-cruise point "d" using cruise and descent angles, @@ -597,9 +603,12 @@ function mission!(pari, parg, parm, para, pare, Ldebug)#, iairf, initeng, ipc1) cosg = cos(gamVcr1) FoW[ip] = Ftotal / (BW * cosg) - DoL - FFC[ip] = Ftotal * TSFC / (W * V * cosg) + gee * mdot_boiloff / (W * cosg * V) #second term accounts for fuel boiloff in cryo tanks - Vgi[ip] = 1.0 / (V * cosg) + mfuel = Ftotal * TSFC / gee + mdot_vent = max(mdot_boiloff * (1 - ρfgas/ρf) - mfuel * ρfgas/ρf, 0) #Vent boiloff gas if excessive + FFC[ip] = (mfuel + mdot_vent) * gee / (W * V * cosg) #second term accounts for fuel venting in cryo tanks + + Vgi[ip] = 1.0 / (V * cosg) ip1 = ipcruise1 ipn = ipcruisen @@ -769,13 +778,13 @@ function mission!(pari, parg, parm, para, pare, Ldebug)#, iairf, initeng, ipc1) # store integrands for Range and Weight integration FoW[ip] = F / (BW * cosg) - DoL - FFC[ip] = F / (W * V * cosg) * TSFC + gee * mdot_boiloff / (W * cosg * V) #second term accounts for fuel boiloff in cryo tanks + Vgi[ip] = 1.0 / (V * cosg) # if F < 0, then TSFC is not valid, so calculate mdot_fuel directly mfuel = pare[ieff, ip] * pare[iemcore, ip] * parg[igneng] - FFC[ip] = gee * mfuel / (W * cosg * V) + gee * mdot_boiloff / (W * cosg * V) #second term accounts for fuel boiloff in cryo tanks - + mdot_vent = max(mdot_boiloff * (1 - ρfgas/ρf) - mfuel * ρfgas/ρf, 0) #Vent boiloff gas if excessive + FFC[ip] = (mfuel + mdot_vent) * gee / (W * V * cosg) #second term accounts for fuel venting in cryo tanks if (ip > ipdescent1) # corrector integration step, approximate trapezoidal diff --git a/src/structures/tankWmech.jl b/src/structures/tankWmech.jl index 652f36f0..8240e332 100644 --- a/src/structures/tankWmech.jl +++ b/src/structures/tankWmech.jl @@ -43,6 +43,7 @@ function size_inner_tank(fuse_tank, t_cond::Vector{Float64}) Wfuel = fuse_tank.Wfuelintank ρfuel = fuse_tank.rhofuel + ρfgas = fuse_tank.rhofuelgas ftankadd = fuse_tank.ftankadd Δp = fuse_tank.ptank sigskin = fuse_tank.inner_material.UTS @@ -72,8 +73,9 @@ function size_inner_tank(fuse_tank, t_cond::Vector{Float64}) #--- Calculate length of cylindrical portion Wfuel_tot = Wfuel #Wfuel already includes the amount that boils off - Vfuel = Wfuel_tot / (gee * ρfuel) - Vinternal = (1 + ullage_frac)*Vfuel # required interal volume + ρfmix = (1 - ullage_frac) * ρfuel + ullage_frac * ρfgas #Density of saturated mixture in tank + Vfuel = Wfuel_tot / (gee * ρfmix) #Total tank volume taken by saturated mixture + Vinternal = Vfuel # required internal volume V_ellipsoid = 2π * (Rtank^3 / AR) / 3 # Half the vol of std ellipsoid = 1/2×(4π/3 ×(abc)) where a,b,c are the semi-axes length. Here a = R/AR, b=c=R # Also see: https://neutrium.net/equipment/volume-and-wetted-area-of-partially-filled-horizontal-vessels/ V_cylinder = Vinternal - 2*V_ellipsoid @@ -439,6 +441,8 @@ function insulation_density_calc(material::String) ρ = 32.0 #kg/m^3 elseif lowercase(material) == "polyurethane35" ρ = 35.0 #kg/m^3 + elseif lowercase(material) == "mylar" + ρ = 1390 #kg/m^3, https://www.matweb.com/search/datasheet_print.aspx?matguid=981d85aa72b0419bb4b26a3c06cb284d elseif lowercase(material) == "vacuum" ρ = 0 #kg/m^3 else diff --git a/src/structures/tankWthermal.jl b/src/structures/tankWthermal.jl index 2d5304cb..cd4b6198 100644 --- a/src/structures/tankWthermal.jl +++ b/src/structures/tankWthermal.jl @@ -38,6 +38,7 @@ function tankWthermal(fuse_tank, z::Float64, Mair::Float64, xftank::Float64, tim t_cond = fuse_tank.t_insul Tfuel = fuse_tank.Tfuel qfac = fuse_tank.qfac + h_v = fuse_tank.hvap #heat of vaporization Wtank_total, l_cyl, tskin, r_tank, Vfuel, Wtank, Wfuel_tot, Winsul_sum, t_head, Whead, Wcyl, Wstiff, Winsul, Sinternal, Shead, l_tank = size_inner_tank(fuse_tank, fuse_tank.t_insul) @@ -73,8 +74,6 @@ function tankWthermal(fuse_tank, z::Float64, Mair::Float64, xftank::Float64, tim end sol = nlsolve(fun, guess, xtol = 1e-7, ftol = 1e-6) #Solve non-linear problem with NLsolve.jl - _, h_v = tank_heat_coeffs(Tfuel, ifuel, Tfuel, l_tank) #Liquid heat of vaporization - Q = qfac * sol.zero[1] # Heat rate from ambient to cryo fuel, including extra heat leak from valves etc as in eq 3.20 by Verstraete mdot_boiloff = Q / h_v # Boil-off rate equals the heat rate divided by heat of vaporization m_boiloff = mdot_boiloff * time_flight # Boil-off mass calculation @@ -144,7 +143,7 @@ function residuals_Q(x::Vector{Float64}, p, mode::String) Rair_conv_rad = 1 / (h_air * (2π * (r_tank + thickness) * l_cyl + 2*Shead[end])) # thermal resistance of ambient air (incl. conv and rad) S_int = (2π * (r_inner) * l_cyl) + 2*Shead[1] #liquid side surface area - h_liq, _ = tank_heat_coeffs(T_w, ifuel, Tfuel, l_tank) #Find liquid-side heat transfer coefficient + h_liq = tank_heat_coeff(T_w, ifuel, Tfuel, l_tank) #Find liquid-side heat transfer coefficient R_liq = 1 / (h_liq * S_int) #Liquid-side thermal resistance N = length(t_cond) # Number of layers in insulation @@ -205,18 +204,20 @@ This function calculates the thermal conductivity of different insulation materi - `k::Float64`: thermal conductivity (W/(m K)). """ function insulation_conductivity_calc(T::Float64, material::String) - if material == "rohacell41S" + if lowercase(material) == "rohacell41s" #Note: Brewer (1991) only had one point for rohacell thermal conductivity. They assumed the same curve as PVC k = 0.001579 + 1.283e-4 * T - 3.353e-7*T^2 + 8.487e-10 * T^3 # W/(m K), polynomial fit to Fig. 4.78 in Brewer (1991) between 20 and 320 K - elseif material == "polyurethane27" #polyurethane with density of 27 kg/m^3 + elseif lowercase(material) == "polyurethane27" #polyurethane with density of 27 kg/m^3 k = 2.114e-13 * T^5 - 1.639e-10 *T^4 + 4.438e-8 * T^3 - 5.222e-6*T^2 + 3.740e-4*T - 2.192e-3 # W/(m K), polynomial fit to Fig. 4.78 in Brewer (1991) between 20 and 320 K - elseif material == "polyurethane32" #polyurethane with density of 32 kg/m^3 + elseif lowercase(material) == "polyurethane32" #polyurethane with density of 32 kg/m^3 k = 2.179E-13 * T^5 - 1.683E-10* T^4 + 4.542E-08* T^3 - 5.341E-06* T^2 + 3.816E-04* T - 2.367E-03 # W/(m K), polynomial fit to Fig. 4.78 in Brewer (1991) between 20 and 320 K - elseif material == "polyurethane35" #polyurethane with density of 35 kg/m^3 + elseif lowercase(material) == "polyurethane35" #polyurethane with density of 35 kg/m^3 k = 2.104E-13* T^5 - 1.695E-10* T^4 + 4.746E-08* T^3 - 5.662E-06* T^2 + 3.970E-04* T - 2.575E-03 # W/(m K), polynomial fit to Fig. 4.78 in Brewer (1991) between 20 and 320 K + elseif lowercase(material) == "mylar" + k = 0.155# W/(m K) at 298 K, https://www.matweb.com/search/datasheet_print.aspx?matguid=981d85aa72b0419bb4b26a3c06cb284d end return k end @@ -259,9 +260,9 @@ mutable struct thermal_params end """ - tank_heat_coeffs(T_w, ifuel, Tfuel, ltank) + tank_heat_coeff(T_w, ifuel, Tfuel, ltank) -This function calculates the latent heat of vaporization of a liquid fuel and its liquid-side heat transfer coefficient. +This function calculates the liquid-side heat transfer coefficient of a cryogenic fuel in a tank. !!! details "🔃 Inputs and Outputs" **Inputs:** @@ -272,19 +273,16 @@ This function calculates the latent heat of vaporization of a liquid fuel and it **Outputs:** - `h_liq::Float64`: liquid-side heat tansfer coefficient (W/m^2/K). - - `h_v::Float64`: liquid's enthalpy of vaporization (J/kg). """ -function tank_heat_coeffs(T_w::Float64, ifuel::Int64, Tfuel::Float64, ltank::Float64) +function tank_heat_coeff(T_w::Float64, ifuel::Int64, Tfuel::Float64, ltank::Float64) #Get fuel properties if ifuel == 11 #CH4 - h_v = 510e3 #J/kg, latent heat of vaporozation Pr_l = 2.0 #Prandtl number, from NIST at 2atm and 120 K β = 3.5e-3 #K^(-1), thermal expansion coefficient https://aiche.onlinelibrary.wiley.com/doi/10.1002/aic.14254 ν_l = 2.4e-7 #m^2/s, kinematic viscosity from NIST at 2atm and 120K k = 0.17185 #W/m/K, thermal conductivity from NIST at 2atm and 120K elseif ifuel == 40 #LH2 - h_v = 447e3 Pr_l = 1.3 #from NIST at 2atm and 20K β = 15e-3 #https://nvlpubs.nist.gov/nistpubs/jres/089/jresv89n4p317_A1b.pdf ν_l = 2.0e-7 #m^2/s, from NIST at 2atm and 20K @@ -297,7 +295,7 @@ function tank_heat_coeffs(T_w::Float64, ifuel::Int64, Tfuel::Float64, ltank::Flo h_liq = Nu_l * k / ltank #heat transfer coefficient for liquid - return h_liq, h_v + return h_liq end """ diff --git a/src/structures/tanksize.jl b/src/structures/tanksize.jl index eb76ee78..1b4a1f03 100644 --- a/src/structures/tanksize.jl +++ b/src/structures/tanksize.jl @@ -126,6 +126,7 @@ function res_MLI_thick(x::Vector{Float64}, fuse_tank, z::Float64, Mair::Float64, iinsuldes = fuse_tank.iinsuldes Tfuel = fuse_tank.Tfuel Wfuel = fuse_tank.Wfuelintank + h_v = fuse_tank.hvap #heat of vaporization # Extract states Δt = x[1] @@ -141,8 +142,6 @@ function res_MLI_thick(x::Vector{Float64}, fuse_tank, z::Float64, Mair::Float64, Wtank_total, l_cyl, tskin, r_tank, Vfuel, Wtank, Wfuel_tot, Winsul_sum, t_head, Whead, Wcyl, Wstiff, Winsul, Sinternal, Shead, l_tank = size_inner_tank(fuse_tank, t_all) - _, h_v = tank_heat_coeffs(Tfuel, ifuel, Tfuel, l_tank) #Liquid heat of vaporizatio - mdot_boiloff = boiloff_percent * Wfuel / (gee * 100) / 3600 # Boil-off rate equals the heat rate divided by heat of vaporization Q_net = mdot_boiloff * h_v # Heat rate from ambient to cryo fuel, including extra heat leak from valves etc as in eq 3.20 by Verstraete diff --git a/test/unit_test_fueltank.jl b/test/unit_test_fueltank.jl index a9a83f25..537f5dd6 100644 --- a/test/unit_test_fueltank.jl +++ b/test/unit_test_fueltank.jl @@ -32,8 +32,12 @@ fuse_tank.outer_material = TASOPT.StructuralAlloy("Al-2219-T87") fuse_tank.theta_outer = [1.0, 2.0] fuse_tank.Ninterm = 1.0 -fuse_tank.rhofuel = 70.0 -fuse_tank.Tfuel = 20.0 +Tsat, ρl, ρg, hv = TASOPT.cryo_fuel_properties("LH2", fuse_tank.ptank) + +fuse_tank.rhofuel = ρl +fuse_tank.rhofuelgas = ρg +fuse_tank.Tfuel = Tsat +fuse_tank.hvap = hv fuse_tank.Wfuelintank = 1e5 @testset "Fuselage tank" begin @@ -41,7 +45,7 @@ fuse_tank.Wfuelintank = 1e5 outputs_size = TASOPT.structures.tanksize!(fuse_tank, z, Mair, xftank, time_flight, ifuel) - outputs_size_check = (0.004247366632687734, 145.623998835008, 1.9102299646350587, 33161.58141606815, 15.64798968801848, 55769.4869666662) + outputs_size_check = (0.004247366632687734, 166.77327116515787, 1.8582041735501131, 39458.39577725139, 17.14040378395799, 62696.45129193434) for i in 1:length(outputs_size) @test outputs_size[i] ≈ outputs_size_check[i] @@ -49,18 +53,18 @@ fuse_tank.Wfuelintank = 1e5 outputs_mech = TASOPT.structures.size_inner_tank(fuse_tank, fuse_tank.t_insul) - outputs_mech_check = (155769.48696666618, 12.754661771033133, 0.0035645268979941384, 1.9102299646350587, 145.623998835008, 55769.4869666662, 100000.0, 33161.58141606815, 0.003561204259200333, 1565.0564653871652, 15934.428537303356, 1488.0999415568983, [10913.366525915251, 12048.59580855156, 10199.619081601335], 192.00106755661065, [15.844936883719953, 19.260709211746004, 23.014988363554576, 27.106823359071385], 15.64798968801848) + outputs_mech_check = (162696.45129193435, 14.195146976056426, 0.003467445742769384, 1.8582041735501131, 166.77327116515787, 62696.45129193434, 100000.0, 39458.39577725139, 0.003464213597222433, 1440.632953061823, 16803.608249490673, 1440.6308577338177, [12862.27440187786, 14345.62808332926, 12250.493292044272], 203.79565769806428, [14.993605136580461, 18.69431805929931, 22.809097377939224, 27.336717237299393], 17.14040378395799) for i in 1:length(outputs_mech) @test outputs_mech[i] ≈ outputs_mech_check[i] end - + l_cyl = outputs_mech_check[2] l_tank = outputs_mech_check[16] r_tank = outputs_mech_check[4] Shead = outputs_mech_check[15] outputs_thermal = TASOPT.structures.tankWthermal(fuse_tank, z, Mair, xftank, time_flight, ifuel) - outputs_thermal_check = (107.03363914390675, 0.004247366632694712) + outputs_thermal_check = (107.0336391439565, 0.0042473666326966865) for i in 1:length(outputs_thermal) @test outputs_thermal[i] ≈ outputs_thermal_check[i] @@ -76,14 +80,14 @@ fuse_tank.Wfuelintank = 1e5 outputs_vac_size = TASOPT.structures.tanksize!(fuse_tank, z, Mair, xftank, time_flight, ifuel) - outputs_vac_size_check = (0.0020612680724211444, 145.623998835008, 2.4, 0.0, 10.488026633984164, 111294.54283358467) + outputs_vac_size_check = (0.0022038332318087754, 166.77327116515787, 2.4, 0.0, 10.884763251319967, 114902.51923675802) for i in 1:length(outputs_vac_size) @test outputs_vac_size[i] ≈ outputs_vac_size_check[i] end outputs_vac_mech = TASOPT.structures.size_inner_tank(fuse_tank, fuse_tank.t_insul) - outputs_vac_mech_check = (124335.29733608443, 7.833183975356766, 0.00435715618585557, 2.335, 145.623998835008, 24335.29733608443, 100000.0, 0.0, 0.004353094705443698, 2858.4704423225353, 14495.185531216492, 1910.8711623970096, [0.0], 166.54673585196167, [23.675162947566548, 25.28189535258484], 10.302533008581799) + outputs_vac_mech_check = (125124.15434149424, 8.219175642639724, 0.00435715618585557, 2.335, 166.77327116515787, 25124.15434149424, 100000.0, 0.0, 0.004353094705443698, 2858.4704423225353, 15209.457128611886, 1913.7422971923515, [0.0], 172.43073912598143, [23.675162947566548, 25.28189535258484], 10.688524675864757) for i in 1:length(outputs_vac_mech) @test outputs_vac_mech[i] ≈ outputs_vac_mech_check[i] @@ -93,11 +97,11 @@ fuse_tank.Wfuelintank = 1e5 fuse_tank.Ninterm = 1.0 Ninterm = TASOPT.structures.optimize_outer_tank(fuse_tank, Winnertank, l_cyl) - Ninterm_check = 14.63037109375 + Ninterm_check = 15.38916015625 @test Ninterm ≈ Ninterm_check outputs_vac_outer = TASOPT.structures.size_outer_tank(fuse_tank, Winnertank, l_cyl, Ninterm_check) - outputs_vac_outer_check = (85410.79548689285 , 31471.569231536196, 12808.143626608879, 20558.32123060317, 168.1077334143945, 24.993050952321795, 118.12163150975091, 0.009563164226558287 , 0.018394143361195915, 10.269972262079158) + outputs_vac_outer_check = (88153.61722490133, 33033.25329386777, 12808.143626608879, 21490.11147555205, 173.928350628553, 24.993050952321795, 123.9422487239094, 0.009566313929022768, 0.018394143361195915, 10.655963929362116) for i in 1:length(outputs_vac_outer) @test outputs_vac_outer[i] ≈ outputs_vac_outer_check[i] end @@ -128,8 +132,8 @@ fuse_tank.Wfuelintank = 1e5 end @testset "Thermal models" begin - outputs_h = TASOPT.structures.tank_heat_coeffs(21.0, ifuel, fuse_tank.Tfuel, 5.0) - outputs_h_check = (105.81383176854604, 447000.0) + outputs_h = TASOPT.structures.tank_heat_coeff(21.0, ifuel, fuse_tank.Tfuel, 5.0) + outputs_h_check = (133.43555789562762) for i in 1:length(outputs_h) @test outputs_h[i] ≈ outputs_h_check[i] end @@ -142,9 +146,37 @@ fuse_tank.Wfuelintank = 1e5 Taw = outputs_free_check[3] Rvac = TASOPT.structures.vacuum_resistance(fuse_tank.Tfuel, Taw, 90.0, 100.0) - Rvac_check = 0.5500815195658247 + Rvac_check = 0.543113131740449 @test Rvac ≈ Rvac_check end + @testset "Fuel properties" begin + ps = [1, 2, 5] * 101325.0 #pressures to test at + + fuel = "LH2" + outputs_lh2_check = [20.369 70.848 1.3322 448.71e3; 22.965 67.639 2.5133 431.474e3; 27.314 60.711 6.1715 371.36e3] + + for (i, p) in enumerate(ps) + outputs_lh2 = TASOPT.cryo_fuel_properties(fuel, p) + + outputs_check = outputs_lh2_check[i, :] + for j in 1:length(outputs_check) + @test outputs_lh2[j] ≈ outputs_check[j] rtol = 1e-2 #Only require 1% diff as function uses fits + end + end + + fuel = "CH4" + outputs_ch4_check = [111.67 422.36 1.8164 510.83e3; 120.81 409.78 3.4382 492.921e3; 135.59 384.63 8.1024 457.614e3] + + for (i, p) in enumerate(ps) + outputs_ch4 = TASOPT.cryo_fuel_properties(fuel, p) + + outputs_check = outputs_ch4_check[i, :] + for j in 1:length(outputs_check) + @test outputs_ch4[j] ≈ outputs_check[j] rtol = 1e-2 #Only require 1% diff as function uses fits + end + end + + end end \ No newline at end of file