Line data Source code
1 : //
2 : // This SCP model implements the method described in the following paper
3 : //
4 : // .. bibliography::
5 : // :list: none
6 : // :filter: False
7 : //
8 : // meier2024diffuse
9 : //
10 : // This method imitates inheritance :code:`Propellant::Propellant` using static polymorphism
11 : //
12 :
13 : #ifndef MODEL_PROPELLANT_FULLFEEDBACK_H
14 : #define MODEL_PROPELLANT_FULLFEEDBACK_H
15 :
16 : #include "Propellant.H"
17 :
18 : namespace Model
19 : {
20 : namespace Propellant
21 : {
22 : class FullFeedback
23 : {
24 : public:
25 : static constexpr const char* name = "fullfeedback";
26 :
27 3 : FullFeedback() = default;
28 : FullFeedback(IO::ParmParse &pp, std::string name)
29 : {
30 : pp_queryclass(name,*this);
31 : }
32 :
33 2550 : void set_pressure(const Set::Scalar a_P)
34 : {
35 2550 : this->Pbar = a_P/Pref;
36 :
37 2550 : zeta_2 = zeta_fit_a + Pbar * zeta_fit_b;
38 :
39 2550 : if (arrhenius_dependency == 1) zeta_1 = zeta_2;
40 0 : else zeta_1 = zeta_0;
41 :
42 : // Pressure-dependent post calculation variables
43 2550 : k1 = a1 * Pbar + b1 - zeta_1 / zeta;
44 2550 : k2 = a2 * Pbar + b2 - zeta_1 / zeta;
45 2550 : k3 = 4.0 * log((c1*Pbar*Pbar + a3*Pbar + b3) - k1 / 2.0 - k2 / 2.0);
46 :
47 2550 : return;
48 : }
49 :
50 : AMREX_FORCE_INLINE
51 : Set::Scalar get_qdot(const Set::Scalar mdot, const Set::Scalar phi)
52 : {
53 1747840 : Set::Scalar qflux = k1 * phi +
54 1747840 : k2 * (1.0 - phi) +
55 1747840 : (zeta_1 / zeta) * exp(k3 * phi * (1.0 - phi));
56 1747840 : Set::Scalar mlocal = (mlocal_ap) * phi + (mlocal_htpb) * (1.0 - phi) + mlocal_comb * phi * (1.0 - phi);
57 1747840 : Set::Scalar mdota = fabs(mdot);
58 1747840 : Set::Scalar mbase = tanh(4.0 * mdota / (mlocal));
59 :
60 1747840 : return mbase * qflux;
61 : }
62 :
63 :
64 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
65 : Set::Scalar get_K (const Set::Scalar phi)
66 : {
67 1747840 : return k_ap * phi + k_htpb * (1.0 - phi);
68 : }
69 :
70 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
71 : Set::Scalar get_rho (const Set::Scalar phi)
72 : {
73 1747840 : return rho_ap * phi + rho_htpb * (1.0 - phi);
74 : }
75 :
76 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
77 : Set::Scalar get_cp (const Set::Scalar phi)
78 : {
79 1747840 : return cp_ap * phi + cp_htpb * (1.0 - phi);
80 : }
81 :
82 : AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE
83 : Set::Scalar get_L (Set::Scalar phi, Set::Scalar T)
84 : {
85 1747840 : Set::Scalar L = NAN;
86 1747840 : if (mob_ap == 1) L = m_ap * Pbar * exp(-E_ap / T) * phi;
87 1747840 : else L = m_ap * exp(-E_ap / T) * phi;
88 1747840 : L += m_htpb * exp(-E_htpb / T) * (1.0 - phi);
89 1747840 : if (T <= bound) L = 0.0;
90 1747840 : return L;
91 : }
92 :
93 :
94 1 : static void Parse(FullFeedback& value, IO::ParmParse& pp)
95 : {
96 : // AP/HTPB interface length
97 4 : pp.query_default("phi.zeta", value.zeta, "10.0_um", Unit::Length());
98 : // Reference interface length for heat integration
99 4 : pp.query_default("phi.zeta_0", value.zeta_0, "10.0_um", Unit::Length());
100 :
101 : // Constant offset for $\zeta$ adjustment with pressure
102 4 : pp.query_default("phi.zeta_fit_a",value.zeta_fit_a, "45.0_um", Unit::Length());
103 : // Scaling factors for $\zeta$ adjustment with pressure
104 4 : pp.query_default("phi.zeta_fit_b",value.zeta_fit_b, "-6.42_um", Unit::Length());
105 :
106 : // Surrogate heat flux model paramater - AP
107 2 : pp_query_required("a1", value.a1, Unit::Less());
108 : // Surrogate heat flux model paramater - HTPB
109 2 : pp_query_required("a2", value.a2,Unit::Less());
110 : // Surrogate heat flux model paramater - Total
111 2 : pp_query_required("a3", value.a3,Unit::Less());
112 : // Surrogate heat flux model paramater - AP
113 2 : pp_query_required("b1", value.b1, Unit::Less());
114 : // Surrogate heat flux model paramater - HTPB
115 2 : pp_query_required("b2", value.b2, Unit::Less());
116 : // Surrogate heat flux model paramater - Total
117 2 : pp_query_required("b3", value.b3,Unit::Less());
118 : // Surrogate heat flux model paramater - Total
119 2 : pp_query_required("c1", value.c1, Unit::Less());
120 :
121 : // Whether to use pressure to determined the reference Zeta
122 2 : pp_query_default("pressure.dependency", value.arrhenius_dependency, true);
123 :
124 : // Reference pressure for pressure nondimensionalization
125 4 : pp_query_default("Pref",value.Pref,"1_MPa",Unit::Pressure());
126 :
127 : // AP volumetric mass flux reference value
128 4 : pp.query_default("mlocal_ap", value.mlocal_ap, "0.0", Unit::Density()/Unit::Time());
129 : // combined volumetric mass flux reference value
130 4 : pp_query_default("mlocal_comb", value.mlocal_comb, "0.0", Unit::Density()/Unit::Time());
131 : // HTPB volumetric mass flux reference value
132 4 : pp.query_default("mlocal_htpb",value.mlocal_htpb, "5000_kg/m^3/s", Unit::Density()/Unit::Time());
133 :
134 : // Set::Scalar massfraction = 0.8;
135 : // value.mlocal_htpb = 685000.0 - 850e3 * massfraction; // Replicating but needs to be revised!
136 : // Util::Message(INFO,value.mlocal_htpb);
137 :
138 : // AP Density
139 2 : pp_query_required("rho_ap", value.rho_ap, Unit::Density());
140 : // HTPB Density
141 2 : pp_query_required("rho_htpb", value.rho_htpb, Unit::Density());
142 : // AP Thermal Conductivity
143 2 : pp_query_required("k_ap", value.k_ap, Unit::ThermalConductivity());
144 : // HTPB Thermal Conductivity
145 2 : pp_query_required("k_htpb", value.k_htpb, Unit::ThermalConductivity());
146 : // AP Specific Heat
147 2 : pp_query_required("cp_ap", value.cp_ap,Unit::SpecificHeatCapacity());
148 : //HTPB Specific Heat
149 2 : pp_query_required("cp_htpb", value.cp_htpb,Unit::SpecificHeatCapacity());
150 :
151 : // AP Pre-exponential factor for Arrhenius Law
152 2 : pp_query_required("m_ap", value.m_ap,Unit::Velocity());
153 : // HTPB Pre-exponential factor for Arrhenius Law
154 2 : pp_query_required("m_htpb", value.m_htpb,Unit::Velocity());
155 : // COMBINED AP Activation Energy, divided by Rg, for Arrhenius Law (has units of temperature)
156 2 : pp_query_required("E_ap", value.E_ap,Unit::Temperature());
157 : // COMBINED HTPB Activation Energy, divided by Rg, for Arrhenius Law(has units of temperature)
158 2 : pp_query_required("E_htpb", value.E_htpb,Unit::Temperature());
159 : // Whether to include pressure to the arrhenius law [??]n
160 2 : pp_query_default("mob_ap", value.mob_ap,false);
161 : // HTPB Activation Energy for Arrhenius Law
162 4 : pp_query_default("bound", value.bound,"0.0", Unit::Temperature());
163 1 : }
164 :
165 : // IO Variables
166 : Set::Scalar zeta = NAN, zeta_0 = NAN;
167 : Set::Scalar a1 = NAN, a2 = NAN, a3 = NAN;
168 : Set::Scalar b1 = NAN, b2 = NAN, b3 = NAN;
169 : Set::Scalar c1 = NAN;
170 : //Set::Scalar h1 = NAN, h2 = NAN;
171 : bool arrhenius_dependency = false;
172 : Set::Scalar mlocal_ap = NAN, mlocal_comb = NAN, mlocal_htpb = NAN;
173 :
174 : //Thermal properties
175 : Set::Scalar k_ap = NAN, k_htpb = NAN;
176 : Set::Scalar rho_ap = NAN, rho_htpb = NAN;
177 : Set::Scalar cp_ap = NAN, cp_htpb = NAN;
178 :
179 : // Regression variables
180 : int mob_ap = -1;
181 : Set::Scalar m_ap = NAN, m_htpb = NAN, m_comb = NAN;
182 : Set::Scalar E_ap = NAN, E_htpb = NAN;
183 : Set::Scalar bound = NAN;
184 : Set::Scalar Pref = NAN;
185 :
186 : // Calculated variables
187 : Set::Scalar Pbar = NAN;
188 : Set::Scalar zeta_1 = NAN, zeta_2 = NAN;
189 : Set::Scalar zeta_fit_a = NAN, zeta_fit_b = NAN;
190 : Set::Scalar k1=NAN, k2=NAN, k3=NAN, k4=NAN;
191 : };
192 :
193 : }
194 : }
195 :
196 :
197 :
198 : #endif
199 :
200 :
|