Numerical approach of a hybrid rocket engine behaviour

(1)

DEGREE PROJECT REPORT

Numerical approach of a hybrid rocket engine behaviour

Modelling the liquid oxidizer injection using a Lagrangian solver

Gustave S

PORSCHILL

Master’s Program in Aerospace Engineering

Examiner: Elena GÛTIERREZ FÂREWIK Supervisors: Stefan WÂLLIN

Jean-Yves LESTRADE

Multi-Physics for Energetics Department

Master’s thesis carried out from 06/02/2017 to 07/07/2017

UNCLASSIFIED

(2)

(3)

Abstract

To access and operate in space, a wide range of propulsion systems has been developed, from high-thrust chemical propulsion to low-thrust electrical propulsion, and new kind of systems are considered, such as solar sails and nuclear propulsion. Recently, interest in hybrid rocket engines has been renewed due to their attractive features (safe, cheap, flexible) and they are now investigated and developed by research laboratories such as ONERA.

This master’s thesis work is in line with their development at ONERA and aims at finding a methodology to study numerically the liquid oxidizer injection using a Lagrangian solver for the liquid phase. For this reason, it first introduces a model for liquid atomiser developed for aeronautical applications, the FIMUR model, and then focuses on its application to a hybrid rocket engine configuration.

The FIMUR model and the Sparte solver have proven to work fine with high mass flow rates on coarse grids. The rocket engine simulations have pointed out the need of an initialisation of the flow field. The methodology study has proven that starting with a reduced liquid mass flow rate is preferable to a simulation with a reduced relaxation between the coupled solvers. The former could not be brought to conclusion due to lack of time but gives an encouraging path to further investigate.

Key-words: Hybrid propulsion, combustion, numerical simulations, two-phase flow, La- grangian approach

(4)

(5)

v

List of Figures

1.1 Ariane 5 at launch (Credit: ESA – CNES) . . . 5

1.2 Liquid propulsion . . . 6

1.3 Solid propulsion . . . 6

1.4 Hybrid propulsion . . . 6

1.5 Studied configuration . . . 7

1.6 Diffusion flame inside a hybrid rocket . . . 7

1.7 Relation between Isp and oxidizer-to-fuel ratio O/F [3] . . . . 9

2.1 Break-up of a jet flow [27] . . . 11

2.2 Primary break-up regimes [22] . . . 13

2.3 Regime domains [24] . . . 13

2.4 Secondary break-up regimes [20] . . . 14

2.5 Jet flow atomisation with DNS [15] . . . 14

2.6 Comparison between ELSA model and DNS (2D cut) [15] . . . 15

2.7 Particle velocity in a hollow-cone atomisation using the FIMUR model (2D cut) . . 15

3.1 Geometry of the atomiser modelled by the FIMUR model [25, 26] . . . 20

3.2 The vertical position shift ∆y of the atomiser . . . 21

4.1 The test case mesh . . . 25

4.2 Gas and droplet temperature at t = 0.1 s . . . 26

4.3 Evaporation of H2O2 in case 1 without reactions at t = 0.1 s . . . 27

4.4 Mass fraction distribution in case 2 with H2O2 decomposition reaction at t = 0.1 s . 27 4.5 Profiles in both test cases at x = 0.15 m . . . 28

(7)

vii

5.1 The rocket engine geometries . . . 29

5.2 Mesh details . . . 31

5.3 Pressure waves due to combustion . . . 31

5.4 Gas injection with a turbulence model . . . 32

5.5 Liquid injection with reverse gas flow . . . 34

5.6 Liquid injection with the narrow atomiser . . . 35

5.7 Liquid injection with the narrow atomiser and the nozzle . . . 36

5.8 Working liquid injection with ˙m = 30 g s⁻¹. . . 38

A.1 Density for H2O2 and H2O . . . 44

A.2 Vapour pressure for H2O2 and H2O . . . 45

A.3 Heat of vaporization for H2O2 and H2O . . . 46

A.4 Heat Capacity for H2O2 and H2O . . . 47

A.5 LiquidViscosity for H2O2 and H2O . . . 48

A.6 Thermal conductivity for H2O2 and H2O . . . 49

A.7 Surface tension for H2O2and H2O . . . 50

B.1 Delavan WDA atomiser specification sheet (extract) . . . 51

List of Tables

1.1 Possible propellants for hybrid propulsion . . . 8

3.1 The atomiser parameters . . . 22

5.1 Dimensions of the geometries . . . 30

A.1 Constant physical properties . . . 43

A.2 Sparte data for density . . . 44

A.3 Sparte data for vapour pressure . . . 45

A.4 Sparte data for heat of vaporization . . . 46

A.5 Sparte data for heat capacity . . . 47

A.6 Sparte data for liquid viscosity . . . 48

A.7 Sparte data for thermal conductivity . . . 49

A.8 Sparte data for surface tension . . . 50

(8)

viii

List of Symbols

Variables

d Droplet diameter [m]

F Thrust [N]

Isp Specific impulse [m s⁻¹] or [s]

k Rate of reaction [mol m⁻³s⁻¹]

m Mass [kg]

˙m Mass flow rate [kg s⁻¹]

O/F Oxidizer-to-fuel ratio [–]

Oh Ohnesorge number [–]

R Radius [m]

Re_l Liquid Reynolds number [–]

T Temperature [K]

U Velocity [m s⁻¹]

U_e Effective exhaust velocity [m s⁻¹]

u, v, w Velocity axial, radial and tangential components [m s⁻¹]

W e Weber number [–]

X Air core ratio [–]

[X] Molar concentration of species X [mol m⁻³]

Y_X Mass fraction of species X [–]

µ Dynamic viscosity [Pa s]

ρ Density [kg m⁻³]

σ Surface Tension [N m⁻¹]

θS Half spray-angle [°]

Constants

g0 Earth’s gravitational acceleration at sea level [m²s⁻¹]

R Universal gas constant [J mol⁻¹K⁻¹]

(9)

ix

Subscripts

0 At the atomiser’s nozzle

g Gas phase

l Liquid phase

Acronyms

CEDRE CFD software for Energetics developed by ONERA CFD Computational Fluid Dynamics

DNS Direct Numerical Simulation DPS Discrete Particle Simulation

ELSA Eulerian-Lagrangian Spray Atomiser

FIMUR Fuel Injection Model by Upstream Reconstruction PIV Particle Image Velocimetry

(10)

(11)

1

Acknowledgement

First of all, I would like to thank Jérôme Anthoine, head of the Propulsion Laboratory unit at ONERA, for letting me carry out my master’s thesis in his unit and giving me the opportunity to witness the rocket propulsion research environment.

