Aerothermal analysis of Themis T3: influence of nozzle cluster geometry on the base heating and heat loads on the grid fins

AEROTHERMAL ANALYSIS OF THEMIS T3: INFLUENCE OF NOZZLE CLUSTER GEOMETRY ON THE BASE HEATING AND HEAT LOADS ON THE GRID FINS M. Laureti * and S. Karl * * DLR - Institute of Aerodynamics and Flow Technology - Spacecraft Bunsenstr. 10, 37073 G¨ ottingen, Germany [email protected] - [email protected] ABSTRACT The Horizon Europe project SALTO (reusable Strategic Space Launcher Technologies & Operations) supports ESA’s Themis program by testing key technologies for reusable launchers. In this framework aerothermal CFD (Computational Fluid Dynamic) analysis have been performed to assess the heat loads acting on the surface of Themis T3: a reference configuration for a future demonstrator, during its entire flight trajectory. This paper shows representative flow field solutions and heating patterns for critical phases of the ascent and re-entry flight with a focus on the base area and aerodynamic control surfaces. The impact of nozzle cluster geometry and gas generator exhaust, used to feed the turbopump and subsequently injected in dedicated nozzles, on the base heating has been also evaluated. Index Terms—Reusable Launch Vehicle, Retro-Propulsion, Aerothermal Loads, CFD, RANS 1. INTRODUCTION As space exploration progresses, the focus is on increasing efficiency, sustainability, and reducing costs, leading to the development of reusable space launch systems. This transition is driven by the need for more environmentally responsible and economically viable access to space, with reusable launch vehicles playing a central role in meeting these needs. Within Europe, a confluence of ambitious programs and initiatives, led by the EU, government agencies and private companies, is driving the advancement of Vertical Take-off Vertical Landing (VTVL) Reusable Launch Vehicles (RLVs). The European knowledge and expertise in that field is growing thanks to systematic studies and detailed investigations of system components [1, 2], a notable example is the Horizon 2020 RETALT project [3]. A valuable contribution to the research and development of reusable launchers has been made also by CALLISTO (Cooperative Action Leading to Launcher Innovation in Stage Toss back Operations) project, born from the collaboration between DLR, CNES and JAXA [4, 5]. Through initiatives such as the Future Launchers Preparatory Programme (FLPP), ESA aims at developing the necessary know-how for future launch vehicles to be able to guarantee Europe’s needs for a secured and autonomous access to space. Part of the strategy consists in promoting the development of cutting-edge technologies and concepts for next-generation launch systems. The RETPRO (Validation of Wind Tunnel Test and CFD Techniques for Retro-propulsion) project, part of ESA’s FLPP, aims to prepare and validate the necessary tools for reliable design and simulation of a descent trajectory of a proposed in-line VTVL launcher configuration using retro-propulsion and aerodynamic control surface to perform re-entry and landing maneuvers [6]. The Horizon Europe project SALTO (reusable Strategic Space Launcher Technologies & Operations) directly supports ESA’s Themis program by testing key technologies for reusable launchers. Themis is an ESA rocket prototype with a height of c.a. 29 meters and a diameter of 3.5 meters, it can be used as a first stage rocket, a modular launcher or a liquid booster [7]. Themis is powered by the Prometheus engine, which uses liquid methane and oxygen as propellants. The SALTO project involves multiple fly/recover/re-fly cycles. Hop tests, where the rocket will take off, fly, and land vertically are planned for the last quarter of 2025 and will take place at the Esrange Space Center in Sweden. As part of the SALTO project, a reference configuration for a future demonstrator, called Themis T3, designed to perform suborbital flights, is also studied. This configuration is equipped with three Prometheus engines and is supposed to reach an altitude of 100 km, typical for a rocket first-stage mission. Current aerothermal studies in the literature have focused on 9-engine cluster configurations [8, 9, 10], with a layout similar to Space X’s FALCON 9, or single-engine configurations such as the Callisto demonstrator [11]. Themis T3 differs from the mentioned layouts by the presence of 3 engines. The number of engines forming the cluster is an important parameter to consider when evaluating the thermal loads acting on the region base. In fact, the geometric configuration strongly influences the topology of the flow field that develops downstream or upstream (depending on the direction of the flight) of the nozzles exit section. The size and shape of the recirculation region that typically develops around the engine plumes, as well as the intensity of plume-to-plume interactions, are strongly affected by the free area on the base plate in relation to its diameter. Furthermore the gas generator exhaust, used to feed the turbopump and subsequently injected in dedicated nozzles, plays a role on the baseplate heating. This work shows the results of the aerothermal CFD (Computational Fluid Dynamic) analyses of Themis T3 configuration, performed in the framework of the SALTO project, focusing on the typical flow field phenomena occurring during the flight trajectory (ascent and re-entry) especially in the plumes region and the resulting heating pattern on the rocket rear-part. Heat loads on other important structural parts like aerodynamic control surfaces, i.e. grid fins, and landing legs will be presented as well. 2. MATHEMATICAL & NUMERICAL MODEL The heat flux data at different flight regimes and surface temperatures, used for creating aero-thermal data bases, are provided by an extensive CFD simulations campaign. These CFD analyses are performed with the hybrid structuredunstructured DLR Navier–Stokes solver TAU [12].The TAU code is a second order finite-volume solver for the Euler and Navier-Stokes equations which includes a comprehensive range of RANS-based or scale resolving turbulence models and extensions for chemically reacting high enthalpy flows. It uses hybrid structured/unstructured computational grids to facilitate the analysis of complex geometries and is highly optimized for the application on massively parallel HPC systems. TAU has been successfully applied to a wide range of subto hypersonic flow problems, both in scientific and industrial applications, including the analysis of re-usable launcher configurations. The baseline set of numerical models which have been applied for the present investigations provides accurate and robust treatment for all flow conditions as defined in the CFD test matrix on a wide range of Mach numbers. The calculation of the inviscid fluxes in the finite volume framework is based on the application of the AUSMDV flux vector splitting scheme [13], together with MUSCL gradient reconstruction to achieve second order spatial accuracy [14]. Viscous fluxes are treated with a low-dissipation central discretization scheme. Turbulence was modeled with the Spalart-Allmaras oneequation RANS model in the low Reynolds formulation. This model provides a good compromise between numerical efficiency and accuracy and is particularly applicable to flows with strong shocks. It has been shown in previous studies [15] that the RANS plume structure of a Spalart-Allmaras-Model compares well to LES reference data. The model completely resolves the structure of the turbulent boundary layer including the laminar sub-layer. The vehicle surfaces are treated as isothermal walls with a constant temperature of 200K and 600K. Interpolation algorithms, implemented in the aerothermal database, allow to evaluate the heat flux at each point of the vehicle surface as a function of the surface temperature. Nozzle exhaust of Prometheus engine, where applicable, is modeled by imposing the exit flow conditions of the thrust nozzles at their respective exit planes in the computational domain of the launcher. The radial distribution of fluid dynamic variables, as well as the exhaust chemical composition, is determined by separate 2D-axisymmetric simulations of the flow inside the thrust nozzle, which provides an accurate representation of the exhaust jets. The LOX/(L)CH4 skeletal reaction mechanism of Slavinskaya et al.[16] is used to evaluate chemical non-equilibrium effects. The thermo-chemical modeling is generally based on a reacting mixture of thermally perfect gases and the solution of a transport equation for each component. The properties of the individual species are computed by either spectroscopic constants using partition functions that include an accurate representation of high temperature effects such as anharmoniccorrections and coupling of rotational and vibrational degrees of freedoms for molecules [17] or from NASA-Polynomials [18]. Appropriate mixture rules are applied to compute thermodynamic properties depending on the local gas composition, pressure and density. Species viscosities are calculated through Blottner or Sutherland curve fits. Chemical non-equilibrium is assumed to account for the additional heat release due to post-combustion of fuel-rich exhaust gases in the flow field around the rocket configuration and endothermic dissociation in the stagnation region. The global 3-step mechanism of Westbrook and Dryer is applied [19]. Thermal conductivity is estimated from the Euckencorrection. Wilke ´ s rule gives the thermal conductivity of the mixture. Turbulent heat conductivity is obtained from a constant ratio of laminar to turbulent Prandtl numbers of 0.8. The turbulent diffusivity is calculated using a constant turbulent Schmidt number of 0.7. 3. THEMIS T3 ASCENT AND RE-ENTRY CONFIGURATION The significant size of the rocket results in large Re-numbers during flight phases with high dynamic pressure associated with thin boundary layers, which need thin prismatic mesh layers in the computational grid. During high altitude flight and retro-maneuvers, boundary layers are much thicker. Furthermore, the thickness of the boundary layer depends on the wall temperature. Hence grids with different boundary layer resolutions are generated to ensure accurate resolution and computational efficiency. Generally, the prismatic sub-layers close to the wall are constructed with a first dimensionless wall spacing which ensures a y+ in the order of one, and a wall normal stretching ratio of grid cells lower than 1.3. As- suming a zero angle of attack, the use of a quarter computational domain allows for significant savings in computational resources. Surface meshing used a hybrid approach. Unstructured triangular meshes are applied to complex structures such as landing legs or grid fins. Smooth regions on the surface are meshed with structured quadrilateral meshes which slightly improves the accuracy of the flow solver due to the availability of line detection algorithms for gradient reconstruction. The overall sizes of the computational grids are around 2.4M and 16.6M points for the ascent and re-entry configuration, respectively. Surface meshes are shown in Fig.1. Ducts are not included in the ascent configuration, while they are retained in the re-entry one. The resolution of the flow field around and inside the grid fin requires a large number of grid points, which is the main reason for the size difference between the two configurations. (a) Ascent configuration (b) Descent configuration Fig. 1: Surface meshes 4. CFD REPRESENTATIVE RESULTS Due to the sensitive nature of the data, both inputs and outputs are presented in dimensionless form. Fluid dynamic variables are made dimensionless using the free stream conditions, denoted by the subscript ∞. Altitude and time are made dimensionless by dividing by the maximum altitude and the time in which the vehicle reaches the maximum altitude. The pressure distribution on the vehicle is given by the pressure coefficient, cp, calculated according to Eq.1. cp=p−p∞ 1 2ρ∞u∞ 2(1) The heat fluxes are scaled by dividing by a reference value corresponding to the maximum heat flux occurring at the nose of the launcher during the ascent phase. 4.1. Ascent phase This study shows that the heat loads characterizing the ascent phase are substantially lower than those on reentry. During the ascent phase it is possible to identify two characteristic base flow modes depending on altitude and strength of plumes interaction. Fig.2 and Fig.3 show the flow field topology associated to low and high altitude mode, respectively. Fig.2a and Fig.3a illustrate the pressure coefficient distribution on the surface. The Mach number on the symmetry plane is plotted in purpleorange scale. The green iso-surface represents the boundary of 50% nitrogen mass fraction and indicates the geometrical extent of the engine plumes. Fig.2b and Fig.3b show the scaled heat flux on the surface and the dimensionless temperature field on the symmetry plane. This time the white iso-surface depicts the region of flow reversal and recirculation with axial velocity u = 100 m/s. The base flow topology for the two characteristic flow modes is shown in Fig.2c and Fig.3c. Contrary to what was observed in other projects [9, 20], where significant base heating occurs only during high altitude flight and in small regions between the nozzle clusters, here the highest thermal loads on the base plate occur during the low altitude mode and involve a large area. At low altitudes due to the large atmospheric ambient pressure, exhaust plumes are confined and the intensity of plume-plume interactions is relatively limited. The heat flux is dominated by flow recirculation in the near field, as can be seen in Fig. Fig.2c. The forebody turbulent boundary layer separates at the base plate shoulder, hits the hot nozzle exhaust and brings heated gas back in the bubble. The heat flux streaks at the baseplate shoulder are placed in a region where the flow goes from the bottom to the top. In the low altitude mode, due to the high dynamic pressure and density, significant thermal loads can be observed downstream of the shock forming in front of the folded grid fins and on the fin surfaces perpendicular to the flow direction, see Fig.3b (bottom left frame). (a) Mach number on the symmetry plane, cpsurface distribution and iso-surface 50% N2 mass fraction (b) T/T∞on the symmetry plane, q/qnose max surface distribution and iso-surface of flow reversal (c) Base flow topology: T/T∞on the symmetry plane, q/qnose max surface distribution and streamlines Fig. 2: Ascent flow field structure: low altitude mode, h/hmax= 0.114, Mach =1.1) (a) Mach number on the symmetry plane, cpsurface distribution and iso-surface 50% N2 mass fraction (b) T/T∞on the symmetry plane, q/qnose max surface distribution and iso-surface of flow reversal (c) Base flow topology: T/T∞on the symmetry plane, q/qnose max surface distribution and streamlines Fig. 3: Ascent flow field structure: high altitude mode, h/hmax=0.57 , Mach =2.8) At higher altitudes (h/hmax >0.25) the ambient pressure is lower and consequently exhaust plumes significantly spread out. Stronger plume-plume interactions occur in the vicinity of the launcher base leading to a recirculation of hot exhaust gases between the nozzle clusters. The exhaust back flow leaves a mark on the heat flux pattern on the base plate. Remarkable changes in the base flow topology can be noticed by comparing Fig.2c and Fig.Fig.3c, the most relevant consists in the absence of the recirculation bubble generated by the detachment of the turbulent boundary layer at the base plate shoulder. The temporal evolution of the average heat flux on the most thermally stressed components is now analyzed and the effect of the presence of gas generators is evaluated. As expected, the time history of the average heat flux on the nose is not affected by the presence of gas generators, as shown in Fig.4a. In the time window from take-off to engines shutdown (time/timeh max = 0.56), heating is dominated by the difference between wall temperature and ambient temperature. After engines shutdown, heat fluxes gradually become zero. (a) Nose (b) Baseplate (c) Nozzle wall Fig. 4: Ascent time history of average heat flux (scaled by the reference value) on different components) Regarding the baseplate and the nozzle walls, the discussion is much more complex. Unlike in the RETALT project [9], the transition between low and high altitude mode is gradual, so there is no sharp jump in the average heat flux. At low altitude, the heat flux on the base plate is dominated by flow recirculation in the base region (time/timeh max <0.42), see Fig.4b. As the altitude increases, plume to plume interactions occur as an additional driving phenomenon. Moreover, it is necessary to consider the deep impact provided by the presence of gas generator plumes. Fig.4b shows that the presence of gas generators does not significantly affect the average heat flux on the baseplate calculated for Tw=200K, while for Tw=600K the difference is significant. Gas generator plumes have a large effect on the average heat flux observed on the nozzle walls, and this is true for both the considered wall temperatures, see Fig.4c. Again, after the engines’ shutdown, average heat fluxes become zero. Fig. 5: Peak heat flux position For components where the heat peak always occurs at the same point, see Fig.5, it is possible to monitor its evolution over time. In the case of the folded grid fin, the maximum value is at the center of the surface belonging to the folding mechanism and perpendicular to the flow. Regarding the central nozzle, in the case without gas generators, the maximum heat flux occurs on the symmetry plane close the exit. In the case with gas generators, the point is always on the symmetry plane but not far from the base plate. (a) Folded grid fin (b) Central nozzle wall Fig. 6: Ascent time history of peak heat flux on different components A complex system of interacting shocks and expansions is generated upstream of the folded grid fins. The location where the heat peak occurs does not vary with time however very high values are found at the corners of the fin surfaces perpendicular to the flow. Just like the nose, the heating is influenced by the difference between wall temperature and ambient temperature, the maximum value occurs at time/timeh max = 0.42, as shown in Fig.6a. The presence of the gas generators obviously does not affect the flow field around the grid fins. The time evolution of the peak heat flux on the surface of the central nozzle is qualitatively similar to the mean value. When the exhaust of the gas generators is considered, the peak value is significantly higher than in the case without gas generators. This effect is particularly pronounced at low altitudes and for a wall temperature Tw=200K, see Fig.6b. (a) Gas generators off (b) Gas generator on Fig. 7: Low altitude mode: heat flux distribution obtained with and without gas generators (a) Gas generators off (b) Gas generator on Fig. 8: High altitude mode: heat flux distribution obtained with and without gas generators The presence of gas generators significantly changes the topology of the flow field in the base region and consequently the heat flux pattern on the base plate and nozzle walls as can be seen in Fig.7 and Fig.8. One of the reasons behind these changes lies in the interaction of the gas generator plumes with the recirculation bubble that forms at low altitude close to the base plate. The presence of the gas generators results in the additional introduction of gases with momentum and different chemical composition that can potentially give rise to local post-combustion effects (it is necessary to consider the presence of oxygen from the free stream in the recirculation bubble). Furthermore, the plumes of gas generators impinge on the nozzle walls and interact with the main plumes. 4.2. Descent phase Aerothermal load predictions for Vertical Take-Off Vertical Landing concepts impose a comprehensive set of challenges which are related to the complex flow structure during retroburn maneuvers. A wide range of complex gas dynamic and aerothermal effects occurs. They require adequate modeling and are difficult to identically replicate in wind tunnel experiments. Most of these effects are related to the interaction of the propulsive jets with the supersonic free stream. The flow field topology at the beginning of retro-propulsion is shown in Fig.9. In particular, Fig.9a shows the Mach field on the symmetry planes (purple-green scale), the sonic line (solid black line) and the heat flux on the rocket surface (rainbow scale). Fig.9b depicts the temperature field on the symmetry planes (yellow-red scale), streamlines colored with the axial velocity and the heat flux distribution on the vehicle structure (rainbow scale). Finally, Fig. 9c shows the pressure field on the symmetry planes and corresponding iso-lines (blue-yellow scale) and the pressure coefficient on the surface (blue-red scale). The iso-surface of 30% N2 mass fraction represents the geometrical extent of the plume. Note that, again, the heat flux and fluid dynamic variables are presented in dimensionless form. This flow field topology has never been observed before and is the result of the choice to start the retro-propulsion at relatively low altitude and to use only the lateral engines. Fig.9a shows that along the centerline the free stream is decelerated by the bow shock. The flow downstream of the bow shock is subsonic, slightly reaccelerates toward the central nozzle, but is stagnated by a second shock that forms in front of it. The engine jet is laterally recompressed by a barrel shock and a contact surface separates free stream air from engine exhaust. Due to the high pressure on the axis of symmetry and in front of the central nozzle, the engine plume is bent outward as shown in Fig.9c, therefore engine plumes do not merge. Instead of the classical recirculation zone that surrounds the plume and shields the launcher, Fig.9b highlights the presence of a vortex relatively far from the launcher surface that provides no protection. In fact, it is possible to observe streamlines lapping the landing leg cover in the positive y-direction and enveloping the rocket from the outside inward. The impact of streamlines on the base plate leaves a trace on the mechanical and thermal loads. The high temperatures of the mixture of exhaust gases and air present between the nozzles and the low flow velocity cause high heat flux values on the nozzles wall. (a) Mach number and sonic line on the symmetry plane, q/qnose max surface distribution (b) T/T∞on the symmetry plane, q/qnose max surface distribution, streamlines colored according to u/u∞ (c) p/p∞with iso-lines on the symmetry plane, cpsurface distribution and iso-surface 30% N2 mass fraction Fig. 9: Re-entry flow field structure: beginning of retropropulsion, h/hmax= 0.168, Mach=3) During the retro-propulsion phase, especially at the beginning of it, large heat loads are observed at the leading edge of the grid fins, as well as their support structures, which are directly exposed to the hot exhaust gases. The values found on the leading edge of the grid fin are in line with those found in RETPRO project [16] for similar free stream conditions, although the launcher configuration is different. Except for the leading edge where q/qnose max = 62, the average value on the grid fin remains relatively low (q/qnose max ≈10). High peaks are also found on the grid fin support (q/qnose max = 43.5). The heat flux distribution on the grid fin is given in Fig.10a. The Mach field and temperature iso-lines on a slice taken at the central section of the grid fin are shown in Fig.10b, a cut-off at Mach=1.9 has been applied. On the same slice, the divergence of velocity field normalized by its maximum value and the sonic line are displayed in Fig.10c. This quantity is extremely useful for identifying waves present in the flow field. The divergence of velocity is negative in the presence of isentropic compressions or shocks and positive in the case of expansions. (a) Heat flux surface distribution (b) Mach number contour, temperature iso-lines and heat flux surface distribution (c) Divergence of velocity on the slice, heat flux surface distribution and sonic line (d) N2+O2 mass fraction with iso-lines on the slice and heat flux surface distribution Fig. 10: Grid fin analxsis: beginning of retro-propulsion, h/h˙max= 0.168, Mach =3) As shown Fig.10a, the heat flux distribution on the grid fin is not uniform, in fact the lower part undergoes higher thermal loads and this can be explained observing the flow field in Fig.10b and Fig.10c. The lower part is in the wake of the landing leg cover and protrusions. The flow is subsonic (see sonic line in Fig.10c) and characterized by high temperatures. An isentropic compression occurs in front of the grid fin (total pressure is conserved). Moving away from the base the temperature decreases. At the top of the grid fin the temperature is lower by a factor of 2. The negative value of the divergence of velocity highlights the formation of a shock in front of the grid fin (upper part), the sonic line shows that the grid fin operates in chocked mode under the current flight conditions. In addition to thermal loads, the chemical composition of the atmosphere in which grid fins are immersed is an important parameter in the choice of the manufacturing material. This allows to account for wall chemical reactions and possible deposition effects. Since the engine jets do not contain N2, the mass fraction of molecular nitrogen can be used to determine in which regions of the field and to what extent the exhaust gases are present. The sum of N2 and O2 mass fractions in a slice passing through the center of the grid fin indicates that the atmosphere in which it is immersed is predominantly air, see Fig.10d. The time evolution during the re-entry phase of average dimensionless heat fluxes occurring on the most thermally stressed components for the two analyzed wall temperatures is shown in Fig.11. (a) Twall = 200 K (b) Twall = 600 K Fig. 11: Re-entry time history of average heat flux on different components Regarding the T=200K case the convective heat fluxes are