CC BY
30
UDC 621.454.3
DOI: 10.18698/0236-3941-2021-3-98-107
NUMERICAL COMPUTATION OF SPECIFIC IMPULSE AND INTERNAL FLOW PARAMETERS IN SOLID FUEL ROCKET MOTORS WITH TWO-PHASE COMBUSTION PRODUCTS
T.S. Sultanov G.A. Glebov
[email protected] [email protected]
Kazan National Research Technical University named after A.N. Tupolev-KAI, Kazan, Russian Federation
Abstract
Eulerian — Lagrangian method was used in the Fluent computational fluid dynamics system to calculate motion of the two-phase combustion products in the solid fuel rocket motor combustion chamber and nozzle. Condensed phase is assumed to consist of spherical particles with the same diameter, which dimensions are not changing along the motion trajectory. Flows with particle diameters of 3, 5, 7, 9, and 11 ^m were investigated. Four versions of the engine combustion chamber configuration were examined: with slotted and smooth cylindrical charge channels, each with external and submerged nozzles. Gas flow and particle trajectories were calculated starting from the solid fuel surface and to the nozzle exit. Volumetric fields of particle concentrations, condensed phase velocities and temperatures, as well as turbulence degree in the solid propellant rocket engine flow duct were obtained. Values of particles velocity and temperature lag from the gas phase along the nozzle length were received. Influence of the charge channel shape, degree of the nozzle submersion and of the condensate particles size on the solid propellant rocket engine specific impulse were determined, and losses were estimated in comparison with the case of ideal flow
Keywords
Solid propellant rocket engines, numerical methods, two-phase flows, specific impulse, supersonic nozzle
Received 04.11.2020 Accepted 16.12.2020 © Author(s), 2021
Introduction. Due to real thermal and gas-dynamic processes in the solid fuel rocket motors (SRM), thrust and specific impulse differ from their theoretical (ideal) values obtained by thermodynamics calculation or using gas-dynamic relations [1-4]. The difference between the Isp specific impulse actual value from its Isp.id ideal value is usually considered the specific impulse loss factor [1]:
Cc = l^f^L = £ Q, (1)
Isp.id i
where n is the number of types of losses; Q is the coefficient of losses related to various physical processes (including those caused by flow scattering at the nozzle exit C,sc; friction inside the nozzle C, fr; chemical nonequilibrium C,n; multiphase effect C,s and other losses C,otk).
Losses due to submerged nozzle, nozzle flare, etc. are considered as other losses. Work [5] considers five types of losses, work [6] — eight, and work [7] — ten types of losses. Semi-empirical or empirical methods were developed to calculate each type of loss. Friction losses can be determined using the boundary layer methods [8-10], methods for determining two-phase losses can be found in [1, 2] and the method for determining losses due to chemical nonequilibrium can be found in [2]. Empirical formulas for the C,rec, C,fr, C,sc, C,s calculation obtained on the experiment basis are provided in [11].
Disadvantages of this approach in identifying the real specific impulse Isp are the following. First, each type of loss under consideration is determined without taking into account their mutual influence on each other. Second, summation of a large number of the Q coefficients leads to a significant error.
According to [12], the AC,c total root-mean-square error could be determined as following:
ACc (AC; )2. (2)
If it is assumed that the AQ error in determining each type of loss is 10 %, then for n = 6 the total error will be AC,c = 25 % (for AQ = 15 % and n = 6, we'll get A^c = 37.5 %). The error is rather significant.
The same can be said about the calculation of the discharge coefficient = m / mdd, where mid is the flow rate determined in an approximation of the equilibrium flow in the nozzle, and m is the actual flow rate. In [11], this
7
coefficient is defined as ^ = , where are the discharge coefficient
1
components which depend on the nozzle geometric parameters and other factors. If we'd assume = 10 % as an error in determining each coefficient,
[7 2
then the total error would be A^c = () « 26.5 %. It seems reasonable to
determine the and Isp flow parameters directly by numerical simulation taking into consideration real physical processes in the SRM flow duct.
Methods and assumptions. SRM combustion products flow was simulated for the following flow duct configurations: smooth cylindrical channel and an external nozzle; cylindrical channel and nozzle submerged by 30 % of the nozzle length; channel with 8 slots and an external nozzle, channel with 8 slots and nozzle recessed by 30 %. Considered configurations are shown in Fig. 1.
Fig. 1. Considered versions of the engine channel geometry
Metallized composite fuel was used as the charge. The fuel has the following composition: NH4ClO4 60 %, Al 15 %, rubber 20 %, other components 5 %. Combustion products composition and properties in the engine chamber were determined based on thermodynamics calculation at pressure in the chamber P0c = 9 MPa using ASTRA software [13]. AhO3 condensate mass fraction was z = 0.3 of the combustion products total mass. Burning surface and burning rate for both configurations of the charge channel were the same. The nozzle was of conical shape with the 38 % opening angle and geometric expansion ratio of Fa = 8.4. Minimum nozzle section diameter was 200 mm. Burning surface area and burn rate were the same in all cases. Nozzle and combustion chamber walls were considered adiabatic. When the particle trajectory hits the wall, it was assumed that the particle gets "frozen" on the wall [11, 14], and its trajectory is terminated.
Polydisperse flow is replaced by monodisperse flow with the equivalent particle diameter. This work represents the condensed phase by particles of the same size corresponding to the d43 mass-average diameter values [2]. Particle shape is assumed to be spherical. Calculation was performed for particles with the following diameters: d43 = 3, 5, 7, 9 and 11 ^m. Particle sizes are not changing along their trajectory, i.e., it is assumed that condensation process occurs instantly near the solid fuel surface. The process of particles breaking and coagulation is not taken into account in this work.
Working fluid is considered consisting of two phases: continuous gas and dispersed condensed phases. Flow parameters at the nozzle inlet are determined as a result of numerical calculation of the combustion chamber flow.
The well-known condition of supersonic output is used as a boundary condition at the combustion chamber exit (nozzle exit section). The wall boundary condition is specified by the standard wall function. The mesh is structured, the mesh cell size is 5 mm, the total number of cells is about 1.5 million.
To calculate the gas phase flow, the Favre averaged Navier — Stokes equations are used:
dp d (— . n ~dt {m ) = 0;
dt
dp d
dxi dxj
d
(ц + ^t)
dxj
( düi dü j
v dx j dxi
— 5,j
дщ
ЛЛ
3 4 dxi
J J
closed by the Realizable k-s turbulence model:
о
(X + Xt)
dT
dx
+ Ss
j J
i(Pk ) + l^(Püjk)= 6
f
dt
dxj
dxj
ш Л dk W +—1
Ok J dxj
+ Ц -pe;
¿(P8) + l^(PSüj )= 0
dt
dxj
dxj.
+ щ ^ ds ÜR J dx j
л
+ pQSs-pC2
j
k + Vvs
where
Ci = max 0.43; ^ ; ц = S —; S = J2SijSij;
l Л+1) e
k2
Mt = PC|i-; C|l='
Ao As
kU
U — ^ jj + OjOij ; Oij — Qjj 2Sxjk®k; ^ij ~ ^ij &ijk®k;
1 S'S S
A0 = 4.04; AS = V6cos9; 9 = -arccos(V6W); W = ij j ki
3 S
Si = SijSij; C2 = 1.9; Ok = 1.0; oe= 1.2, and the ideal gas state equation:
p = pRT.
(3)
(4)
(5)
DPM method is used to simulate the condensed phase motion, its essence is that for groups of particles emitted from one section of the input boundary
(fuel surface), motion trajectory of each group is constructed based on the calculated gas phase flow field, after that the flow field recalculation is performed, but with a consideration for forces and thermal fluxes to and from particles. This process is iteratively repeated until convergence is achieved. Trajectories are calculated by integrating the force balance equation for each group of particles [15]
dus _ ug - Us (6)
dt psdj 24
18^ Cd Re
where Ug is the gas phase velocity; us is the particle speed; is the gas viscosity; ps is the particles density; ds is the particles diameter; Re is the Reynolds number for the particle diameter; Cd is the particle drag coefficient
(Cd =a1 +—- + —3r, where the ai, a2, a3 coefficients are taken according Re Re2
to [16]).
Heat exchange between particles and gas is calculated as follows:
msCs^y- = hAs (T«- Ts), (7)
dt v 7
where As is the particle surface area; h is the heat transfer coefficient obtained from the Ranz and Marshall correlation [17]
Nu = ^ = 2 + 0.6Red/2 Pr1/3 .
X d
To improve numerical stability and accelerate convergence, the method of averaging the particles position relative to the mesh cells based on the Gaussian distribution was used, where the free distribution parameter a = 0.5 [18].
Results and discussion. Calculations showed significant nonuniformity in the particle concentration field over the nozzle cross section for the slotted channel (Fig. 2). Significant nonuniformity was also observed for the gas and particles velocity field (Figs. 3 and 4). These results are in qualitative agreement with work [19], which also notes the flow significant inhomogeneity in the channel-slot charge, but calculations in that work are based on the Euler equations. White fields near the nozzle walls indicate that concentrations of particles near the wall are tending towards zero in the considered configuration of the engine.
Velocity field of particles lag behind the gas phase at the nozzle exit is presented in Fig. 3 on a logarithmic scale. Large values of the particles lag along
Slotted channel Smooth channel
Fig. 2. Distribution of particles concentration at the nozzle exit in slotted and smooth cylindrical charge channel (di3 = 7 ^.m)
the nozzle periphery are caused by the absence of particles in this area. In such cases, the condensed phase velocity is interpreted as zero. It could be seen that significant inhomogeneity in the velocity lag is observed for the slotted channel.
Fig. 5 presents profiles of particles lagging behind the gas phase along the nozzle length. Results are in satisfactory agreement with work [1]. Difference between the Ug - uz curves from [1] and our data could be explained by the difference in the nozzle profiles.
Data were obtained on the Tz - Tg particles temperature lag. For particles with the 11 ^.m diameter, the lag is 900-1,050 K, for 3 ^.m particles, it is 100-150 K, and for 1 ^.m particles, the lag is close to zero.
Slotted channel Smooth channel
Fig. 4. Particles temperature lag field at the nozzle exit (d43 = 7 ^.m)
Fig. 5. Particles velocity lag along the nozzle length (d43 = 5 ^.m)
Calculations demonstrated that significant turbulence energy is generated in the combustion chamber axial region, and it is carried away into supersonic part of the nozzle. However, its relative value or turbulence degree (Vk / u) in the nozzle supersonic region tends to zero, which is associated with the flow significant acceleration. In the nozzle supersonic part near the wall layer, the turbulence degree value is ~ 10 % (Fig. 6).
0.01 0.00
Fig. 6. Turbulence intensity
Results of calculating the specific impulse are presented in Fig. 7 as the ratio of the ISpp calculated value to its ideal value without considering the losses. The inpid = 2.864 m/s value was obtained based on thermodynamic calcula-
тП /тП
1sp'1sp.id
Fig. 7. Influence of the condensate particles diameter on specific impulse
tion [13]. It is known from [20] that particles with the 1 ^.m diameter practically follow the gas flow. Therefore, the difference between ISp and iSpp id for d43 =
= 1 ^.m is mainly caused by friction and scattering losses. For a nozzle with a = 38°, according to [8], the scattering loss is 2.7 %. Then, according to calculations conducted, the friction losses will be ~ 0.5-0.6 %. For the investigated nozzle with d * = 200 mm and average particle size d43 = 5 ^.m, the two-phase losses will be ~ 3 %, which corresponds to the known data by other authors [4], for particles d43 = 10 ^.m, and two-phase losses will be ~ 5.5 %.
Conclusion. Based on the numerical investiagetion carried out, it could be concluded that when determining the SRM specific impulse, the need to calculate separately each type of losses and their effect on specific impulse is eliminated.
It was found out that both nozzle submerging and complicating the charge channel shape cause an increase in the specific impulse losses.
Translated by K. Zykova
REFERENCES
[1] Alemasov V.E., Dregalin A.F., Tishin A.P. Teoriya raketnykh dvigateley [Theory of rocket engines]. Moscow, Mashinostroenie Publ., 1969.
[2] Glushko VP., ed. Termodinamicheskie i teplofizicheskie svoystva produktov sgoraniya. T. 1 [Thermodynamical and thermophysical properties of combustion products. Vol. 1]. Moscow, AN SSSR Publ., 1971.
[3] Sorkin R.E., ed. Gazotermodinamika raketnykh dvigateley na tverdom toplive [Thermogasdynamics of solid fuel rocket engines]. Moscow, Nauka Publ., 1967.
[4] Timnat Y.M. Advanced chemical rocket propulsion. New York, Academic Press, 1987.
[5] Abugov D.I., Bobylev V.M. Teoriya i raschet raketnykh dvigateley tverdogo topliva [Theory and calculation of solid fuel rocket engines]. Moscow, Mashinostroenie Publ., 1987.
[6] Shishkov A.A., Panin S.D., Rumyantsev B.V. Rabochie protsessy v RDTT [Solid fuel rocket engines internal processes]. Moscow, Mashinostroenie Publ., 1989.
[7] Erokhin B.T. Teoreticheskie osnovy proektirovaniya RDTT [Theoretical foundation of solid fuel rocket engines design]. Moscow, Mashinostroenie Publ., 1982.
[8] Dobrovol'skiy M.V. Zhidkostnye raketnye dvigateli [Liquid fuel rocket engines]. Moscow, Bauman MSTU Publ., 2005.
[9] Kutateladze S.S., Leont'yev A.I. Teplomassobmen i trenie v turbulentnom pogranich-nom sloe [Heat transfer and friction in turbulent boundary layer]. Moscow, Ener-goatomizdat Publ., 1985.
[10] Avduevskiy V.S., Koshkin V.K., eds. Osnovy teploperedachi v aviatsionnoy i raketno-kosmicheskoy tekhnike [Foundations of heat transfer in airspace engineering]. Moscow, Mashinostroenie Publ., 1992.
[11] Koroteev A.S., ed. Gazodinamicheskie i teplofizicheskie protsessy v raketnykh dvigatelyakh tverdogo topliva [Gas dynamics and heat transfer processes in solid fuel rocket engines]. Moscow, Mashinostroenie Publ., 2004.
[12] Taylor J.R. Introduction to error analysis. CA, Mill Valley, University Science Books, 1982.
[13] Trusov B.G. Modelirovanie khimicheskikh i fazovykh ravnovesiy pri vysokikh tem-peraturakh [Modelling of chemical and phase equilibrium at high temperatures]. Moscow, Bauman MSTU Publ., 1991.
[14] Savel'yev S.K., Emel'yanov V.N., Benderskiy B.Ya. Eksperimental'nye metody issledo-vaniya gazodinamiki RDTT [Experimental methods of solid fuel rocket motors gas dynamics research]. St. Petersburg, Nedra Publ., 2007.
[15] Volkov K.N., Emel'yanov V.N. Techeniya gaza s chastitsami [Flow of gas with particles]. Moscow, FIZMATLIT Publ., 2008.
[16] Morsi S.A., Alexander A.J. An investigation of particle trajectories in two-phase flow systems. J. Fluid Mech., 1972, vol. 55, no. 2, pp. 193-208.
DOI: https://doi.org/10.1017/S0022112072001806
[17] Ranz W.E., Marshall W.R. Evaportaion from drops. Chem. Eng. Prog., 1952, vol. 48, no. 3, pp. 141-146.
[18] Apte S.V., Mahesh K., Lundgren T. Accounting for finite-size effects in simulations of disperse particle-laden flows. Int. J. Multiph. Flow, 2008, vol. 34, no. 3, pp. 260-271. DOI: https://doi.org/10.1016/j.ijmultiphaseflow.2007.10.005
[19] Shimada T., Daimon Y., Sekino N. Numerical simulation of flow inside a solid rocket motor by eulerian-hybrid approach with relation to nozzle inlet ablation. Proc. 8th Int. Symp. on Experimental and Computational Aerothermodynamics of Internal Flows, 2007, paper ISAIF8-00109.
[20] Soskin M.S., ed. Lazernaya anemometriya, distantsionnaya spektroskopiya i interfero-metriya [Laser anemometry, ranged spectroscopy and interferometry]. Kiev, Naukova dumka Publ., 1985.
Sultanov T.S. — Post-Graduate Student, Department of Reaction Engines and Power Plants, Kazan National Research Technical University named after A.N. Tupolev-KAI (K. Marksa ul. 10, Kazan, 420111 Russian Federation).
Glebov G.A. — Dr. Sc. (Eng.), Professor, Department of Reaction Engines and Power Plants, Kazan National Research Technical University named after A.N. Tupolev-KAI (K. Marksa ul. 10, Kazan, 420111 Russian Federation).
Please cite this article as:
Sultanov T.S., Glebov G.A. Numerical computation of specific impulse and internal flow parameters in solid fuel rocket motors with two-phase combustion products. Herald of the Bauman Moscow State Technical University, Series Mechanical Engineering, 2021, no. 3 (138), pp. 98-107. DOI: https://doi.org/10.18698/0236-3941-2021-3-98-107