I am grateful to Jean-Yves Lestrade, my supervisor, and Jérôme Messineo, both researchers of the unit, for their feedbacks and suggestions to progress in this thesis work.

I also thank Jouke Hijlkema and Olivier Rouzaud, senior users and former developers of Sparte, for their help with this solver, and more broadly Jean-Mathieu Senoner, head deve- loper of Sparte, and the CEDRE software support team for their unlocking explanations and solutions.

Finally, I would like to thank Jean-Étienne Durand and Quentin Levard, PhD students, Duc Minh Le, postdoc researcher, and everyone else at the Propulsion Laboratory unit for welcoming me among them and the nice atmosphere they provided.

(12)

(13)

3

Introduction

Nowadays, space has become an important scientific, economic and strategic location. A crowd of satellites orbits the Earth to monitor its atmosphere and oceans, predict the weather, broadcast data such as global positioning to users all over the world, and observe closely war zones. During the last century, space propulsion has been developed to access these orbits and even space beyond, through chemical propulsion such as liquid and solid rocket engines and electrical propulsion. As improvements of the existing systems are investigated, new solutions like solar sails and nuclear propulsion are considered.

Hybrid chemical propulsion, using both liquid and solid propellants, has recently drawn back the attention on it with Virgin Galactic’s SpaceShipOne successful flight. Especially interested in its safety and flexibility features, aerospace laboratories have decided to speed up its development. ONERA, the French aerospace laboratory, studies numerically and experimentally various types of hybrid rocket engines in its Propulsion Laboratory unit, part of the Multi- Physics for Energetics Department on the site of Le Fauga-Mauzac.

Carried out there, this master’s thesis work investigates the methodology to follow in order to run simulations of a hybrid rocket engine with a liquid oxidizer injection, using a Lagrangian solver and the FIMUR model for the atomiser. First, the context of space propulsion is introduced, highlighting the characteristics of the hybrid propulsion system, and the liquid injection is presented. Then the tools and global models for the simulations are described.

Finally, the FIMUR model is tested, followed by the analysis of the methodology and of the difficulties encountered on the rocket geometry.

(14)

(15)

5

Chapter 1

Context of the Master’s Thesis

1.1 Space propulsion

Figure 1.1: Ariane 5 at launch (Credit: ESA – CNES) Space propulsion systems basic goal is to accele-

rate a spacecraft. Since the spacecraft is isolated when in space, its motion can still be altered by ejecting mass according to the action-reaction principle. Derived from Newton’s law of motion, the thrust vector F applied on the spacecraft is:

F = mdU

dt = − ˙m Ue (1.1) where Ue is the velocity vector of the expelled mass. The propulsion system performances are then evaluated with the specific impulse Isp:

Isp= F

˙m ≡ U_e [m s⁻¹] (1.2) It represents the propulsion efficiency and corresponds to the effective exhaust velocity of the propellant. A high Isp is desirable, since it means low propellant consumption. Another widely used definition is given in (1.3), which can be interpreted

as the length of time for which one kilogram of propellant produces a thrust equivalent to a one-kilogram mass (i.e. a force of about 9.81 N) in Earth’s gravitational field g0.

I_sp = F

˙m g0 [s] (1.3)

Space missions are of various types, such as launching payloads into orbit (e.g. Ariane 5 in Figure 1.1) or orienting a spacecraft on its orbit, so the propulsion systems are varied as well. More particularly, for high-thrust applications, chemical propulsion systems are divided into three categories:

– Liquid propulsion (Figure 1.2): the propellants are stored separately as liquid in tanks.

They are mixed in the combustion chamber before being accelerated through the nozzle.

(16)

6 CHAPTER 1. CONTEXT OF THE MASTER’S THESIS Liquid rockets are characterized by the highest Ispfor chemical propulsion (up to 465 s for the Vinci engine developed by ASL [7]).

– Solid propulsion (Figure 1.3): both oxidizer and fuel are premixed together in one solid phase propellant, the grain, usually with a metallic additive (e.g. aluminium) to increase the thrust performances. Producing the highest thrust amongst the chemical rockets (14 · 10⁶N for the Space Shuttle SRB) but with low Isp, it is mainly used on launcher boosters at take-off. Its design is quite simple but it is impossible to throttle or to shut down. Moreover, the solid propellant is explosive and highly sensitive to manufacturing imperfections.

– Hybrid propulsion (Figure 1.4): the propellants are stored in two different states, a solid one and a liquid one, generally a solid fuel with a liquid oxidizer. Its theoretical performances lie between the solid and liquid propulsion performances, with a vacuum Isp

that can exceed 360 s.

Figure 1.2: Liquid propulsion Figure 1.3: Solid propulsion

Figure 1.4: Hybrid propulsion

Out of these three chemical propulsion systems, the first two have been widely studied and developed due to fast breakthroughs and improvements of their performances. The interest in the last one has been recently renewed thanks to Virgin Galactic’s SpaceShipOne, but also thanks to its promises: safer, cheaper, more flexible.

1.2 Hybrid propulsion main concepts

The hybrid propulsion can be designed according to different configurations. The traditional one, studied here and shown in Figure 1.5, consists in a cylindrical combustion chamber, also called a port, inside a tubular solid fuel, as in solid propulsion systems. To this port are added a pre-combustion chamber, to insure the heating and the vaporisation of the oxidizer before combustion, and a post-combustion chamber, to increase the combustion advancement before

(17)

1.2. HYBRID PROPULSION MAIN CONCEPTS 7 exiting the gas through the nozzle. The oxidizer is either injected directly as liquid inside the rocket engine, or through a catalyst that heats and vaporises the oxidizer before it enters the pre-chamber.

Pre-combustion chamber

Fuel grain Port

Post-combustion chamber

Nozzle Liquid

Oxidizer

Figure 1.5: Studied configuration

One of the main differences between the hybrid propulsion and the two other chemical propulsions presented earlier lies in the combustion. Indeed, in solid and liquid propulsions, the fuel and the oxidizer are closely mixed in an almost homogeneous phase in the combustion chamber, as fine solid grains or as small injected droplets. This leads to a premixed flame for the solid rocket and a short diffusion flame with a large flame area for the liquid rocket, both with a controlled oxidizer-to-fuel ratio set to the stoichiometric ratio.

In hybrid propulsion, the gaseous fuel is produced at the surface of the solid fuel by pyrolysis (i.e. decomposition of the solid phase directly into fuel gas due to high heat fluxes), while the oxidizer is injected at the entrance of the combustion chamber. This results in a long diffusion flame inside the boundary layer, represented in Figure 1.6, where the propellants must diffuse towards one another in order to burn. The diffusion flame is less efficient here since the fuel is injected on the whole length of the combustion chamber: the fuel blown from the grain end has less time to correctly burn, and the mixing area between fuel and oxidizer is smaller than the one in the liquid rocket, where the flame appears all around each droplet. Adding the post- combustion chamber enlarges the residence time of the propellants and the flame area, thereby helping gain in performance. Besides, the oxidizer-to-fuel ratio is not constant and decreases along the port due the consumption of the oxidizer.

Figure 1.6: Diffusion flame inside a hybrid rocket

Concerning the choice of propellants, several combinations have been investigated. Table 1.1 gives examples of fuels and oxidizers used in hybrid rockets. HTPB, HDPE and PP correspond to common rubbers, while paraffins are liquefiable waxes that will produce a thin liquid layer between the solid grain and the gaseous reactants.

(18)

8 CHAPTER 1. CONTEXT OF THE MASTER’S THESIS

Table 1.1: Possible propellants for hybrid propulsion

Fuels Oxidizers

Hydroxyl-terminated Polybutadiene (HTPB) Liquid Oxygen (LOx) High Density Polyethylene (HDPE) Hydrogen Peroxide (H2O2)

Polypropylene (PP) Nitrogen Peroxide (N2O4)

Paraffins Nitrous Oxide (N2O)

Cellulose (Wood) Nitric Acid (NO3H)

. . . .

1.3 Hybrid propulsion characteristics

As mentioned in Section 1.1, interest in hybrid rockets is related to its many advantages [4]:

(+) Safety: The propellants are separated by distance and by phase, with an inert fuel such as rubber, so that a spark cannot light the engine on. This allows for safe manufacturing, safe transportation, safe assembly and safe storage. The launchpad operations are also safer and the risks of catastrophic failure are reduced, while the hybrid rocket can be shut-down (cf. Flexibility). The working pressures, between 10 bar to 50 bar, are lower than in solid and liquid propulsion, from 50 bar up to 150 bar.

(+) Reliability: The mechanical system is less complex than the ones in liquid propellant rockets, since it uses half of the storage and feed system. In the same time, the solid fuel is insensitive to cracks and imperfections, contrary to solid oxidizer-fuel grains in solid rockets.

(+) Flexibility: The hybrid rocket can be throttled, shut-down and restarted, by adjusting the oxidizer flow inside the combustion chamber. As seen in Table 1.1, hybrid propulsion also offers a large choice of propellants. Besides, it can be designed for large boosters as well as for satellite manoeuvre thrusters.

(+) Environmental friendliness: The hybrid rocket has an environmentally clean exhaust, without hydrogen chloride or aluminium oxide, which impact the environment. Actually, for LOx/HTPB or H2O2/HDPE rockets, exhausts are only composed of H2O, CO, CO2

and of unreacted propellants.

(+) Low cost: The safety of the materials and the simplicity of the systems reduce the cost.

However, improvements have still to be investigated to correct or compensate some drawbacks:

(–) Slow regression rate¹: Compared to solid propulsion, the fuel grain regression rate is low, making it difficult to obtain sufficient mass flow rate of pyrolysed fuel vapour in order to achieve high thrust levels. But recent studies have shown possible improvements through innovative fuels and injection methods [4, 18].

(–) Mixing and combustion inefficiencies: The diffusion flame inside the hybrid rocket combustion chamber is less efficient than the flames found in the other propulsion systems.

Moreover, the fuel is “injected” all along the combustion chamber, thus a part of it might leave the rocket unburned.

1The regression rate corresponds to the thickness of solid fuel that is vaporised per unit time.

(19)

1.3. HYBRID PROPULSION CHARACTERISTICS 9 (–) Low volumetric loading: To enhance the combustion efficiency, a post-combustion chamber is added on hybrid rockets, decreasing further the low volumetric loading due to slow regression rates.

(–) Oxidizer-to-fuel ratio shift: The regression of solid fuel leads to an increasing combus- tion chamber diameter and thus a shift in the mixing ratio over time. Noting that the ideal ratio is the stoichiometric ratio, this shift directly impacts the performances of the hybrid rocket, as shown in Figure 1.7. Therefore, current studies investigate solid fuel geometry optimisation.

Figure 1.7: Relation between I_sp and oxidizer-to-fuel ratio O/F [3]

(20)

(21)

11

Chapter 2

Liquid injection

2.1 Two-phase flow

Liquid injection in a combustion chamber lets two phases interface inside the flow, the liquid droplets and the carrying gas stream. This corresponds to one of the different regimes of multiphase flows, classified by Ishii [13] as dispersed flow.

Actually, close to the injector, the considered flow consists in a jet flow, i.e. a separated flow without mixing of the two phases. Different phenomena occur then: atomisation, which breaks up the jet flow into droplets and thus into the dispersed flow, heating and evaporation of these droplets, but also collisions between droplets and with the wall, or further fragmentation.

Figure 2.1 shows a jet flow going through atomisation, resulting in a dispersed flow.

Figure 2.1: Break-up of a jet flow [27]

2.2 Atomisation

Atomisation is an important process that transforms a liquid jet into a diluted spray of small droplets. Its required characteristics depend on the application, which can be chemistry, agriculture or motor design.

In rocket propulsion, the injector vaporizes the propellant to enable the combustion between gas. The atomisation phenomenon has to produce the smallest drops possible so that the interface between the two phases is increased, leading to a fast evaporation and a reduced combustion chamber in length. For a classical hybrid rocket configuration, the fine atomisation

(22)

12 CHAPTER 2. LIQUID INJECTION influences the length of the pre-combustion chamber. Indeed, in an optimized pre-combustion chamber the residence time and the evaporation time of a droplet are equal.

The process can be split into two steps. The primary atomisation consists in the first break-up. The high differential velocity between the jet and the gas produces an important shear stress on the phases interface, leading to instabilities that detach fragments from the jet.

The primary atomisation depends on the considered liquid streams, since the exact mechanisms may differ between cases (e.g. liquid jet, liquid sheet, swirl jet,...).

The secondary atomisation corresponds to the second break-up: large drops into a spray of droplets. Carried downstream by the gas and in transverse directions due to turbulence, the large drops of the primary break-up are destabilized and break into smaller drops. This step is mostly controlled by the ratio between the destabilizing aerodynamic forces and the stabilizing surface tension forces acting on a liquid fragment, which is expressed through the Weber number Weas:

We= ρg(Ug− U_l)²d

σ_l (2.1)

The destabilization regimes can also be characterized by other non-dimensional numbers taking in account the properties of the gas, of the liquid and of the geometry. For example, the liquid Reynolds Rel corresponds to the ratio between aerodynamic and viscosity forces in the liquid drop, while the Ohnesorge number Oh is the ratio between viscosity and surface tension forces acting on a drop:

Re_l= ρ_lU_ld

µ_l (2.2)

Oh= √µl

ρ_lσ_ld (2.3)

2.2.1 Primary atomisation

The primary atomisation is an important process, since it gives the secondary atomisation its initial conditions, i.e. the large drops and liquid fragments distribution. Therefore, the choice of the injector is also important. There are three main types of injectors:

– Pressure atomiser: the liquid is injected at high velocity, through a small nozzle, thanks to a high pressure difference between the flows inside and outside the injector. This is the one used in the studied case, and usually used inside propulsion systems.

– Rotary atomiser: the liquid is injected radially from a spinning surface such as a disc. The centrifugal energy helps achieve high relative velocity between the liquid stream and the carrying gas, leading to primary break-up.

– Twin-fluid atomiser: the liquid is injected at low velocity inside a coaxial high velocity gas stream, usually air. This injector, also called airblast or air-assist atomiser, uses the kinetic energy of the airstream to shatter the liquid stream.

Depending on the flow parameters through the previously given non-dimensional numbers We, Rel and Oh, experiments have shown distinctions between several break-up regimes, which have been classified by Reitz [22] into four main regimes, represented in Figures 2.2 and 2.3:

– the Rayleigh mechanism: at low We and Rel, instabilities in the jet stream are due to aerodynamic forces from the gas relative flow, resulting in drops at least as large as the nozzle exit diameter, far from the nozzle.

(23)

2.2. ATOMISATION 13 – the first wind induced break-up: the aerodynamics forces apply a torsion to the liquid

jet stream, leading to drop of the size of the injector nozzle diameter, far from it.

– the second wind induced break-up: surface tension instabilities produced by the aerodynamic forces are added to the torsion of the previous case. The drops are smaller than the nozzle diameter.

– the atomisation: for large We or Rel, the jet stream is shattered as soon as it exits the injector nozzle into really small droplets compared to the nozzle diameter.

Figure 2.2: Primary break-up regimes [22] Figure 2.3: Regime domains [24]

For the liquid injection in rocket propulsion, as well as in all propulsion applications, the last regime is the required one. The droplets being smaller will let them evaporate faster and mix more efficiently with the gaseous fuel in the combustion chamber.

2.2.2 Secondary atomisation

Following the primary atomisation, the liquid fragments and drops are broken up into smaller droplets, due to shear stress resulting amongst others from their velocity differential with respect to the carrying gas flow. The secondary break-up is thus mainly characterised by the Weber number, based on the primary drops and fragments length scale. Once again, several regimes have been identified [8, 23], starting from a critical Weber number Wecusually set to 12. They are shown in Figure 2.4.

– for Wec ≤ We ≤ 100, the bag break-up regime occurs. In this first regime, the drop flattens and grows hollow near the stagnation point, under the effect of the dynamic pressure. It starts to break up at the bottom of the “bag”. On the second half of the We-range starts a transitional regime: the bag and stamen break-up.

– for 100 ≤ We ≤ 350, the shear break-up (or sheet stripping) regime occurs. Unlike the previous regime, the flattened drop is deformed into a ligament and broken up starting from the edges. This phenomenon could be explained by destabilizing capillary waves or by the stripping of the drop boundary layer due to shear stress.

– for 350 ≤ We, the catastrophic break-up regime occurs. The large relative velocity between the drop and the carrying gas of the last regime leads to small-wavelength

(24)

14 CHAPTER 2. LIQUID INJECTION perturbations at the surface of the drop. The resulting small ligaments on the drop then break up into droplets.

Figure 2.4: Secondary break-up regimes [20]

2.2.3 Numerical description

To model accurately the whole atomisation process is quite difficult. It can be done using a Direct Numerical Simulation (DNS) as proposes Lebas et al. [15], shown in Figure 2.5.

This method corresponds to a numerical experiment and yields detailed results that replace experimental data, hard to obtain because the sensors are not adapted to the dense spray of the primary break-up. It is however limited to short simulation times due to today’s limited processing performances.

Figure 2.5: Jet flow atomisation with DNS [15]

Another way is to use a complex model developed for the atomisation process, such as the ELSA model (Eulerian-Lagrangian Spray atomisation) based on the work of Vallet and Borghi [28]. It describes the dense zone using an eulerian approach of a two-phase fluid, i.e. the fluid has a variable density and two species. The model also calculates the mean area of the liquid- gas interface, in order to switch to Lagrangian approach for the liquid phase (see below) in the dilute zone. The good accuracy in the near nozzle region comes with some drawbacks, amongst which the limitation to a 3D mesh for a correct description [12]. Figure 2.6 compares the fields calculated with the ELSA model to DNS.

(25)

2.2. ATOMISATION 15

(a) Liquid volume fraction (b) Liquid-gas surface density

Figure 2.6: Comparison between ELSA model and DNS (2D cut) [15]

In the present study, to avoid this still heavy modelling, the FIMUR model (Fuel Injection Model by Upstream Reconstruction) introduced by Sanjosé [25, 26, 2, 11] will be used. This model, developed for spray simulations in aircraft motors, directly produces the droplets distribution of a swirl pressure atomiser after the primary and secondary break-ups, based on the atomiser and the spray characteristics. As shown in Figure 2.7, the resulting spray forms a hollow cone, due to the swirl motion given to the injected liquid. The model is detailed in Section 3.2.

Figure 2.7: Particle velocity in a hollow-cone atomisation using the FIMUR model (2D cut)

Finally, to model the dispersed phase, the large number of particles (liquid in the current case) makes it difficult to reproduce them numerically. Therefore, three different formulations for CFD software have been developed [19]:

– DPS (Discrete Particle Simulation): this is the direct approach, where all the physical particles are tracked individually. The same constraints as for the DNS apply: calculations are extremely heavy and can only be carried out on small parts during short simulation times.

– Lagrangian approach: simplification of the DPS, it solves for numerical particles, also called parcels, containing several “real” particles with the same characteristics (size, velocity,...) and located by their center of gravity. Quite accurate, it still can be heavy

(26)

16 CHAPTER 2. LIQUID INJECTION in calculation depending on the number of parcels to compute and is difficult to speed up through parallel processing. It has to be coupled with an eulerian solver for the gas flow.

The spray distribution in Figure 2.7 illustrates the Lagrangian approach.

– Eulerian approach: similar to gas solvers, it solves the flow fields based on conservative equations, generally on a finite volume grid. More suitable for implicit algorithms and parallel processing, it is confronted to numerical diffusion and precision issues since the dispersed phase is modelled by average per cell. Moreover, it cannot explicitly take into account some particle phenomena, such as collisions.

2.3 Evaporation

In the considered hybrid rocket configuration, the oxidizer is injected as liquid droplets at storage temperature (approximatively 300 K) inside the combustion chamber. There, the droplets are heated and evaporated.

While the liquid droplet is cold, only little vapour is present at its surface, which leads to low evaporation rate and mass transfer. Therefore, the received heat is mainly used to increase the droplet temperature. This increase produces vapour, decreasing the heat flux at the liquid surface and slowing down the temperature increase. The temperature inside the droplet becomes uniform and the evaporation is accelerated.

The evaporation rate depends on the carrying gas (pressure, temperature and properties), on the liquid phase (temperature, diameter and properties) and on the relative velocity of the droplet with respect to the carrying flow.

The evaporation is generally not modelled globally for the whole dispersed phase, but independently for each droplet. Three main models have been developed [1]:

– the D² model: the simplest model, the temperature is supposed uniform in the droplet and constant over time. All the received heat is used to evaporate, with a law proportional to the square of the drop diameter.

– the infinite conductivity: the conductivity of the droplet is supposed infinite, so the temperature is uniform in the liquid, but it varies with respect to time.

– the effective conductivity: in this model, the heating of the droplet is more detailed, leading to a non-uniform temperature inside the droplet.

Using one of these models implies some simplification hypotheses, among which:

– the process is quasi-stationary

– the considered drop is spherical and isolated – the liquid only contains one species

– the phase change is way faster than the vapour transport in the ambient air (the vapour is produced at the surface temperature of the drop and its corresponding vapour pressure) – the heat transfer by radiation is negligible

2.4 Objectives of the thesis

The objectives of the current master’s thesis is to find a methodology to get simulations of liquid H2O2injection inside a hybrid rocket to run, using a Lagrangian description.

(27)

2.4. OBJECTIVES OF THE THESIS 17 Such simulations have already been carried out, especially by Lazzarin et al. [14] with N2O. At ONERA, simulations of hybrid rocket engines have been first investigated in 2009 by Lestrade [16], and highlighted some difficulties with the Lagrangian solver and its evaporation models. Until now, these simulations have therefore been performed with catalysed injection, i.e. the oxidizer enters the pre-chamber already as a hot gas. This work focuses back on the liquid injection, as the solvers and their models have been further developed and improved.

When Lestrade introduced the liquid spray inside the rocket engine, no atomiser models were implemented in the software. The droplets were thus injected from a large distribution of point injectors, each having unvarying characteristics, to produce the injection spray based on experimental data. Since then, the FIMUR model was developed and should simplify the injection settings while yielding a better droplets distribution in the spray. As the atomiser model has been implemented for aeronautic propulsion systems, the thesis investigates its application with the Lagrangian solver to rocket propulsion systems with much larger mass flow rates (100 g s⁻¹ instead of usually around 2 g s⁻¹).

(28)

(29)

19

Chapter 3

General settings

3.1 CEDRE

All simulations have been carried out on CEDRE¹, the CFD software developed by ONERA [21]. It is dedicated to numerical simulations for energetics and propulsion, and is redistributed internally but also externally to aerospace industrials. Contrary to commercial CFD codes, CEDRE is meant for research applications and focuses therefore more on the correctness of the physics and less on the robustness, by including regularly new models developed by researchers and PhD students.

The software is divided into several solvers interacting with one another in order to consider the different physical phenomena occurring in this kind of flow, for example reactive compressible flow, but also dispersed flows, thermal conduction in the walls or heat radiation.

In the present study, only the Lagrangian solver for dispersed flow Sparte is used in addition to Charme, the main solver for reactive compressible flows. The gas phase is modelled as a thermal perfect gas while the heat radiation is neglected, though it can play a significant role in rocket combustion chambers, depending on the chemical species involved and on the gas conditions (temperature, pressure). The input grid is designed with Gmsh, an open source 3D finite element grid generator [9].

Both solvers Charme and Sparte can be used in steady or unsteady flow. In the case of Charme, a steady simulation corresponds to a simulation using local time step size, which can be specified through several models. When no local time step model is chosen, the simulation is unsteady. For Sparte, the steady simulation consists in finding the trajectory of an injected particle until either it evaporates or exits the domain. However, this kind of simulation is reserved for dilute sprays where the particles do not have important influence on the carrying gas. Therefore, when introducing droplets, both solvers are set to unsteady simulation and the convergence is assessed by looking to the variation of parameters such as the pressure, the temperature and the mass fractions.

Concerning evaporation, the Sparte solver proposes only two of the presented models in Section 2.3: infinite conductivity and effective conductivity. To keep the model simple and fast to solve, the first one is used in all simulations.

The interaction between both solvers is managed by CEDRE as follow:

– One-Way coupling: Only Charme affects Sparte. The gas is solved as if no particles were injected, while the particles follow the gas flow and feed on its energy. It corresponds

1Calcul d’Écoulements Diphasiques Réactifs pour l’Énergétique: Two-phase Reactive Flow Calculation for Energetics, version 6.1.1

(30)

20 CHAPTER 3. GENERAL SETTINGS to a correct approximation for highly dilute dispersed flows in terms of volume fraction, where the influence of the particles on the carrying gas is negligible (e.g. particles in PIV).

– Two-Way coupling: Charme and Sparte affect each other. The Lagrangian solver adds source terms in the gas flow, corresponding to the exchanges in mass, momentum and energy between the particles and the gas. To avoid large gradients when starting the liquid injection, a relaxation parameter moderates the interaction of Sparte on Charme and can be set between 0 (i.e. One-Way coupling) and 1 (i.e. 100% Two-Way coupling).

The software do not include the Four-Way coupling, which also solves the interactions between each particle in very dense dispersed flows.

3.2 FIMUR model

3.2.1 General description

As indicated in Section 2.2.3, the FIMUR model defines the distribution of the dispersed phase after the atomisation of a pressure swirl atomiser (shown in Figure 3.1), to avoid its complex modelling. It also greatly simplifies previous “manual” technique described in Section 2.4.

Figure 3.1: Geometry of the atomiser modelled by the FIMUR model [25, 26]

Given the input parameters for the atomiser and the spray, the model yields the droplet distribution for position, velocity and diameter downstream the atomisation process, and re- constructs the path of the droplets upstream to the atomiser nozzle to insure continuity of the injection. The atomiser is characterized by the mass flow-rate of the injected liquid ˙m, the orifice radius R0and the air core radius Rarelated to the air core ratio X, defined in Equation (3.1), and the spray by its half spray-angle θS and the characteristic diameter of the droplet distribution dafter atomisation.

X = A_a

A₀ =R_a R₀

2

= sin²θ_S

1 + cos²θ_S (3.1)

3.2.2 Implementation in SPARTE

In addition to the parameters mentioned above, the Lagrangian solver requires the position, the orientation and the injection period of the atomiser, which determines the number of new numerical particles injected per time step.

(31)

3.2. FIMUR MODEL 21 The position of the atomiser and the time step for the Sparte solver are constraints by the inputs parameters. Indeed, for each time step, the atomiser model places the new numerical particles inside a cone, determined by the axial velocity of the particles and the half spray-angle, as shown in Figure 3.2. If this cone goes over the limits of the mesh, some particles (in red in Figure 3.2a) will be statistically placed outside the mesh and lead to errors.

θ_S

Δy = ul Δttanθ_S ul Δt

(a) Bad placement

θ_S

Δy = ul ΔttanθS

ul Δt

(b) Correct placement

Figure 3.2: The vertical position shift ∆y of the atomiser

Due to that, the atomiser cannot be placed on the symmetry axis in a 2D axisymmetric mesh but have to be slightly shifted inside, and the closer from the axis it is, the smaller the time step must be, thus increasing the processing time.

To avoid the shift, one could think of inclining the injection cone of θS, while reducing the half spray-angle to the spray half thickness angle δθS, to reproduce the same injection. However, the velocity distribution (3.2) inside the spray is defined according to the injection axis: this will therefore lead to a wrong distribution.











u(x0, r, φ) = ˙m ρlπR₀²(1 − X) v(x0, r, φ) = 0 w(x0, r, φ) = ˙m

ρ_lAp

R0+ Ra

2r

(3.2)

3.2.3 The atomiser in the simulations

Using Sparte and FIMUR, the modelled injector is a swirl pressure atomiser, which parameters are given in Table 3.1. These parameters are based on the Delavan WDA nozzle, used in the experimental settings and presented in Appendix B.

This corresponds to a dense hollow cone injection, with a distribution of the particles within

±5° of the half spray-angle. The initial temperature of the droplets is set to 300 K. The droplet diameter distribution is given by a log-normal distribution, defined by its probability density function p(d) in (3.3):

p(d) = Log-N (ln(dm), σ²) = √ 1

2πdσexp −(ln(d) − ln(dm))² 2σ²

!

(3.3)

where σ = 0.4 is the standard deviation.

(32)

22 CHAPTER 3. GENERAL SETTINGS

Table 3.1: The atomiser parameters

Parameter Value

Mass flow rate ˙m 100 g s⁻¹ Orifice radius R0 0.89 mm

Air core ratio X 0.14 Half spray-angle θS 30°

Droplet mean diameter dm 50 µm

According to Equation (3.2), the droplets at the injector nozzle have an axial velocity of 32.4 m s⁻¹. This leads to a vertical position shift of the atomiser of about 5 · 10⁻⁵m using a time step ∆t = 10⁻⁶s (as it will be used in the following chapters). Considering a margin and noting that the first cell height at the injection point is of 2 · 10⁻⁴m in the rocket grids (see Section 5.1), the chosen shift is ∆y = 10⁻⁴s. This offset works fine and no particles are lost.

3.3 Chemistry

In the considered hybrid rocket design, the propellants are stored as liquid hydrogen peroxide H2O2 and solid HDPE (high-density polyethelene). In the combustion chamber, due to high temperature and heat fluxes, the H2O2 evaporates and decomposes itself into O2 and H2O following Reaction (3.4), while the HDPE yields mainly C2H4 by pyrolysis.

H2O2 ka

H2O + 1

2O2 (3.4)

No catalyst is used, so that the H2O2 decomposition is thermally driven. The kinetics is given by a first order reaction rate according to Giguère and Liu [10], modelled with an Arrhenius equation:

k_a= Aa exp−E_a,a RT

[H2O2] (3.5)

where R = 8.314 J mol⁻¹K⁻¹ is the universal gas constant and

A_a= 1 · 10¹³s⁻¹ E_a,a= 200 832 J mol⁻¹

The hydrogen peroxide is naturally unstable and is difficult to sustain above a concentration of 99% in an aqueous solution. Practical and economical considerations therefore limit its concentration to 98% in industrial applications [5]. This leads to a two-step evaporation, the water evaporating at a lower temperature than H2O2 for a given pressure (cf. Appendix A.2).

Although the Sparte solver manages multi-species particles, it cannot yet evaporate each species independently. For that reason, and to keep the model simple, the hydrogen peroxide is assumed to be pure in the following simulations. For simplification still, the fuel pyrolysis is assumed to yield only gaseous C2H4.

The combustion occurs between the gaseous fuel C2H4 and the decomposition product O2, which is the real oxidizer of the combustion. It is described according to the two-step reaction model (3.6), using the coefficients established by Westbrook and Dryer [29].

C2H4+ 2 O2 kb

2 CO + 2 H2O CO + 1

2O2 kc

k−c CO2

(3.6)

(33)

3.3. CHEMISTRY 23

The rates of reaction are modelled with the following Arrhenius equations:

k_b = Ab exp−E_a,b RT

[C2H4]^0.1[O2]^1.65 (3.7)

k_c = Ac exp−E_a,c RT

[CO][H2O]^0.5[O2]^0.25 (3.8)

k−c= A−c exp−E_a,−c RT

[CO2] (3.9)

where

A_b = 7.589 · 10⁷m^2.25mol^−0.75s⁻¹ E_a,b = 125 520 J mol⁻¹ A_c = 1.259 · 10¹⁰m^2.25mol^−0.75s⁻¹ E_a,c = 167 360 J mol⁻¹ A−c = 5 · 10⁸s⁻¹ E_a,−c= 167 360 J mol⁻¹

Finally, the propellants are considered to be in stoichiometric proportions, corresponding to the optimum oxidizer-to-fuel ratio O/F = 2.5 for the couple O2/C2H4 [4]. The H2O2 mass flow rate is set to 100 g s⁻¹, i.e. an injected O2 mass flow rate of 47 g s⁻¹, this leads to a fuel mass flow rate of 18.8 g s⁻¹.

(34)

(35)

25

Chapter 4

FIMUR test cases

4.1 Test settings

Several test cases have been carried out to investigate the general behaviour of the model with high mass flow rate and evaporation.

For these quick simulations, a simple 2D-axisymmetric coarse grid (0.2 × 0.15 m) has been used with the time step ∆t = 10⁻⁶s. It is shown in Figure 4.1:

X Y

AXIS Z

WALL

INLET

OUTLET

ATOMISER

Figure 4.1: The test case mesh

It contains 2000 elements and its boundary conditions are defined as follows:

– on the left stands an adiabatic wall, – below is the symmetry axis,

– on the top, the boundary condition is set to a pressure condition, acting as inlet as default, – on the right, an outlet with fixed pressure,

– the atomiser is placed in the refined lower left corner.

After some test cases with the default species, i.e. air as a species for the gas phase and water for the liquid phase, an input file for the Sparte solver, defining the liquid hydrogen peroxide properties, has been written (cf. Appendix A) and two cases have been tested:

(36)

26 CHAPTER 4. FIMUR TEST CASES

1. The liquid H2O2 simply evaporates into the gas phase without any chemical reactions, 2. The liquid H2O2evaporates into the gas phase and the gaseous H2O2decomposes according

to Reaction (3.4).

The grid is initialized with water vapour at 30 bar, 2000 K. The atomiser is placed on the axis and injects liquid hydrogen peroxide at 300 K with the mass flow rate ˙m = 100 g s⁻¹.

4.2 Results and discussion

In the steady state, the injection cone is deformed by the flow it has created and tends to spread toward the axis. This flow also cools down the inner part of the cone to 300 K by convecting downstream the heat transfer with the droplets exiting the atomiser (see Figure 4.2).

(a) Case 1 (b) Case 2

Figure 4.2: Gas and droplet temperature at t= 0.1 s

As it can be seen in Figure 4.3a for case 1, the smallest droplets are driven closest to the axis due to their low inertia and end up in the cold region Therefore they cannot evaporate.

The heaviest droplets are little affected by the flow on the grid length scale and go through the hot gas. They need more time to heat, so they evaporate little in the grid. Finally, the gaseous H2O2forms in-between these two groups (Figure 4.3b).

In case2, the H2O2 forms a line between the cold, where the droplets do not evaporate, and the hot regions, where it decomposes quickly and completely (see Figure 4.4). Reaction (3.4) is exothermic and thus releases heat, which affects the temperature profile and causes the H2O2

to decompose faster. Consequently, a sharp increase in temperature parallel to a sharp drop in peroxide mass fraction can be observed in Figure 4.5. Also, this sustained high temperature in the spray increases the evaporation rate of the droplet compared with the first case : in Figure 4.2b, the droplets almost disappear inside the light green zone, whereas they can be seen down to the cold border in Figure 4.2a.

(37)

4.2. RESULTS AND DISCUSSION 27

(a) Droplets diameter and streamlines (b) Gas phase H2O2mass fraction

Figure 4.3: Evaporation of H₂O₂ in case 1 without reactions at t= 0.1 s

(a) Gas phase H2O2 (b) O2

Figure 4.4: Mass fraction distribution in case 2 with H₂O₂decomposition reaction at t= 0.1 s

At this point, no specific issues are encountered with the FIMUR model, or more broadly with Sparte, despite the dense spray and the importance of the mass flow rate, and can therefore be applied to a more complex geometry. However, the resulting spray will have to be compared to coming experimental measurements to insure the accuracy of the model. If necessary, the FIMUR model could be fitted to the experimental observations in the range of the rocket propulsion conditions (high mass flow rate, high pressure, high temperature).

(38)

28 CHAPTER 4. FIMUR TEST CASES

(a) Temperature (b) Mass fractions

Figure 4.5: Profiles in both test cases at x= 0.15 m

(39)

29

Chapter 5

Simulations on the rocket geometries

5.1 Geometry

The simulations are processed on two different 2D-axisymmetric geometries. The first one (Figure 5.1a) represents only the pre-chamber and the port. It is a simplified case of the hybrid rocket, containing only the parts of interest. Indeed, the focus is put on the atomiser and on the region close to it.

1

3 2

5

4 x y

L_prechamber L_port

R_prechamber R_port

(a) Without nozzle

1

3 2

5

x y

Lprechamber Lport

R_prechamber Rport

Lnozzle

2

4

(b) With nozzle

Figure 5.1: The rocket engine geometries

The second geometry (Figure 5.1b) adds a short post-chamber and a nozzle, to represent the whole internal flow from the injection to the nozzle throat. This sets more correctly the outlet condition and solves some issues explained in the following sections, but it also adds complexity and new difficulties. The throat radius R^∗ has been calculated using (5.1), derived from the mass flow rate definition with the isentropic relations.

A^∗= π(R^∗)²= ˙m P0

s RT₀ Mgasγ

γ+ 1 2

^(γ+1)/2(γ−1)

(5.1) where P0 and T0 are the pressure and temperature conditions inside the combustion chamber.

Mgas and γ are the molecular mass and the heat capacities ratio for the mixed gas, and are

(40)

30 CHAPTER 5. SIMULATIONS ON THE ROCKET GEOMETRIES based on the exit gas of simulations on the first geometry.

The dimensions are given in Table 5.1. They are chosen approximatively according to usual dimensions and do not represent an experimental set-up. Therefore they are not optimized yet for perfect working conditions.

Table 5.1: Dimensions of the geometries Part Length Radius

[mm] [mm]

Pre-chamber 50 37.5

Port 230 12.5

Nozzle 75 3.8¹

On these geometries, fives boundary conditions are defined:

1. The inlet (in green on the figures). During initialisation (cf. Section 5.2) and until the injection is fully liquid, this boundary injects decomposed H2O2, as if catalysed upstream.

Hence it injects O2 and H2O at 1225 K, with the mass fraction YO₂ = 0.47. The radius of the inlet is chosen at the beginning of the simulation, but can be changed between two different cases. Once the atomiser injects the whole oxidizer in liquid phase, the inlet is turned off and becomes an adiabatic wall.

2. Adiabatic walls (in black).

3. The fuel pyrolysis (in red). It corresponds to a wall that injects C2H4 at ˙mfuel = 18.8 g s⁻¹ when its temperature exceeds 900 K.

4. The outlet (in blue). In the first geometry, it is a simple subsonic outlet set to 30 bar.

In the second geometry, it is set to supersonic outlet once the nozzle throat is choked. In order to reach that condition, the outlet is first set to a low pressure condition (≈ 1 bar).

5. The symmetry axis (in yellow).

Due to the simplicity of the first geometry, the mesh can be regularly structured with rectangle elements and refinement close to the boundaries. The corner at the port entry is also refined, to better capture the gradients where the injection flow separates into a recirculating flow inside the pre-chamber and the port flow (Figure 5.2a).

The second geometry takes the same structured grid in the pre-chamber and the port.

However, to avoid too many cells at the nozzle throat, the post-chamber is partially unstructured, as shown in Figure 5.2b. The post-chamber and nozzle regions are not of interest in this study and are only used here to force the exit condition. Therefore, the slight irregularities at the borders between the structured and the unstructured regions are not important here.

The resulting grids consist in between 30 000 and 40 000 elements. These numbers are based on previous works at ONERA on a similar geometry to obtain a sufficient precision, enough here to define a methodology. Nevertheless, a grid refinement study should be carried out to insure the accuracy of the processed result fields.

1This corresponds to the throat radius, the pre-chamber and post-chamber radii being equal.

(41)

5.2. INITIALISATION 31

X Y

(a) The injection and port entrance mesh Z

X Z Y

(b) The nozzle mesh

Figure 5.2: Mesh details

5.2 Initialisation

Before introducing liquid H2O2droplets inside the grid, an initialisation of the flow must be carried out.

Indeed, a direct simulation with the particles leads to divergence in less than 20 iterations.

The cause for that is the sudden and important gradient in momentum at the injection point of the pre-chamber. Another downside with this method, the cold droplets will not be able to start the combustion. In a real case, the combustion would actually be started with a pyrotechnic igniter.

Therefore, the rocket is first ignited with a gas phase injection, corresponding to the case of a hybrid rocket with catalysed injection. The objective here is to obtain a converged velocity field, with a stable flame. Once this stage is reached, the injection is to be progressively transferred to liquid droplets.

The time step size is driven by the chemistry. The chemical characteristic time is hard to assess, because it depends on many local variables, such as the local concentration or the local temperature. The time step has been reduced until the disappearance of the numerical instabilities due to the combustion, which provoked amplifying pressure waves with high gradients at the entrance edge of the port (Figure 5.3). The resulting value is ∆t = 10⁻⁶s, corresponding to a maximum CF L number below 0.5 in the mesh.

Figure 5.3: Pressure waves due to combustion

(42)

32 CHAPTER 5. SIMULATIONS ON THE ROCKET GEOMETRIES

5.3 Catalysed gas injection

The injected gas density being three order of magnitude lower than the liquid phase, the inlet diameter is chosen equal to the port diameter. The resulting injection velocity for ˙m = 100 g s⁻¹, a bit less than 30 m s⁻¹, is acceptable regarding to the liquid injection velocity.

The flow inside the rocket engine is turbulent, so a turbulence model is added to yield results closer to reality. However, the input turbulent parameters are hard to estimate, as no measurements inside the combustion chamber are available. The chosen input values are thus theoretical and let the numerical convergence correct them.

The turbulence is modelled with the two-equation k-ω SST model, based on Messineo’s work [17]. In the input file, CEDRE allows the user to enter the two scalars needed as turbulence level T u and turbulent length scale `, and converts them automatically to k and ω. For a cylindrical combustion chamber, these values are usually given by (5.2):

T u ≈0.05 and ` ≈0.07D = 0.0018 m (5.2) The gas-initialised solution displayed in Figure 5.4 does not exactly correspond to the steady state. The O2is still convected inside the pre-chamber recirculation zone and will carry on until the mass fraction is homogeneous inside of it, but its field in the port is converged. Moreover, the velocity field has converged and the diffusion flame is well defined and stable, which corresponds to the goals of this initialisation. Liquid injection can therefore start progressively.

(a) Mass fraction of O₂and streamlines

(b) Gas temperature, showing the diffusion flame

Figure 5.4: Gas injection with a turbulence model

(43)

5.4. STARTING THE LIQUID INJECTION 33

5.4 Starting the liquid injection

The liquid particles are added in Two-Way coupling, using the FIMUR atomiser described in Section 3.2. Starting directly with a full coupling leads to fast divergence, due to the large source terms in the dense region. Especially two source terms are in cause.

The first one is the momentum source term. In this case, the calculation states that there is a cell with negative density or negative pressure in the near region of the atomiser. Due to the transfer of momentum, the gas is highly accelerated in the first cells of the mesh and might end up empty of gas. This is actually one of the drawbacks of the FIMUR model: the CFD code does not take into account the separate phase flow in the vicinity of the atomiser.

The second source term causing divergence is the energy exchange. This can be explained by the coupling mechanism. Each time step is divided into two sub-steps for each solver. First, the Charme solver processes the gas phase based on the previous time step solution. At the end of this sub-step, each cell possesses a determined amount of energy in form of heat. Then the Sparte solver moves the particles inside the gas and calculates the new source terms to send to the gas solver. Especially, if the gas temperature in a cell is higher than the temperature of the particles that it contains, heat is transferred from the cell to the particles. However, each particle absorbs an amount of heat as if it was alone inside the cell, while the gas of the cell is not refreshed and keeps its heat until the end of the time step. Indeed, the coupling does not take into account interactions between particles, whatever their proximity, and the evaporation model considers the droplets as isolated (Section 2.3). Also, note that the particles are located with their center of gravity, without consideration of their spatial extend: a droplet astride several cells only produces source terms in one of them. As a result, in the case of high evaporation rate of a dense spray, more heat is absorbed by the particles than available inside the cell, causing its absolute temperature to become negative. Since such a behaviour is not physical, the calculation stops.

Three approaches were considered to solve this issue. The first one consists in using larger cells, and thereby reducing the accuracy of the simulation, so that each cell contains more gas and thus more heat to exchange for a given time step. The coupling actually works fine on a coarser grid, as it has been demonstrated with the FIMUR test cases (Chapter 4), but in the long run the objectives are to obtain a good representation of the flow inside the combustion chamber, so the cells will have to be kept fine. Besides, having coarse cells close to the injector and fine elsewhere affects the grid quality as well as the accuracy of the flow field in the pre-chamber.

The time step can also be reduced to let the evaporation occurs more progressively. The time step was therefore set to 10⁻⁷s, because reducing further the time step increases harshly the computation time.

The last method consists in reducing the influence of the particles on the gas by reducing the relaxation parameter. This parameter needs then to be progressively increased as the flow field converges towards its steady solution. To get the energy issue corrected and the simulation run, the relaxation parameter had to be reduced below 5% on the whole grid, in a first phase. To that, a mask in the vicinity of the atomiser was added. This mask allows to choose a different relaxation parameter inside a rectangular area, but only one can be defined and its limits are sharp, leading to a jump in the relaxation parameter value. In the current cases, it was set below 0.1%, i.e. almost One-Way coupling, to hide the atomisation process region.

Numerical approach of a hybrid rocket engine behaviour

DEGREE PROJECT REPORT