B Fromentin-Denozière

H-K Phan

A Cheraly

A Guy

D Gueyffier

Ch Erades

J Labaune

New developments in numerical modelling of the radar signature of rocket exhaust plumes

Keywords: plume, radar signature, rocket,LES, turbulence

In this paper, we present the new developments of the multiphysics numerical models designed for predicting the time-dependent electromagnetic signature of rocket engine plumes. We first describe the physical and numerical model for aerothermochemical physics and for electromagnetic applications. A mock up rocket engine (MoSER) has been designed and its signature has been measured in UHF and X band. The results were compared to simulation.

I. INTRODUCTION

The plume flow that occurs behind a solid propellant rocket engine undergoes postcombustion leading to the formation of a flame. As the propellant contains alkali impurities (K, Na) with low ionization potentials, electrons are release in the postcombustion flame, leading to an interaction with a given incoming radar wave. Hence, it is of practical interest to predict the radar signature of a rocket engine and this has been studied in literature. More accurate information is given by the frequency spectrum of the reflected radar wave. This spectrum is usually modeled using basic turbulence models, hence relying on rather crude assumptions. In this work, we present an original approach to study the role of the large eddies on the spectral signature of the engine.

The numerical approach results from the coupling of two models: i-An aerothermochemical numerical code called CEDRE [START_REF] Refloch | CEDRE software[END_REF] which simulates the fluid flow and the dynamics of the charged species in the rocket plume using LES(large eddy simulation ) model for the turbulence.

ii-An electromagnetic numerical code called Maxwell3D [START_REF] Levadoux | Recent advances in the pre-conditioning of integral equations of electromagnetism[END_REF] [START_REF] Levadoux | New trends in the preconditioning of integral equations of electromagnetism[END_REF] which uses the electron density and the electronneutral collision frequency calculated by CEDRE.

Earlier work has been conducted using a 2D-axisymetric approach [START_REF] Fromentin-Denoziere | Numerical modelling of the radar signature of rocket exhaust plumes[END_REF][2] using a RANS turbulence model. The new development concern a more precise model of turbulence in order to emphasize the role of the large eddies on the spectral signature of the engine.

A mock up of solid propellant motor, MoSER, has been realized. Its signature has been measured in the UHF band and the frequency spectrum has been compared to the simulated one.

II. MODEL DESCRIPTION

The fluid mixture in the plume is assumed to be weakly ionized; we also assume that for the radar frequencies under consideration the electromagnetic field interacts only with the plume electrons and not with the heavier ions.

A. Simulation of ionized rocket plumes

The method exhibited in this paper for calculating the radar dynamic radar signature relies on the 3D unsteady flowfields of the plume undergoing postcombustion. In this section we briefly present the various models used to simulate the plume flow, the transport of charged species and the chemical reactions acting upon them.

Equilibrium composition in the chamber. Electrons and ions are formed due to high temperatures in the rocket combustion chamber. These charged species result from the thermal ionization of small amounts of alkali metals contained as impurities in the propellant. In order to obtain the density of these species in the chamber we use the chemical equilibrium code EQUIL contained in the CHEMKIN II package.

Aerothermochemical simulations. Simulations are carried out using CEDRE, a numerical tool able to simulate supersonic reactive gas flows coupled to a disperse phase [START_REF] Fromentin-Denoziere | Numerical modelling of the radar signature of rocket exhaust plumes[END_REF][START_REF] Gueyffier | Numerical Simulation of Ionized Rocket Plumes[END_REF][START_REF] Refloch | CEDRE software[END_REF][START_REF] Murrone | Numerical modeling of dispersed twophase flows[END_REF]. This tool has been used to simulate a number of problems involving liquid and solid propulsion systems. The solver performs calculations on unstructured meshes with cells of arbitrary geometry.

Ambipolar diffusion approximation.

The Debye length in the plume is sufficiently small to model the plasma as an ambipolar gaz mixture. Because the flow is highly turbulent, it is not necessary to account for ambipolar diffusion, hence all aspects of ionized gas physics are removed from the model, except ionization reaction which are dealt with as any other chemical reactions.

Ionized chemistry. We have used a set of chemical reactions describing the thermal ionization of alkali metals. These reactions have been obtained from previous investigations of ionization of alkali metals in rocket plume [START_REF] Smoot | Prediction of microwave attenuation characteristics of Rocket Exhaust[END_REF][START_REF] Jensen | Reaction rate coefficients for flame calculations[END_REF]. Coherence of chemical constants and thermodynamic parameters has been thoroughly checked. The ionization mechanism involves the following reactions for the chemistry of potassium : K + + e -+ M ↔ K + M [START_REF] Smoot | Prediction of microwave attenuation characteristics of Rocket Exhaust[END_REF], K + + Cl -↔ K + Cl [START_REF] Refloch | CEDRE software[END_REF], K + + Cl -+ M ↔ KCl + M [START_REF] Smoot | Prediction of microwave attenuation characteristics of Rocket Exhaust[END_REF], Cl + e -+ M ↔ Cl -+ M [START_REF] Smoot | Prediction of microwave attenuation characteristics of Rocket Exhaust[END_REF], HCl + e -↔ Cl -+ H [START_REF] Smoot | Prediction of microwave attenuation characteristics of Rocket Exhaust[END_REF], K + HCl ↔ KCl + H [START_REF] Refloch | CEDRE software[END_REF] These reactions supplement another set of chemical reactions used to simulate afterburning. The latter mechanism contains 17 reactions involving 12 species (N 2 , H 2 , O 2, O, H 2 O, OH, H, CO, CO 2 , HCl, Cl 2 , Cl). We use the Arrhenius law to describe the temperature dependence of chemical reaction rates.

Simulation of turbulence. The flow in the plume is highly turbulent. The purpose of this study is to provide the structure and the dynamics of the large structures of the flow. Of particular interest is the dynamics of ionization, closely related to the flame dynamics. Hence we perform a three-dimensional, unsteady simulation of the rocket exhaust plume. The physics include mixing of exhaust gases with ambient air, postcombustion and ionization of exhaust gases. A MiLES simulation method is exploited. As a first step, a simulation was conducted on a 1 million cells structured mesh.

B. Simulation of time domain electromagnetic response

The wave propagates in a heterogeneous medium whose properties (e.g. temperature, electronic density) fluctuate in time.

We assume that for a given wavelength, the time t (=c/L, c: speed of light and L: plume length) for the wave to travel across the medium is much lower than the characteristic time T of fluctuations of space properties of the medium. In our case the maximum fluctuation frequency f max is around 300 MHz and the time t is then much lower than T(=1/f max ). Moreover the radar wavelength λ i is much larger than the Debye length l D and than the distance between species. With this assumption, we can use a quasi-static approximation. The macroscopic Maxwell equations in harmonic time domain can be used and the influence of time can be taken into account by this decomposition:

frequency netic electromag time ) 2 exp( )) ( , , ( ) , ( =
= = ∑ e n n e n r e n e f t t f j t f t R f t S π ε (1)
The sample R(t n ) is evaluated by the harmonic solver Maxwell3D with the various dielectric permittivity ) ( n r t ε evaluated starting from the results of Cedre at each step in time Diffusion will be primarily controlled by the behaviour of the free electrons with number density N e and temperature T e .

The relative dielectric permittivity ε r is obtained from the following formula

( ) ( ) 2
j f f j f N q f m f N ω ω ω ε ω ω ω ω ω ω ε - = - = - ⎛ ⎞ + + ⎜ ⎟ ⎝ ⎠ = = = ≈ π = (2)
with ω the radar pulsation, m e electron mass, ε 0 the dielectric permittivity in free space.

The collision frequency f c reads [START_REF] Smoot | Prediction of microwave attenuation characteristics of Rocket Exhaust[END_REF] ( )

1 2 8 , B e c s s e k T e s f N Q r t m π ⎛ ⎞ = ⎜ ⎟ ⎝ ⎠ ∑ (3)
with k B the Boltzmann constant, N s the number density for species colliding with electrons and Q es electron collision cross section with the species. We formulate Maxwell's equations for 3D geometries [3][4]. Simulation cases consist of perfectly conducting bodies and of dielectric isotropic materials with or without losses (magnetic and/or electric). The problem is solved using an integral equation formulation on the boundary coupled to a differential equation formulation in the bulk. Different solvers and preconditioner can be used. In that application a MGCR_MPI (Message Passing Interface) solver was used with FMM(Fast Multipole Method) matrixvector product accelerator.

The frequency spectrum is obtained by Fourier transform of the time domain response S .

III. VALIDATION : A PROPELLANT MOCK UP MOSER

The combined solvers have been tested on a solid propellant MoSER and the results have been compared to measurement.

The simulation has been realized on an under dense meshing that does not take into account small turbulence.

A. Description of the propellant

Eight reduced scale mockup of solid propellant engine were built and fired in Onera pyrotechnic facilities.

The shot that has been the subject of the present work is a 68 % ammonium perchlorate -18 % aluminum propellant, To emphasize the physics of ionized species, 300ppm of alkali impurities were added to the composition as alkali salts (NaCl, KCl).

B. Measurement facilities

The engines have been ground fired in Onera's pyrotechnics areas.

C. Numerical results

The combined solvers have been used to simulate the far field time domain response of the mock-up at 720 MHz. The aspect angles correspond to those used in the measurement campaign. Time observation was 6.10 -3 s with time step of 2.10 - 5 s.

1) Calculation of electronic densities in the rocket plume.

Electronic and ion densities in the combustion chamber. The mass fractions of the combustion products, including ionized species, are computed as the equilibrium state obtained after the propellant has burned.

Electronic densities in the plume. We force the charged species mass fraction at nozzle entrance to be equal to those in the combustion chamber. A flux of charged species enters the computational domain trough the nozzle. An instantaneous map of the electrons mass fractions calculated by CEDRE is shown on figure 1. In the axial direction we observe that electronic densities decrease trough the nozzle but increase in the afterburning flame region. In this region the rate of thermal ionization reactions is increased by the high temperatures created by chemical reactions driving the postcombustion process. Number densities of charged species highly correlate with the temperature field. It must be noted that the ionized species create only a small perturbation with respect to the unperturbed combustion problem. In the radial direction we observe that ionization takes place in the mixing layer between the plume and fresh air flowing around the rocket. The spectrum of the electromagnetic response versus angle with 25 dB of dynamic is shown below. The frequency domain is quite large and is generated by the turbulence with more energy in the negative part due to the velocity of the particles The shapes are closed with same frequency domain. High frequencies do not appear in simulation and may be due to the under dense meshing.

D. Future work

We envision the following future work to improve our models:

• Simulation of MoSER mock up using refined meshing to better represent the turbulent eddies, • Simulation and comparison to all measurements of the MoSER campaign, • Sensitivity study to obtain a simplified model.

E. Conclusion

Our multiphysics model developed to simulate the radar response of a rocket plume has been improved to simulate the influence of the fluctuations on the electromagnetic signature dynamics. This tool combining two solvers, CEDRE and Maxwell3D, has been used to calculate the radar signature of the MoSER mock-up. The first results are quite close to the measurement in terms of signal shape. Consolidation of these results is in progress.

Number equations consecutively. Equation numbers, within parentheses, are to position flush right, as in (1), using a right tab stop. To make your equations more compact, you may use the solidus ( / ), the exp function, or appropriate exponents. Italicize Roman symbols for quantities and variables, but not Greek symbols. Use a long dash rather than a hyphen for a minus sign. Punctuate equations with commas or periods when they are part of a sentence, as in

a + b = γ (1) α + β = χ. (1) (1)
Note that the equation is centered using a center tab stop. Be sure that the symbols in your equation have been defined before or immediately following the equation. Use "(1)," not "Eq. ( 1)" or "equation [START_REF] Fromentin-Denoziere | Numerical modelling of the radar signature of rocket exhaust plumes[END_REF]," except at the beginning of a sentence: "Equation ( 1) is ..."

F. Some Common Mistakes

• The word "data" is plural, not singular.

• The subscript for the permeability of vacuum μ 0 , and other common scientific constants, is zero with subscript formatting, not a lowercase letter "o." • In American English, commas, semi-/colons, periods, question and exclamation marks are located within quotation marks only when a complete thought or name is cited, such as a title or full quotation. When quotation marks are used, instead of a bold or italic typeface, to highlight a word or phrase, punctuation should appear outside of the quotation marks. A parenthetical phrase or statement at the end of a sentence is punctuated outside of the closing parenthesis (like this). (A parenthetical sentence is punctuated within the parentheses.) • A graph within a graph is an "inset," not an "insert."

The word alternatively is preferred to the word "alternately" (unless you really mean something that alternates). • Do not use the word "essentially" to mean "approximately" or "effectively." • In your paper title, if the words "that uses" can accurately replace the word using, capitalize the "u"; if not, keep using lower-cased. • Be aware of the different meanings of the homophones "affect" and "effect," "complement" and "compliment," "discreet" and "discrete," "principal" and "principle." • Do not confuse "imply" and "infer."

• The prefix "non" is not a word; it should be joined to the word it modifies, usually without a hyphen. • There is no period after the "et" in the Latin abbreviation "et al." • The abbreviation "i.e." means "that is," and the abbreviation "e.g." means "for example."

An excellent style manual for science writers is [START_REF] Murrone | Numerical modeling of dispersed twophase flows[END_REF].

IV. USING THE TEMPLATE

After the text edit has been completed, the paper is ready for the template. Duplicate the template file by using the Save As command, and use the naming convention prescribed by your conference for the name of your paper. In this newly created file, highlight all of the contents and import your prepared text file. You are now ready to style your paper; use the scroll down window on the left of the MS Word Formatting toolbar.

A. Authors and Affiliations

The template is designed so that author affiliations are not repeated each time for multiple authors of the same affiliation. Please keep your affiliations as succinct as possible (for example, do not differentiate among departments of the same organization). This template was designed for two affiliations. e) Reassign number of columns: Place your cursor to the right of the last character of the last affiliation line of an even numbered affiliation (e.g., if there are five affiliations, place your cursor at end of fourth affiliation). Drag the cursor up to highlight all of the above author and affiliation lines. Go to Column icon and select "2 Columns". If you have an odd number of affiliations, the final affiliation will be centered on the page; all previous will be in two columns.

B. Identify the Headings

Headings, or heads, are organizational devices that guide the reader through your paper. There are two types: component heads and text heads.

Component heads identify the different components of your paper and are not topically subordinate to each other. Examples include ACKNOWLEDGMENTS and REFERENCES, and for these, the correct style to use is "Heading 5." Use "figure caption" for your Figure captions, and "table head" for your table title. Run-in heads, such as "Abstract," will require you to apply a style (in this case, italic) in addition to the style provided by the drop down menu to differentiate the head from the text.

Text heads organize the topics on a relational, hierarchical basis. For example, the paper title is the primary text head because all subsequent material relates and elaborates on this one topic. If there are two or more sub-topics, the next level head (uppercase Roman numerals) should be used and, conversely, if there are not at least two sub-topics, then no subheads should be introduced. Styles named "Heading 1," "Heading 2," "Heading 3," and "Heading 4" are prescribed.

We suggest that you use a text box to insert a graphic (which is ideally a 300 dpi resolution TIFF or EPS file with all fonts embedded) because this method is somewhat more stable than directly inserting a picture. As an example, write the quantity "Magnetization," or "Magnetization, M," not just "M." If including units in the label, present them within parentheses. Do not label axes only with units. In the example, write "Magnetization (A/m)" or "Magnetization (A ( m(1)," not just "A/m." Do not label axes with a ratio of quantities and units. For example, write "Temperature (K)," not "Temperature/K."

ACKNOWLEDGMENT (Heading 5)

The preferred spelling of the word "acknowledgment" in America is without an "e" after the "g." Avoid the stilted expression "one of us (R. B. G.) thanks ...". Instead, try "R. B. G. thanks...". Put sponsor acknowledgments in the unnumbered footnote on the first page.

The template will number citations consecutively within brackets [START_REF] Fromentin-Denoziere | Numerical modelling of the radar signature of rocket exhaust plumes[END_REF]. The sentence punctuation follows the bracket [START_REF] Gueyffier | Numerical Simulation of Ionized Rocket Plumes[END_REF]. Refer simply to the reference number, as in [START_REF] Levadoux | Recent advances in the pre-conditioning of integral equations of electromagnetism[END_REF]-do not use "Ref. [START_REF] Levadoux | Recent advances in the pre-conditioning of integral equations of electromagnetism[END_REF]" or "reference [START_REF] Levadoux | Recent advances in the pre-conditioning of integral equations of electromagnetism[END_REF]" except at the beginning of a sentence: "Reference [START_REF] Levadoux | Recent advances in the pre-conditioning of integral equations of electromagnetism[END_REF] was the first ..." Number footnotes separately in superscripts. Place the actual footnote at the bottom of the column in which it was cited. Do not put footnotes in the reference list. Use letters for table footnotes.

Unless there are six authors or more give all authors' names; do not use "et al.". Papers that have not been published, even if they have been submitted for publication, should be cited as "unpublished" [START_REF] Levadoux | New trends in the preconditioning of integral equations of electromagnetism[END_REF]. Papers that have been accepted for publication should be cited as "in press" [START_REF] Smoot | Prediction of microwave attenuation characteristics of Rocket Exhaust[END_REF]. Capitalize only the first word in a paper title, except for proper nouns and element symbols.

For papers published in translation journals, please give the English citation first, followed by the original foreign-language citation [START_REF] Refloch | CEDRE software[END_REF].

  Facilities have been used to measure: • Pressure • Number density Ne by electromagnetic transmission at 10GHz • Electromagnetic forward and backward signature in frequency domain [830MHz-880MHz] • Frequency spectrum of the radar signature between 700-760MHz Eight reduced scale mockup of solid propellant engine were built and fired in Onera pyrotechnic facilities.

Fig 1 :

 1 Fig 1: Electron mass fraction in the plume 2) Electromagnetic simulation in UHF band.

Fig 2 :

 2 Fig 2: simulation in UHF band

Fig 3 :

 3 Fig 3: dielectric permittivity at a time step

Fig 4 :

 4 Fig 4: Frequency spectrum of the radar energy for different angle aspect 1) Comparison of the results. The histograms of the spectra with a dynamic of 25dB are presented below:

Fig 5 :

 5 Fig 5: frequency spectra of energy in UHF band

1 )

 1 For author/s of only one affiliation (Heading 3): To change the default, adjust the template as follows. a) Selection (Heading 4): Highlight all author and affiliation lines. b) Change number of columns: Select the Columns icon from the MS Word Standard toolbar and then select "1 Column" from the selection palette. c) Deletion: Delete the author and affiliation lines for the second affiliation. 2) For author/s of more than two affiliations: To change the default, adjust the template as follows. a) Selection: Highlight all author and affiliation lines. b) Change number of columns: Select the "Columns" icon from the MS Word Standard toolbar and then select "1 Column" from the selection palette. c) Highlight author and affiliation lines of affiliation 1 and copy this selection. d) Formatting: Insert one hard return immediately after the last character of the last affiliation line. Then paste down the copy of affiliation 1. Repeat as necessary for each additional affiliation.

Fig. 1 .

 1 Fig. 1. Example of a figure caption. (figure caption) Figure Labels: Use 8 point Times New Roman for Figure labels. Use words rather than symbols or abbreviations when writing Figure axis labels to avoid confusing the reader.As an example, write the quantity "Magnetization," or "Magnetization, M," not just "M." If including units in the label, present them within parentheses. Do not label axes only with units. In the example, write "Magnetization (A/m)" or "Magnetization (A ( m(1)," not just "A/m." Do not label axes

  To have non-visible rules on your frame, use the MSWord "Format" pull-down menu, select Text Box > Colors and Lines to choose No Fill and No Line.C. Figures and Tables 1) Positioning Figures and Tables: Place figures and tables at the top and bottom of columns. Avoid placing them in the middle of columns. Large figures and tables may span across both columns. Figure captions should be below the figures; table heads should appear above the tables. Insertfigures and tables after they are cited in the text. Use the abbreviation "Fig.1," even at the beginning of a sentence.

TABLE I .TABLE STYLES Table Column Head Table Head Table column subhead Subhead Subhead copy

 ISTYLES

More table copy a a. Sample of a Table footnote. (Table footnote) b.


Daniel Hirschkoff

Damien Pous

A Distribution Law for CCS and a New Congruence Result for the π-calculus

We give an axiomatisation of strong bisimilarity on a small fragment of CCS that does not feature the sum operator. This axiomatisation is then used to derive congruence of strong bisimilarity in the finite π-calculus in absence of sum. To our knowledge, this is the only nontrivial subcalculus of the π-calculus that includes the full output prefix and for which strong bisimilarity is a congruence.

Introduction

In this paper, we study strong bisimilarity on two process calculi. We first focus on microCCS (µCCS), the very restricted fragment of CCS that only features prefix and parallel composition. Our main result on µCCS is that adding the following distribution law to the laws of an abelian monoid for parallel composition yields a complete axiomatisation of strong bisimilarity (in the law above, η is a CCS prefix, of the form a or a, and P is any CCS process -the same number of copies of P appear on both sides of the equation).

The distribution law is not new: it is mentioned -among other 'mixed equations' relating prefixed terms and parallel compositions -in a study of bisimilarity on normed PA processes [START_REF] Hirshfeld | Bisimulation Equivalence is Decidable for Normed Process Algebra[END_REF]. In our setting, this equality can be oriented from left to right to rewrite processes into normal forms, which intuitively exhibit as much concurrency as possible. Strong bisimilarity (∼) between processes is then equivalent to equality of their normal forms. This rewriting phase allows us to actually compute unique decompositions of processes into prime processes, in the sense of [START_REF] Milner | Unique Decomposition of Processes[END_REF]: a process P is prime if P is not bisimilar to the inactive process 0 and if

P ∼ Q | R implies Q ∼ 0 or R ∼ 0.
The distribution law is an equational schema, corresponding to an infinite family of axioms, of the form η.(P | (η.P ) k ) = (η.P ) k+1 , for k ≥ 1 (where Q k denotes the k-fold parallel composition of process Q). We show that although our setting is rather simple, there exists no finite axiomatisation of ∼ on µCCS.

We then move to the study of strong bisimilarity in the π-calculus. Because of the presence of the input prefix, and of the related phenomenon of name-passing, bisimilarity is more complex in the π-calculus than in CCS. In particular, both early and late bisimilarity, that differ in their treatment of name substitution, fail to be congruences in the full π-calculus.

There exist subcalculi of the π-calculus for which strong bisimilarity is a congruence (we discuss these in Section 5). When this is the case, this equivalence coincides with ground bisimilarity (∼ g ), which allows one to consider only one fresh name when inspecting an input transition, instead of the usual quantification involving all free names of the process. Congruence of strong bisimilarity is hence an important property: not only is it necessary in order to reason in a compositional way, but it also helps making bisimulation proofs simpler, by reducing the size of case analyses.

In the full π-calculus, in order to get congruence, one has to work with Sangiorgi's open bisimilarity [START_REF] Sangiorgi | A Theory of Bisimulation for the π-Calculus[END_REF], which has a more involved definition than the early and late variants. Tools like the Mobility Workbench [START_REF] Victor | The Mobility Workbench[END_REF], for instance, have adopted this equivalence on processes.

It is known [START_REF] Sangiorgi | The π-calculus: a Theory of Mobile Processes[END_REF] that bisimilarity in the π-calculus fails to be a congruence as soon as we have prefix, parallel composition, restriction and replication. In this work, we focus on the finite, sum-free π-calculus, that we call π 0 . We rely on the axiomatisation of strong bisimilarity on µCCS to prove that ground bisimilarity (∼ g ) is closed under substitutions in π 0 , i.e., that whenever P ∼ g Q, then P σ ∼ g Qσ for any substitution σ. Closure under substitution of ground bisimilarity entails that on π 0 , ground, early, late and open bisimilarities coincide, and are congruences. The problem of congruence of ∼ g on π 0 is mentioned as an open question in [START_REF] Sangiorgi | The π-calculus: a Theory of Mobile Processes[END_REF]Chapter 5], and is known since at least 1998 [START_REF] Boreale | Some Congruence Properties for π-calculus Bisimilarities[END_REF]. To our knowledge, this is the first congruence result for a subcalculus of the π-calculus that includes the full output prefix (see Section 5 for a discussion on this).

At the heart of our proof of congruence is a notion that we call mutual desynchronisation, and that corresponds to the existence of processes T, T 12 , T 21 such that T η1 -→ η2 -→ T 12 and T η2 -→ η1 -→ T 21 , for two distinct actions η 1 and η 2 , and with T 12 ∼ T 21 . We additionally require in the two sequences of transitions from T to T 12 and T 21 respectively that the second prefix being fired should occur under the first prefix in T . In other words, in such a situation, the process T behaves as if the two actions η 1 , η 2 were offered concurrently, but the simultaneous firing of these actions can only be emulated by triggering consecutive prefixes.

Using our analysis of strong bisimilarity on µCCS, we show that mutual desynchronisations do not exist in µCCS. This is essentially due to the fact that our axiomatisation of ∼ on µCCS does not allow one to relate two distinct prefixes when performed concurrently and sequentially. When moving to the π-calculus, it turns out that substitution closure of ∼ g amounts to observing absence of mutual desynchronisations in π 0 . We exploit a transfer property, that extracts a bisimilarity proof in µCCS from a bisimilarity proof in π 0 , to relate the two calculi and to show that mutual desynchronisations do not exist in π 0 , yielding congruence of ∼ g .

Paper outline. We introduce µCCS and the distribution law in Section 1. Section 2 is devoted to the characterisation of ∼ on µCCS using normal forms. In Section 3, we prove that no finite axiomatisation of ∼ on µCCS exists. Section 4 presents the proof of our congruence result in the π-calculus, and we give concluding remarks in Section 5.

MicroCCS Processes and Normal Forms

We consider an infinite set N of names, ranged over with a, b . . . . We define on top of N the set of processes of µCCS, the finite, public, sum-free CCS calculus, ranged over using P, Q, R . . . , as follows:

η ::= a a , P ::= 0 η.P P 1 | P 2 .

0 is the nil process. η ranges over visible actions and co-actions, called interactions, and we let η stand for the co-action associated to η (we have η = η).

For k > 0, we write P k for the parallel composition of k copies of P , and we write i∈I P i for the parallel composition of all processes P i for i ∈ I. It can be noted that our syntax does not include a construction of the form τ.P -see Remark 2.3 below.

Structural congruence, written ≡, is defined as the smallest congruence satisfying the following laws:

(C 1 ) P | Q ≡ Q | P (C 2 ) P | (Q | R) ≡ (P | Q) | R (C 3 ) P | 0 ≡ P
We introduce a labelled transition system (LTS) for µCCS. Actions labelling transitions, ranged over with µ, are either interactions, or a special silent action, written τ .

Definition 1.1 (Operational semantics and behavioural equivalence).

The LTS for µCCS is given by the following rules:

η.P η -→ P P η -→ P Q η -→ Q P | Q τ -→ P | Q P µ -→ P P | Q µ -→ P | Q P µ -→ P Q | P µ -→ Q | P
A bisimulation is a symmetrical relation R between processes such that whenever P R Q and P µ -→ P , there exists Q such that Q µ -→ Q and P R Q . Bisimilarity, written ∼, is the union of all bisimulations. Definition 1.2 (Size). Given P , #(P ) (called the size of P ) is defined by:

#(0) def = 0 #(P 1 | P 2 ) def = #(P 1 ) + #(P 2 ) #(η.P ) def = 1 + #(P ) . Lemma 1.3. P ≡ Q implies P ∼ Q which in turn implies #(P ) = #(Q).
Proof. The first implication follows by proving that ≡ is a bisimulation. Suppose then by contradiction that there exist P, Q such that P ∼ Q and #(P ) < #(Q); and choose such P with minimal size. Q has at least one prefix: Q η -→ Q and we get P η -→ P with P ∼ Q . Necessarily, we must have #(P ) < #(Q ) and #(P ) < #(P ) which contradicts the minimality hypothesis. We shall use this equality, oriented from left to right, to rewrite processes. We write P P when there exist P 1 , P 2 such that P ≡ P 1 , P 2 ≡ P and P 2 is obtained from P 1 by replacing a sub-term of the form of the left-hand side process with the right-hand side process. Remark 1.1 (On the distribution law and PA). Among the studies about properties of ∼ in process algebras that include parallel composition (see [START_REF] Aceto | Finite Equational Bases in Process Algebra: Results and Open Questions[END_REF] for a recent survey on axiomatisations), some works focus on calculi where parallel composition is treated as a primitive operator (as opposed to being expressible using sum or other constructs like the left merge operator). As mentioned above, particularly relevant to this work is [START_REF] Hirshfeld | Bisimulation Equivalence is Decidable for Normed Process Algebra[END_REF], where Hirshfeld and Jerrum "develop a structure theory for PA that completely classifies the situations in which a sequential composition of two processes can be bisimilar to a parallel composition". [START_REF] Hirshfeld | Bisimulation Equivalence is Decidable for Normed Process Algebra[END_REF] establishes decidability of ∼ for normed PA processes: in that setting, the formal analogue of the distribution law (Def. 1.4) holds with η and P being two processes -the 'dot' operator is a general form of sequential composition. This equality is valid in [START_REF] Hirshfeld | Bisimulation Equivalence is Decidable for Normed Process Algebra[END_REF] whenever η is a 'monomorphic process', meaning that η can only reduce to 0 (which corresponds to µCCS), or to η itself. [START_REF] Fokkink | An ω-complete Equational Specification of Interleaving[END_REF] presents a finite axiomatisation of PA that exploits the operators of sum and left merge.

Lemma 1.5. The relation is strongly normalising and confluent.

Proof. If P P then the weight of P (defined as sum of the depths of all prefixes occuring in P ) is strictly smaller than the weight of P , whence the strong normalisation. We then remark that is locally confluent, and conclude with Newman's Lemma. Thus, for any process P , defines a normal form unique up to ≡, that will be denoted by n(P ). We let A, B, . . . range over normal forms.

The following lemma states that preserves bisimilarity:

Lemma 1.6. If P P , then P ∼ P . For any P , P ∼ n(P ).

Proof. The relation ( ∪ ( ) -1 ∪ ≡) is a bisimulation.

Characterisation of Bisimilarity in MicroCCS

Our characterisation of ∼ on µCCS makes use of the notion of decomposition into prime processes, defined as follows:

Definition 2.

i ∼ Q i for all i ∈ [1.
.n], up to a permutation of the indices.

Proof. Similar to the proof of [START_REF] Moller | Axioms for Concurrency[END_REF]Theorem 4.3.1]: the case of µCCS is not explicitly treated in that work, but the proof can be adapted rather easily.

An immediate consequence of the above result is the following property:

Corollary 2.3 (Cancellation). For all P, Q, R, P | R ∼ Q | R implies P ∼ Q. Note that this is not true in presence of replication: a | !a ∼ 0 | !a, but a ∼ 0.
The characterisation of ∼ using the distribution law follows from the observation that if a normal form is a prefixed process, then it is prime. This idea is used in the proof of Lemma 2.5. We first establish a technical result, that essentially exploits the same argument as the proof of Theorem 4.2 in [START_REF] Hirshfeld | Bisimulation Equivalence is Decidable for Normed Process Algebra[END_REF].

Lemma 2.4. If η.P ∼ Q | Q , with Q, Q ∼ 0, then there exist A and k > 1 such that η.P ∼ (η.A) k and η.A is a normal form. Proof. By Lemma 1.6, we have η.P ∼ n(Q | Q ). Furthermore, we have that n(Q | Q ) ≡ i≤k η i .A i ,
where k > 1 and the processes η i .A i are in normal form.

Since the η prefix must be triggered to answer any challenge from the right hand side, we have η i = η and P ∼ A i | l =i η.A l for all i ≤ k. In particular, when i = j, we have

P ∼ A i | η.A j | l ∈{i,j} η.A l ∼ η.A i | A j | l ∈{i,j} η.A l and hence, by Corollary 2.3, A i | η.A j ∼ η.A i | A j .
By reasoning on the sizes of the parallel components in the prime decompositions of these two terms, we conclude that η.A i ∼ η.A j for all i, j ≤ k.

Hence, we have η.P ∼ (η.A 1 ) k with k > 1 and η.A 1 is a normal form. Proof. We show by induction on n that for all A with #(A) = n, we have

(i) if A is a prefixed process, then A is prime; (ii) for any B, A ∼ B implies A ≡ B.
The case n = 0 is immediate. Suppose that the property holds for all i < n, with n ≥ 1.

(i) We write A = η.A , and suppose by contradiction A ∼ P Theorem 2.6. Let P, Q be two µCCS processes. Then P ∼ Q iff n(P ) ≡ n(Q).

A ≡ B. • Otherwise, A = η 1 .A 1 | . . . | η k .A k with k > 1,
Remark 2.1 (Unique decomposition of processes). Our proof relies on unique decomposition of processes (Prop. 2.2), that first appeared in [START_REF] Milner | Unique Decomposition of Processes[END_REF]. Unique decomposition has been established for a variety of process algebras, and used as a way to prove decidability of behavioural equivalence and to give complexity bounds for the associated decision procedure ([9, 3] cite relevant references).

In the present study, beyond the existence of a unique decomposition, we are interested in a syntactic characterisation of ∼ (which will in particular allow us to derive Lemma 4.6 below). In this sense, our work is close to [START_REF] Corradini | Towards parallelization of concurrent systems[END_REF], where the notion of maximally parallel process in CCS (with choice) is studied. [START_REF] Corradini | Towards parallelization of concurrent systems[END_REF] defines a rewriting process through which maximally parallel normal forms can be computed, and shows that in the case of µCCS, such normal forms are unique. However, no syntactical characterisation of the set of normal forms is presented, and such a characterisation cannot be directly deduced from the (rather involved) definition of the rewriting process for full CCS.

We instead restrict ourselves to µCCS from the start, and rely explicitly on the distribution law in order to 'extract' prime components of processes.

Remark 2.2 (Closure under substitutions). In (full) CCS, two strongly bisimilar processes need not remain bisimilar whenever we apply a substitution that replaces names with names. The standard counterexample is given by a.b + b.a ∼ a | b: when we replace b with a, we obtain two processes that are distinguished by ∼, since the latter can perform a τ transition that cannot be matched by the former. This irregularity is the basis of the standard counterexample showing that strong bisimilarity is not a congruence in the π-calculus.

In µCCS, on the other hand, ∼ is closed under substitutions: the intuitive reason is that two processes related by an instance of the distribution law remain equivalent when a substitution is applied (we can show in particular that for any substitution σ, n(P σ) ≡ n(n(P )σ)). This is not the case for the expansion law, of which the counterexample above is an instance.

Remark 2.3 (τ prefix and weak bisimilarity). We do not address weak bisimilarity in the present work. In µCCS, strong and weak bisimilarity coincide, i.e., the internal transitions of processes are completely determined by the visible actions (interactions). When including τ prefixes in the syntax, it can be proved that adding the law τ.P = P is enough to characterise weak bisimilarity. The τ prefix is usually absent in the π-calculus, to which we shall move in Section 4. Since some results on CCS will be transferred to the π-calculus, we did not include this construct in µCCS.

Nonexistence of a Finite Axiomatisation

We let D stand for the set of equations consisting of the three axioms of structural congruence (C 1 , C 2 , C 3 ), and the infinite family of distribution axioms

(D k ) : η.(P | (η.P ) k ) = (η.P ) k+1 , k ≥ 1 .
We let D k stand for the finite restriction of D where only the first k distribution axioms are included ((D i ) 1≤i≤k ). We shall write E P = Q whenever P = Q can be derived in equational logic using a given set E of axioms, and E P = Q when this is not the case.

(D k ) k≥1 forms an equational schema for the distribution law, and Theorem 2.6 states that D is a complete axiomatisation of strong bisimilarity on µCCS. Using a rather classical approach (i.e., establishing ω-completeness and proving compactness, see [START_REF] Aceto | Finite Equational Bases in Process Algebra: Results and Open Questions[END_REF]), this leads to the nonexistence of a finite axiomatisation of ∼ on µCCS. The lemma below provides the central technical property satisfied by the (D k ) k≥1 which is necessary to derive Theorem 3.2, that says that D is intrinsically infinite. Lemma 3.1. Let a be a name. For any k, there exists n s.t. D k a.a n = a n+1 .

Remember that a n stands for the n-fold parallel composition of a.0, so that the above equality is an instance of axiom (D n ).

Proof. Let n be a number strictly greater than k such that n + 1 is prime, and let θ(P, Q) denote the predicate: "P ∼ Q ∼ a n+1 , P ≡ a.P , and

Q ≡ Q 1 | Q 2 with Q 1 , Q 2 ≡ 0". Suppose now that D k
a.a n = a n+1 , and consider the shortest proof of D k P = Q for some processes P, Q such that either θ(P, Q) or θ(Q, P ). Since θ(a.a n , a n+1 ) holds, such a minimal proof does exist. We reason about the last rule used in the derivation of this proof in equational logic. For syntactic reasons, this cannot be reflexivity, a contextual rule, nor one of the structural congruence axioms. It can be neither symmetry nor transitivity, since otherwise this would give a shorter proof satisfying θ. The only possibility is thus the use of one of the distribution axioms, say D i with 1 ≤ i ≤ k and a n+1 ∼ Q ≡ (a.Q ) i+1 . By Lemma 1.3, since #(a n+1 ) = n + 1, i + 1 has to divide n + 1. This is contradictory, because we have 2 ≤ i + 1 ≤ k + 1 < n + 1, and n + 1 is prime. Theorem 3.2 (No finite axiomatisation of ∼). For any finite set of axioms E, there exist processes P and Q such that P ∼ Q but E P = Q.

Proof. Standard, by proving that D is ω-complete and then using the Compactness Theorem (see [START_REF] Aceto | Finite Equational Bases in Process Algebra: Results and Open Questions[END_REF]).

A New Congruence Result for the π-calculus

4.1 The Finite, Sum-free π-calculus π-calculus processes are built from an infinite set N π of names, ranged over using a, b . . . , m, n . . . , p, q . . . , x, y . . . , according to the following grammar: φ ::= m(x) mn , P ::= 0 φ.P P 1 | P 2 (νp)P .

The input prefix m(x) binds name x in the continuation process, and so does name restriction (νn) in the restricted process. A name that is not bound is said to be free, and we let fn(P ) stand for the free names of P . We assume that any process that we manipulate satisfies a Barendregt convention: every bound name is distinct from the other bound and free names of the process. We shall use a, b, c to range over free names of processes, p, q, r (resp. x, y) to range over names bound by restriction (resp. by input), and m, n to range over any name, free or bound (note that these naming conventions are used in the above grammar). Structural congruence on π 0 , written ≡, is the smallest congruence that is an equivalence relation, contains α-equivalence, and satisfies the following laws:

P | 0 ≡ P P | (Q | R) ≡ (P | Q) | R P | Q ≡ Q | P (νp)0 ≡ 0 (νp)(νq)P ≡ (νq)(νp)P P | (νp)Q ≡ (νp)(P | Q) if p / ∈ fn(P )
We let P [n/x] stand for the capture avoiding substitution of name x with name n in P . We use σ to range over substitutions in π 0 (that simultaneously replace several names).

Definition 4.1 (Late operational semantics and ground bisimilarity).

The late operational semantics of π 0 is given by a transition relation whose set of labels is defined by: µ ::= a(x) ab a(p) τ .

Names x and p are said to be bound in actions a(x) and a(p) respectively, and other names are free. We use bn(µ) (resp. fn(µ)) to denote the set of bound (resp. free) names of action µ.

The late transition relation, written -→ π , is given by the following rules (symmetrical versions of the rules involving parallel composition are omitted):

φ.P φ -→ π P P a(x) ---→ π P Q ab -→ π Q P | Q τ -→ π P [b/x] | Q P ab -→ π P (νb)P a(b) --→ π P a = b P a(x) ---→ π P Q a(p) --→ π Q P | Q τ -→ π (νp)(P [p/x] | Q ) P µ -→ π P P | Q µ -→ π P | Q bn(µ) ∩ fn(Q)=∅ P µ -→ π P (νp)P µ -→ π (νp)P p / ∈ fn(µ)
A ground bisimulation is a symmetric relation R between processes such that whenever P R Q and P

µ -→ π P , there exists Q s.t. Q µ -→ π Q and P R Q .
Ground bisimilarity, written ∼ g , is the union of all ground bisimulations. Note that we do not respect the convention on names in the rule to infer a bound output, precisely because we are transforming a free name (b) into a bound name. Proof. Similar to the proof of Lemma 1.4.13 in [START_REF] Sangiorgi | The π-calculus: a Theory of Mobile Processes[END_REF], where the early transition semantics is treated.

Mutual Desynchronisations

We now introduce the notion of mutual desynchronisation in µCCS, which is defined as the existence of processes obeying certain conditions in the calculus. We shall see that because of τ synchronisations, the absence of mutual desynchronisations is related to substitution closure of ∼.

Definition 4.3 (Mutual desynchronisation in µCCS)

. We say that there exists a mutual desynchronisation in µCCS whenever there are two prefixes η 1 , η 2 , and five µCCS processes P, P , Q, Q , R such that

η 1 = η 2 , P η1 -→ P , Q η2 -→ Q and η 2 .P | Q | R ∼ P | η 1 .Q | R.
The notion of mutual desynchronisation is not specific to µCCS. As explained in the introduction, it corresponds to a situation where three processes T, T 12 , T 21 satisfy:

-T η1 -→ η2 -→ T 12 and T η2 -→ η1 -→ T 21
, where the second prefix being triggered occurs under the first one in both sequences of transitions.

η 1 = η 2 and T 12 ∼ T 21 .

The proofs of Lemmas 4.9 and 4.10 will expose analogous situations in π 0 . Definition 4.4. We define, for any µCCS process P and prefix η, the contribution of P at η, written s η (P ), by

s η (0) def = 0 s η (η .P ) def = 0 if η = η s η (P 1 | P 2 ) def = s η (P 1 ) + s η (P 2 ) s η (η.P ) def = #(η.P )
Intuitively, s η (P ) is the total size of the parallel components of P that start with the prefix η. Proof. Follows from Theorem 2.6 and the observation that the distribution law preserves the contribution of a process at a given interaction prefix.

Lemma 4.6 (No mutual desynchronisation).

There exists no mutual desynchronisation in µCCS.

Proof. Suppose by contradiction that there are processes such that P η1

-→ P , Q η2 -→ Q and η 2 .P | Q | R ∼ P | η 1 .Q | R.
By the cancellation property (Corollary 2.3), we have

η 2 .P | Q ∼ P | η 1 .Q, hence for all η, s η (η 2 .P | Q ) = s η (P | η 1 .Q) (Lemma 4.5). Since s η1 (η 2 .P | Q ) = s η1 (Q ) ≤ #(Q ) and s η1 (P | η 1 .Q)) ≥ s η1 (η 1 .Q) = #(Q ) + 2, by taking η = η 1 we finally get #(Q ) ≥ #(Q ) + 2.
This result will be used to show that a situation corresponding to a mutual desynchronisation cannot arise in π 0 . Notice that the proof depends in an essential way on Lemma 4.5, which in turn relies on the axiomatisation of ∼ in µCCS (Theorem 2.

6).

In what follows, we fix two distinct names a and b, that will occur free in the processes we shall consider. The definitions and results below will depend on a and b, but we avoid making this dependency explicit, in order to ease readability. Names a and b will be fixed in the proof of Lemma 4.11. Definition 4.7 (Erasing a π 0 process). Given a π 0 process P , we define the erasing of P , written E(P ), as follows:

E(P 1 | P 2 ) def = E(P 1 ) | E(P 2 ) E((νp)P ) def = E(P ) E(0) def = 0 E(a(x).P ) def = a.E(P ) E(m(x).P ) def = 0 if m = a E(bn.P ) def = b.E(P ) E(mn.P ) def = 0 if m = b
Note that a and b play different roles in the definition of E(•).

It is immediate from the definition that E(P ) is a µCCS process whose only prefixes are a and b. Intuitively, E(P ) only exhibits the interactions of P at a (in input) and b (in output) that are not guarded by interactions on other names. Lemma 4.8 (Transitions of E(P )). Consider a π 0 process P . We have:

-If P a(x) ---→ π P , then E(P ) a -→ E(P ). -If P bc -→ π P or P b(p) --→ π P , then E(P ) b - → E(P ). -Conversely, if E(P )
If P ∼ g Q in π 0 , then E(P ) ∼ E(Q) in µCCS.
Proof. We reason by induction on the size of P (defined as the number of prefixes in P ). Consider a transition of E(P ); as observed above, it can only be a transition along a or a transition along b.

Suppose E(P ) a -→ P 0 . By Lemma 4.8, P a(x) ---→ π P and P 0 = E(P ). Since

P ∼ g Q, Q a(x) ---→ π Q for some Q such that P ∼ g Q . By induction, the latter relation gives E(P ) ∼ E(Q ), and Q a(x) ---→ π Q gives by Lemma 4.8 E(Q) a -→ E(Q ).
The case E(P ) b -→ P 0 is treated similarly: by Lemma 4.8, there are two cases, according to whether P does a free output or a bound output. Reasoning like above allows us to conclude in both cases.

We can now present our central technical result about π 0 , which comes in two lemmas. Lemma 4.9. If Q ∼ g (ν p)(a(x).P 1 | bc.P 2 | P 3 ), then there exist some

Q 1 , Q 2 , Q 3 , q, such that Q ≡ (ν q)(a(x).Q 1 | bc.Q 2 | Q 3 ) and (ν p)(P 1 | P 2 | P 3 ) ∼ g (ν q)(Q 1 | Q 2 | Q 3 ).
Proof. Let P = (ν p)(a(x).P 1 | bc.P 2 | P 3 ) and P = (ν p)(P 1 | P 2 | P 3 ).

Note that by our conventions on notations, c / ∈ p. Since Q ∼ g P and P can perform two transitions along a(x) and bc respectively, Q can also perform these transitions, which gives

Q ≡ (ν q)(a(x).Q 1 | bc.Q 2 | Q 3 ) for some q, Q 1 , Q 2 , Q 3
, the first (resp. second) component exhibiting the prefix that is triggered to answer the challenge on a(x) (resp. bc).

Consider now the challenge

P bc -→ π a(x) ---→ π P , to which Q answers by per- forming Q bc -→ π a(x) ---→ π Q ba , with P ∼ g Q ba . If Q ba = (ν q)(Q 1 | Q 2 | Q 3 ), that is, if Q triggers the
prefixes on top of its first and second components, then we are done. Similarly, if Q triggers a prefix in Q 3 to answer the second challenge, say

Q 3 = a(x).Q 4 | Q 5 , we can set Q 1 = a(x).Q 4 and Q 3 = Q 1 | Q 5 ,
and the lemma is proved.

The case that remains to be analysed is when

Q 2 a(x) ---→ π Q 2 and Q ba = (ν q)(a(x).Q 1 | Q 2 | Q 3 ) ∼ g (ν p)(P 1 | P 2 | P 3 ).
We then consider the challenge where P fires its two topmost prefixes a(x) and bc in the other sequence, namely P a(x) ---→ π bc -→ π P . By hypothesis, Q triggers the prefix of its first component for the first transition. To perform the second transition, Q can fire the prefix bc either in its second or third component, in which case, as above, we are done, or, and this is the last possibility, the prefix bc occurs in Q 1 . This means

Q ab = (ν q)(Q 1 | bc.Q 2 | Q 3 ) ∼ g (ν p)(P 1 | P 2 | P 3 ), with Q 1 bc -→ π Q 1 .
To sum up, we have

Q ab = (ν q)(Q 1 | bc.Q 2 | Q 3 ) ∼ g (ν q)(a(x).Q 1 | Q 2 | Q 3 ) = Q ba , with Q 1 bc -→ π Q 1 and Q 2 a(x)
---→ π Q 2 : this resembles the mutual desynchronisation of Definition 4.3, translated into the π-calculus.

Indeed, we can construct a mutual desynchronisation in µCCS:

Q ab ∼ g Q ba implies E(Q ab ) ∼ E(Q ba ) by Prop. 4.1, and Q 1 bc -→ π Q 1 (resp. Q 2 a(x) ---→ π Q 2 )
implies by Lemma 4.8

E(Q 1 ) b - → E(Q 1 ) (resp. E(Q 2 ) a -→ E(Q 2 )
). Finally, using Lemma 4.6, we obtain a contradiction, which concludes our proof. Lemma 4.10. If Q ∼ g (νp, p)(a(x).P 1 | bp.P 2 | P 3 ), then there exist some

Q 1 , Q 2 , Q 3 , such that Q ≡ (νp, q)(a(x).Q 1 | bp.Q 2 | Q 3 ) and (ν p)(P 1 | P 2 | P 3 ) ∼ g (ν q)(Q 1 | Q 2 | Q 3 ).
Proof (Hint). The proof follows the same lines as for the previous lemma. The only difference is when analysing the transitions that lead to Q ab : to perform the second transition, Q can either extrude the name called p in the equality Q ≡ (νp, q)(a(x).Q 1 | bp.Q 2 | Q 3 ), or otherwise Q can be α-converted in order to extrude another name. In the case where Q chooses to extrude a different name, we can suppose without loss of generality that the necessary α-conversion is a swapping between name p and a name q 1 ∈ q, which brings us back to the case where name p is the one being extruded.

The presence of a bound output introduces some notational complications when expressing Q ab , but basically it does not affect the proof w.r.t. the proof of Lemma 4.9, because the function E(•) is not sensitive to name permutations that do not involve a or b.

Congruence

Theorem 4.11 (Closure of ∼ g under substitution). If P ∼ g Q then for any substitution σ, P σ ∼ g Qσ.

Proof. We prove that the relation R def = {(P σ, Qσ) | P ∼ g Q} is a ground bisimulation. We consider P , Q such that P ∼ g Q and suppose P σ µ -→ π P 0 . We examine the transitions of P that make it possible for P σ to do a µ-transition to P 0 .

According to Lemma 4.2, there are two possibilities. The first possibility corresponds to the situation where µ comes from an action that P can perform, i.e., P µ -→ π P for some µ , with P σ = P 0 and µ σ = µ (cases 1 and 2a in Lemma 4.2). Since P

∼ g Q, Q µ -→ π Q and P ∼ g Q for some Q . We can prove that Qσ µ -→ Q σ, and since P ∼ g Q we have (P σ, Q σ) ∈ R.
The second possibility (which corresponds to the difficult case) is given by µ = τ , where the synchronisation in P has been made possible by the application of σ. There are in turn two cases, corresponding to whether the synchronisation involves a free or a bound name. In the former case, P a(x) ---→ π P and P bc -→ π P for some a, x, b, c, P , P . This entails P ≡ (ν p)(a(x).P 1 | bc.P 2 | P 3 ) for some p, P 1 , P 2 , P 3 , and, since P ∼ g Q, we conclude by Lemma 4.9 that

Q ≡ (ν q)(a(x).Q 1 | bc.Q 2 | Q 3 ) and (ν p)(P 1 | P 2 | P 3 ) ∼ g (ν q)(Q 1 | Q 2 | Q 3 ) .
By definition of R, this equivalence implies that we can apply any substitution to these two processes to yield processes related by R, and in particular [c/x]σ, which gives:

((ν p)(P 1 | P 2 | P 3 ))[c/x]σ R ((ν q)(Q 1 | Q 2 | Q 3 ))[c/x]σ .
Using the Barendregt convention hypothesis, this amounts to

P 0 ≡ ((ν p)(P 1 [c/x] | P 2 | P 3 ))σ R ((ν q)(Q 1 [c/x] | Q 2 | Q 3 ))σ def = Q 0 .
We can then conclude by checking that Qσ τ -→ π Q 0 .

We reason similarly for the case where the synchronisation involves the transmission of a bound name, using Lemma 4.10 instead of Lemma 4.9. We remark that Lemma 4.10 gives (ν p)(P

1 | P 2 | P 3 ) ∼ g (ν q)(Q 1 | Q 2 | Q 3 ),

and in this case

P σ τ -→ π (νp, p)(P 1 [p/x] | P 2 | P 3 )σ (resp. Qσ τ -→ π (νp, q)(Q 1 [p/x] | Q 2 | Q 3 )σ).
In order to be able to add the restriction on p to the terms given by Lemma 4.10, we rely on the fact that ∼ g is preserved by restriction: P ∼ g Q implies (νp)P ∼ g (νp)Q for any P, Q, p. We can then reason as above to conclude. Corollary 4.12 (Congruence of bisimilarity in π 0 ). In π 0 , ground, early and late bisimilarity coincide and are congruences.

Proof. By a standard argument (see [START_REF] Sangiorgi | The π-calculus: a Theory of Mobile Processes[END_REF]): since ∼ g is closed under substitution, ∼ g is an open bisimulation.

It is known (see [START_REF] Sangiorgi | The π-calculus: a Theory of Mobile Processes[END_REF]) that adding either replication or sum to π 0 yields a calculus where strong bisimilarity fails to be a congruence.

Conclusion

We have presented an axiomatisation of strong bisimilarity on a small subcalculus of CCS, and a new congruence result for the π-calculus.

Technically, the notion of mutual desynchronisation is related to substitution closure of strong bisimilarity, as soon as substitutions can create new interactions by identifying two names.

We have shown in Section 4 that there exists no mutual desynchronisation in π 0 , and that ∼ g is a congruence. In (full) CCS, mutual desynchronisations exist, a simple example being given by a.b+b.a. The latter process is bisimilar to a | b, but the equality fails to hold when b is replaced with a. The same reasoning holds for the π-calculus with choice. It hence appears that in finite calculi, mutual desynchronisations give rise to counterexamples to substitution closure of strong bisimilarity. The situation is less clear when infinite behaviours can be expressed. -→ ≡ P . We do not know at present whether ∼ is substitution-closed in this extension of µCCS (we may remark that the two aforementioned processes remain bisimilar when b is replaced with a).

Some subcalculi of the π-calculus where strong bisimilarity is a congruence are obtained by restricting the output prefix [START_REF] Sangiorgi | The π-calculus: a Theory of Mobile Processes[END_REF]. In the asynchronous π-calculus (Aπ), mutual desynchronisations do not appear, basically because the output action is not a prefix. Strong bisimilarity is a congruence on Aπ. In the private π-calculus (P π), since only private names are emitted, no substitution generated by a synchronisation can identify two previously distinct names. Hence, although mutual desynchronisations exist in P π (due to the presence of the sum operator), strong bisimilarity is not substitution closed, but is a congruence. Indeed, to obtain the latter property, we only need to consider the particular substitutions at work in P π, which cannot identify two names.

The question of substitution closure can also be raised in the framework of location sensitive behavioural equivalences such as distributed bisimilarity (see [START_REF] Castellani | Handbook of Process Algebra, chapter Process Algebras with Localities[END_REF]). Without having a formal proof for this claim, we expect this equivalence to be substitution closed on restriction-free CCS. We believe this should be the case because in absence of restriction, distributed bisimilarity is discriminating enough to analyse the maximum degree of parallelism in processes (in particular, the expansion law is not valid for location sensitive equivalences).

Regarding future extensions of this work, we would like to study whether our approach can be adapted to analyse weak bisimilarity in π 0 (as mentioned in Remark 2.3, strong and weak bisimilarity coincide in µCCS). Another interesting direction, as hinted above, would be to study strong bisimilarity on infinite, restriction-free calculi (in CCS and the π-calculus).

  η.(P | η.P | . . . | η.P ) = η.P | η.P | . . . | η.P

Definition 1 . 4 (

 14 Distribution law). The distribution law is given by the following equation, where the same number of copies of P appears on both sides: η.(P | η.P | . . . | η.P ) = η.P | η.P | . . . | η.P .

Lemma 2 . 5 .

 25 Let A, B be two normal forms, A ∼ B implies A ≡ B.

Lemma 4 . 2 .

 42 Suppose that P σ µ -→ π P . 1. If µ is ab, a(p) or a(x), then P µ -→ π P with µ σ = µ and P σ = P . 2. If µ = τ then one of the three following properties hold, where the input and output actions are offered concurrently by P in the last two cases. (a) P τ -→ π P and P σ = P , (b) P bc -→ π a(x) ---→ π P where σ(a) = σ(b) and P [c/x]σ ∼ P , (c) P b(p) --→ π a(x) ---→ π P where σ(a) = σ(b) and ((νp)P [p/x])σ ∼ P .

Lemma 4 . 5 .

 45 P ∼ Q implies s η (P ) = s η (Q) for all η.

a-→

  P 0 , then there exist x and P such that P 0 = E(P ) and P a(x) ---→ π P . Similarly, if E(P ) b -→ P 0 , there exist c, p, P such that P 0 = E(P ) and either P bc -→ π P or P b(p) --→ π P . Proof. Simple reasoning on the LTSs of µCCS and π 0 .

Proposition 4 . 1 (

 41 Transfer).

  For instance, in the extension of µCCS with replication, the process !a | !b is bisimilar to P def = !a.b | !b.a. Process P leads to a mutual desynchronisation:

  1. A process P is prime if P ∼ 0 and P ∼ P 1 | P 2 implies P 1 ∼ 0 or P 2 ∼ 0. When P ∼ P 1 | . . . | P n where the P i s are prime, we shall call P 1 | . . . | P n a prime decomposition of P . Proposition 2.2 (Unique decomposition). Any process admits a prime decomposition which is unique up to bisimilarity: if P 1 | . . . | P n and Q 1 | . . . | Q m are two prime decompositions of the same process, then n = m and P

  1 | P 2 with P 1 , P 2 ∼ 0. By Lemma 2.4, we have A ∼ (η.B) k with k > 1 and η.B in normal form. By triggering the prefix on the left hand side, we have A ∼ B | (η.B) k-1 . It follows by induction that A ≡ B | (η.B) k-1 (using property (ii)), and hence A ≡ η.(B | (η.B) k-1 , which is in contradiction with the fact that A is in normal form. (ii) Suppose now A ∼ B. • If A is a prefixed process, B is prime by the previous point (#(B) = #(A) by Lemma 1.3). Necessarily, A ≡ η.A and B ≡ η.B with A ∼ B . By induction, this entails A ≡ B , and

Acknowledgements. We are grateful to Arnaud Carayol for interesting discussions at early stages of this work. An anonymous referee provided numerous helpful suggestions, which helped us in particular to improve the proof of Theorem 2.6. We benefited from support by the french initiative "ACI GEOCAL" and from the ANR project "MoDyFiable".


Ghislain Blanchard
email: ghislain.blanchard@onera.fr

Davide Zuzio
email: davide.zuzio@onera.fr

Philippe Villedieu
email: philippe.villedieu@onera.fr

A large scale multi-fluid/dispersed phase approach for spray generation in aeronautical fuel injectors

Keywords: two-fluid model, dispersed phase, multi-scale, atomization modelling

Atomization is usually neglected in industrial numerical simulations of aeronautical combustors, despite the increasing physical fidelity of the LES approaches. The reason comes from the strong multi-scale nature of the problem: the typical space and time lengths are too small for their resolution to be afforded. The influence of the atomization process should however be taken in account when simulating unsteady processes like ignition and combustion instabilities. This work proposes a large scale numerical methodology (LSS, Large Scale Simulation) to simulate the unsteady assisted atomization process in industrial simulations of combustion chambers. The LSS approach involves a two-fluid resolution for the first stages of primary atomization and a dynamic local transfer of the liquid phase to a dispersed phase solver. An atomization model generates volume source terms for the multi-fluid solver which remove the liquid mass, while acting as an advanced adaptive injector for the dispersed phase solver. The methodology has been successfully tested on a planar sheared liquid sheet.

Introduction

In aircraft combustion chambers the injection of fuel is mostly achieved by air-blast atomizers. These devices are meant to pulverize a liquid jet, sheet or thin film with the shearing effect of an incoming high-speed airflow. The numerical simulation of this process is challenging in reason of its strong multi-scale aspect, where several order of magnitude separate the physical phenomena, from the primary atomization of the liquid jet to the formation and evolution of the small droplets. Several studies have been devoted to the development of numerical simulation methods for dealing with this problem. The methods based on the sharp interface methods aim to describe the whole atomization process, from the larger scales to the smaller ones ( [START_REF] Lebas | Numerical simulation of primary break-up and atomization: Dns and modelling study[END_REF], [START_REF] Gorokhovski | Modeling primary atomization[END_REF], [START_REF] Desjardin | An accurate conservative level set/ghost fluid method for simulating turbulent atomization[END_REF]). More ambitious multi-scale approaches include multi-scale adaptive techniques coupled with a local resolved two-phase DNS for the liquid atomization ( [START_REF] Herrmann | The influence of density ratio on the primary atomization of a turbulent liquid jet in crossflow[END_REF], [START_REF] Tomar | Multiscale simulations of primary atomization using gerris[END_REF], [START_REF] Pau | An adaptive mesh refinement algorithm for compressible twophase flow in porous media[END_REF], [START_REF] Sussman | A parallelized, adaptive algorithm for multiphase flows in general geometries[END_REF], [START_REF] Popinet | An accurate adaptive solver for surface tension driven interfacial flows[END_REF], [START_REF] Zuzio | An efficient block parallel amr method for two phase flow simulations[END_REF]). However, their objective remains to directly simulate atomization up to the smallest scale.

In contrast to these works, the objective of the proposed approach is to develop a numerical method allowing the resolution of the larger scale phenomena of the liquid fragmentation mechanism as well as the spray formation and evolution in the context of more industrial-oriented simulations. Three numerical tools are used to achieve this goal:

• a separate phases approach, based on a finite volume discretization of conservation equations for both liquid and gas phases, considered as compressible fluids. This will be referred as the two-fluid model;

• a dispersed phase solver, specifically a Lagrangian solver which allows the simulation of the droplet cloud carried by the gas flow;

• an additional sub-scale atomization model. The model is able to detect zones where the droplets are most likely to appear and to act subsequently as an advanced numerical injector, transferring the resolved liquid to an appropriate distribution of numerical droplets. Taking in account the instabilities of the injection, it gives more accurate initial position and velocity of each particle. The model naturally respects the resolved liquid core and performs the conversion at realistic break-up lengths. The first point is detailed in [START_REF] Blanchard | Numerical simulation of primary atomization of a sheared liquid sheet, 25th ILASS[END_REF]. This paper focuses mainly on the third point, the sub-scale transition model between the two-phase solvers.

Two-fluid model

The two-fluid model considers two non miscible compressible fluids assumed to be simultaneously present at any point in space. The two fluids are supposed to be in local mechanical equilibrium, having locally the same velocity and the same pressure, in a four-equation model:

v = v l = vg (1) p = p l = pg
(2) The mass and momentum balance equations of the model lead to the system:

∂ ρ ∂t + div ( ρ ⊗ v) = 0 (3) ∂ρv ∂t + div (ρv ⊗ v + pI) = div (τ v + τ c) + ρg + sp (4) with ρ = [α l ρ l , αgρg] T = [ ρl , ρg]
T the bulk densities, α l,g the volume fractions of the liquid and the gas, ρ l,g the bulk densities of the liquid and the gas, τ c and τ v respectively the capillary and viscous stress tensor, g the gravity acceleration and sp the two-way coupling source term. The mixture density is ρ = α l ρ l + αgρg. In system (3) + (4) it is worth noting that both fluids play a symmetric role and that the volume fraction of one of the two-phase, e.g. α, is not a primary unknown but rather an outcome of the model given by the closure relationship, resulting from the assumption of local mechanical equilibrium (2):

ρg ρg/α * g = ρ l (ρ l /α * l ) α * l + α * g = 1 (5)
given the equations of state for each phase k:

p k = p 0,k + c 2 k (ρ0 -ρ 0,k ) (6)
The choice of the ( 6) EOS (where c is the speed of sound) is justified by the assumption that density changes are considered as negligible in the low-Mach regimes of the considered applications. The fluid mixture is supposed to behave as a Newtonian fluid for which the general form of the viscous stress tensor reads:

τ v = - 2 3 µ div(v)I + µ ∇v + t ∇v (7)
where µ = α l µ l + αgµg is the dynamic viscosity of the two fluids mixture. The capillary stress tensor is modelled according to the Continuum Surface Stress (CSS) formulation [START_REF] Pau | An adaptive mesh refinement algorithm for compressible twophase flow in porous media[END_REF] which is well suited to diffuse interface method and which is conservative by construction.

τ c = σ I - ∇α l ⊗ ∇α l ∇α l 2 ∇α l (8)
The numerical resolution of the system (3)+( 4) is assured by a partially implicit finite volume scheme on unstructured 3D meshes. An implicit time discretization of the mass conservation equations avoids the acoustic CFL limitation of the time step [START_REF] Grenier | An accurate low-mach scheme for a compressible two fluid model applied to sloshing phenomena[END_REF]. A second-order MUSCL scheme is used to achieve robust second order space accuracy, while a special attention has been paid to the multidimensional slope limitation procedure in order to limit the smearing of the interface [START_REF] Blanchard | Modélisation et simulation multiéchelle de l'atomisation d'une nappe liquide cisaillée[END_REF]..

Dispersed phase model

The spray is modelled with a statistical dispersed phase approach. It is based on the resolution of the scalar function f (p) representing the average liquid droplets density. At a given time, the number of droplets of radius [r, r+δr] moving at a velocity in the interval [v, v+δv] is given by f (p) (t, xp, vp, rp)dxpdvpdrp. The underlying hypothesis on the dispersed phase are the following:

• the droplets are spherical and rigid

• the volume occupied by the droplet is small (diluted spray)

• the influence of the turbulence small scales is neglected The evolution of the function f (p) is given by

∂ ∂t f (p) + ∇x vpf (p) + ∇v Fp mp f (p) = 0 (9)
with Fp the external forces acting on the particles and mp = 4/3 ρpπr 3 p the particle mass (ρp being the liquid density). The external forces Fp acting on the particles are the aerodynamic drag force F p,d and the gravity force Fp,g. The expressions are the following:

F p,d = 1 2 ρ f @p πr 2 p CD v f @p -vp (v f @p -vp) (10)
Fp,g = mp g

where ρ f @p and v f @p are respectively the density and the velocity of the carrier phase at the particle location xp and C d the drag coefficient. From [START_REF] Schiller | A drag coefficient correlation[END_REF], the drag coefficient is :

CD(Rep) = 24(1 + 0.15Re 0.687 p )Re -1 p if Rep < 1000 0.445 otherwise ( 12
)
with

Rep = 2 rpρg v f @p -vp µ -1 g ( 13
)
By defining a characteristic relaxation time in the form

τp = 4ρpdp 3 CD(Rep)ρg v f @p -vp (14)
the drag force acting on the particle can be written as

F p,d = mp (v f @p -vp) τ -1 p ( 15
)
The spray is coupled to the carrier phase in a full two-way coupling. The numerical discretization of equation ( 9) is done by a Lagrangian method, where volume-less numerical parcels carry the information of an arbitrary number of physical droplets. This number is referred to as numerical weight w p,k of the parcel k.

The time evolution of the Np parcels is determined by the following equations:

∀k ∈ [1, Np] ,          d dt (x p,k ) = v p,k d dt (v p,k ) = (v f @p -vp) τ -1 p,k + g w p,k (t) = w p,k (0) (16)
The coupling source term sp appearing in the momentum equation (4) reads

sp = Np k=1 -w p,k F p,k δ (x -x p,k (t)) (17)
A fractional step integration in performed with the systems (3)+( 4) and [START_REF] Déjean | Experimental study on the influence of liquid and air boundary conditions on a planar air-blasted liquid sheet[END_REF]. More details about the numerical integration are available in [START_REF] Blanchard | Modélisation et simulation multiéchelle de l'atomisation d'une nappe liquide cisaillée[END_REF].

Large Scale coupling model

The multi-scale approach proposed in this work is meant to provide a high quality droplet injection obtained by the timeresolved simulation of primary atomization by the multi-fluid solver, and then by employing an atomization model to dynamically generate the droplets. The multi-fluid solver has been successfully in capturing the main features of the primary atomization of a two dimensional planar liquid sheet [START_REF] Blanchard | Numerical simulation of primary atomization of a sheared liquid sheet, 25th ILASS[END_REF].

In this paper, the coupling model will be referred to as the "atomization model" for the sake of simplicity. The model should assure the transfer and the compatibility of the information from one solver (and from one scale) to another. The first issue is that the two-fluid solver is based on a continuous formulation, while the dispersed phase solver is based on a Lagrangian (discrete) formulation. To overcome this issue, the model is based on a set of new variables storing the informations needed for the coupling, which communicate with the two solvers via appropriate source terms. These variables will employ the lower-case [•]a or the upper-case[•] a . This way, the two solvers can coexist in all the domain, the coupling taking action where particular conditions are met.

General model

The system (4)+(3) can be rewritten in a compact form :

∂w ∂t + ∇ • H + ∇ • (C + V) + sg + s d + sa←s (18)
where s d is the particles drag term and sa←s the atomization coupling source term. The system ( 16) is modified by new "injection conditions":

         d dt (x p,k ) = v p,k mp d dt (v p,k ) = F drag + g + injection conditions depending on s d←a (19)
The model variables are defined as watom = t (ρa, (ρv)a), following an evolution in the general form

∂ ∂t (wa) = ss→a + s d→a (20)
with ss→a and s d→a the terms modelling the exchanges from respectively the separate phases and the dispersed phase model towards the atomization model.

Mass source terms, two-fluid towards model

The source term ss→a = t s (ρ)a s→a, s

(ρv)a s→a include the terms for the mass flow rate and the relative momentum which are transferred from the two-fluid model to the atomization model. The mass flux rate is written as:

s (ρ)a s→a = ρl νaRa ( 21
)
where ρl is the liquid mass in the two-fluid model, νa an atomisation frequency (or ν -1 a = τa a characteristic atomization time) and Ra an activation function determining whether the model is active or not. The atomization frequency corresponds to a physical droplet generation rate. No well-established model is able to predict this value. In the proposed strategy the model should use local resolved quantities. The shearing effect on the sheet being the atomization energy source, the local vorticity of the velocity field can be an indicator of the intensity of the breakup:

νa = ∇ × v (22)
The activation term Ra is responsible for the general behaviour of the model, and it should respect several conditions:

• the atomization should occur where the two-fluid solution is too smeared

• the model should not affect resolved liquid structures

• the droplets should not be generated "upwind" of a resolved structure 1 The resulting function has been written as:

Ra = χa(α l , ∇α l ) ψa(∇α l , v) (23)
The χa function depends on the local liquid volume fraction and its gradient, its function is to activate the model by respecting the first and second point above-mentioned.

χa(α l , ∇α l ) = 1 if α l ≤ α a l and ∇α l ≤ ∇α l a 0 otherwise ( 24
)
where α a l and ∇α l a are threshold values for respectively the volume fraction and its gradient in the two-fluid model (they are constant values set at the beginning of the simulation). The constant α a l determines the cells where a small volume of fluid is present: a value of α a l = 0.1 is a good compromise, as it allows to generate the droplets near the break-up point as well as producing a dispersed phase sufficiently diluted. The term ∇α l is directly related to the thickness of the interface:

δ interf ace ≈ ∇α l -1 (25)
Low values of this gradient depict a zone where the interface is smeared, where the liquid is fragmenting into structures too small to be tracked by the two-fluid model. Conversely, high values describe the perimeter of a well-resolved structure even if the local α l is small. Its value being limited for a theoretical sharp interface by ∇α l max • ∆x = 1, the retained value for practical purposes is ∇α l a = C ∇α l ∆x, with C ∇α l = 0.5.

The ψa function contains geometrical considerations. It is defined as:

ψa(∇α l , v) = max (0, cos(θa)) (26)
θa being the angle between the local interface normal (pointing from the liquid to the gas) and the local two-fluid velocity vector:

cos(θa) = - ∇α l • v ∇α l v ( 27
)
Its function is to respect the third point, which is allowing the model to produce droplets mostly where the normal vector is aligned to the velocity vector. Three scenarios can be distinguished by this function:

(a) the interface normal is perpendicular to the local velocity vector. This corresponds to a pure shearing effect of the gas on the liquid. In the targeted applications the pure tearing of droplets is almost never observed, so that the production of droplets is minimized in this case2 .

(b) the interface normal and the local velocity vector are aligned. This case represents the situation where the stretching of the liquid and the ligament mechanism are active, so the droplet production is fully authorized here.

(c) the interface normal and the local velocity vector are opposite. This represents an upwind face of a resolved structure, no atomization is allowed in this case.

Momentum source terms, two-fluid towards model

When the mass is transferred from the two-fluid model to the model, the momentum carried by the liquid is transferred as well in order to assure global conservation. A negative term appears in equation 4 to balance the two-way coupling with the particles. This term depends on the mass source term in equation (21) :

s (ρv)a s→a = s (ρ)a s→a v = ρl νaRav (28)
The velocity used in (28) should be the liquid velocity: however, in the 4 equation model only one velocity is available in each cell. Further developments will be devoted to find a better approximation of v l . Since ρ l ρg, this approximation is acceptable as it implies α l 1/ (1 + ρ l /ρg) ≈ 10 -3 respected by the threshold α a l = 0.1 of equation ( 24). Within the discretization of equations ( 21) and (28), at each time-step each cell cumulates the quantities into the wa vectors following equation (20).

Mass and momentum source terms, model towards twofluid

The source terms seen by the two-fluid solver are defined by global conservation as:

s ( ρl )a a→s = -s (ρ)a s→a (29) s (ρv)a a→s = -s (ρv)a s→a (30)
Equation ( 29) and ( 30) are respectively related to the liquid mass and the global momentum equations.

Pressure conservation

An important step for the removal of the liquid in the twofluid solver is necessary in order to avoid pressure fluctuations. A sudden gas expansion is expected in the newly created empty volume freed by the liquid, as the dispersed phase solver does not take in account the volume of the droplets. A constant pressure condition is imposed, dp l = dpg = 0 ⇐⇒ dρ l = dρg = 0. This implies

dρ l = ρ l dα l (31) dρg = ρgdαg = ρgd (1 -α l ) = -ρgdα l (32)
and then

dρg = - ρg ρ l dρ l (33)
Equation ( 33) means that a mass source term has to be imposed in the gaseous phase equation to "fill" the volume freed by the extracted liquid:

s ( ρg )a a→s = -(ρg/ρ l )s ( ρl )a a→s = (ρg/ρ l )s (ρ)a s→a ( 34
)
Exact conservation for the gas mass is lost; however, the high density ratio between the phases implies a very low added mass to the system. Moreover, the volume ratio between the injected liquid and gas in aeronautical injectors is small. The momentum source term (30) is updated in consequence:

s (ρv) a→s = -(1 -ρg/ρ l ) ρl νaRav (35)

Mass and momentum source terms, dispersed phase towards model

The mass subtracted from the two-fluid model is temporarily stored into the model variables wa. Without any constraint from the discrete model of the dispersed phase solver, at each timestep the mass could be directly injected in form of a numerical droplet:

s d→a = -ss→a (36)
so that wa = const. However, the droplet characteristic formation time (dependent on equation ( 22)) is likely to be larger than the regular time-step. This would imply the generation of a far too large number of droplets, with numerical weights largely inferior to the unity. Instead, using system (20) allows to store the particle mass until the desired numerical weight is attained: the injection characteristic time becomes independent of the two-fluid timestep and the atomization rate. The number of droplet Np,inj|i injected in the cell Ωi at the time-step t n , t n+1 is then defined as the biggest integer number satisfying the following condition:

ma -Np,inj|im (num) p ≥ 0 ( 37
)
where ma is the total liquid mass contained in the model variables:

ma = Ω i ρa t n+1 , x dx (38) and m (num) p
= mpwp the user-defined mass of the numerical particle. In each cell, the atomization model source terms reads:

s (ρ)a d→a = - 1 |Ωi| Np,inj|im (num) p δ t -t n+1 (39)
s (ρv)a d→a = - 1 |Ωi| Np,inj|im (num) p vaδ t -t n+1 (40)
where the model velocity at t = t n+1 is given by

va = (ρv)a/(ρ)a (41)
Equations ( 39) and (40) mean that each time a droplet is injected, its mass and momentum are subtracted from the model variables wa.

Droplets injection

Once the source terms (39) and (40) are defined, the numerical droplets are injected in each affected cell Ωi by respecting:

• the number of injected particles is Np,inj|i

• the initial velocity of each droplet is va

• the initial position is randomly taken inside Ωi

• the properties of the droplets are the same as the liquid

• the diameter of the droplets is a user-defined parameter at the present time

Assisted atomization of a liquid sheet

The multi-scale methodology has been applied to the simulation of the primary atomization of a liquid sheet sheared by two co-flowing air streams. In this configuration, the sheet-shape is unstable: hydrodynamic instabilities are the source of the atomization mechanisms and determine the primary break-up characteristics. Two different moments can be distinguished [START_REF] Stapper | An experimental study of the effects of liquid properties on the breakup of a twodimensional liquid sheet[END_REF], [START_REF] Mansour | Disintegration of liquid sheets[END_REF]. In the primary atomization, the sheet becomes subject to longitudinal instabilities, which are the results of the shearing effect: Kelvin-Helmholtz instabilities perturb the plane sheet, starting a sinusoidal stream-wise oscillation all along the sheet. Then, fully three dimensional instabilities generate transverse modulations. The sheet breaks into smaller liquid packs, ligaments and baglike structures. This continuous fragmentation ends with the formation of the spray of droplets.

Numerical set-up

The considered numerical configuration aims to reproduce a simple atomization device as the one experimentally investigated in [START_REF] Déjean | Experimental study on the influence of liquid and air boundary conditions on a planar air-blasted liquid sheet[END_REF]. The injector consists of a NACA-shaped injector 89 mm long, immersed into a channel airflow, discharging the liquid by a rectangular fence a = 300 µm thick. The liquid used in the reference experience is water. In the present work, the two air ducts are 9 mm thick at each side of the injector. The velocity Vair is the maximum air velocity at this location. Room conditions, T0 ≈ 300 K and P0 = 1 bar, are supposed. The chosen physical parameters are typical of a "stretched ligament" break-up, where the primary atomization manifests in the form of membranes and longitudinal ligaments formation. The momentum ratios involving the "stretched ligament" breakup are comprised between 0.5 < M < 4 (from [START_REF] Mansour | Disintegration of liquid sheets[END_REF] and [START_REF] Lozano | Longitudinal instabilities in a air-blasted liquid sheet[END_REF]). The parameters are summarized in tables 1 and 2 3 . The numerical three-dimensional domain, depicted in figure 3, includes the channel, the injector and a vast chamber in which the sheet atomization occurs without confinement effects. The mesh consists of about 15 M elements extruded in the transverse direction. The mesh includes a Cartesian zone refined around the liquid inlet, where the sheet breaks-up. The cell size is here ∆x f ine = 30 µm. The droplet diameter has been fixed at d drop = 100 µm.

Results

A snapshot of the simulation is depicted in figure 2. The sheet has already attained its steady oscillation regime. The two-fluid model captures the stream-wise oscillation as well as the transverse modulations which induce the break-up. It is able do capture the atomization process up to the formation of the ligaments, as shown in the comparison 2(a). The atomization model is active and correctly performs the coupling after the break-up point, 2(b). The droplets evolves following their initial velocities and the drag force induced by the airflow.

Two macroscopic quantities have been measured, the global oscillation frequency and the mean break-up length. The results are summarized in table 3 in comparison with the experimental values of [START_REF] Déjean | Experimental study on the influence of liquid and air boundary conditions on a planar air-blasted liquid sheet[END_REF]. The results are in quite good agreement with the experiences, despite a slight overestimation of both values.

In order to evaluate the effect of the atomization model, three more two dimensional simulations have been performed on the same configuration, where the air velocity has been set respectively to 30, 50 and 80 m.s -1 . Figure 4 shows the droplet distribution in the three cases, as well as the time averaged droplet volume fraction < αp >. The averaged fields clearly show a decrease in the spray opening with increased airflow, correctly following the different sheet spatial evolution in the difference regimes. A comparison with the experiences is difficult to perform because of the difficulty to employ the same criterion for the visualizations and the simulations.

Conclusion

A new approach is proposed for the simulation of the assisted atomization. This involves coupling of two types of models. The first one, called two-fluid model, is based on the Navier-Stokes equations for two immiscible compressible fluids. It is used to describe the largest scales of the atomization mechanism in the near-injector area. The second one, called spray model, is based on a dispersed phase approach. This model describes the evolution of the droplet cloud produced by the primary fragmentation of liquid jet. The coupling of these two models has been achieved by introducing an atomization model which ensure liquid transfer between the two-fluid model and the spray model. Special care has been taken in order to assure a robust coupling. The approach has been applied to the numerical simulation of a three dimensional sheared liquid sheet: comparison with experimental data have shown that the presented methodology gives promising results.

Figure 1 :

 1 Figure 1: Effect of the ψa function, based on the interface geometry and the local two-fluid velocity.

Figure 2 :

 2 Figure 2: Results of the planar liquid sheet atomization. (a) Visual comparison of the experience (top) and the simulation (bottom) (b) Fully developed atomization with the atomization model in action. In green the two-fluid model resolved liquid (α l = 0.5 iso-contour), in cyan the Lagrangian droplets. From [12].

Figure 3 :

 3 Figure 3: Numerical domain and local mesh.

Table 3 :

 3 Quantitative liquid sheet characteristics, simulation and experience. (a) Vgas = 30 m.s -1 (c) Vgas = 50 m.s -1 (e) Vgas = 80 m.s -1

Figure 4 :

 4 Figure 4: Snapshots of the two-dimensional simulations with different air velocities, atomisation model activated. (Left) Instantaneous fields and (right) mean particle density < αp >.

Table 1 :

 1 Test case parameters, dimensional.

Table 2 :

 2 Test case parameters, non dimensional parameters.

	Non-dimensional parameters
	Air Reynolds (Re δ )	4900
	Water Reynolds (Rea)	366
	Momentum flow rate	1.61
	Weber (W ea)	31

May 22nd -27th 2016, Firenze, Italy ICMF-2016 -9th International Conference on Multiphase Flow May 22nd -27th 2016, Firenze, Italy

As they would travel inside a liquid structure, which would be at odd with the dispersed phase hypothesis

In liquid propulsion injectors, where the gas velocities can attain

m.s -1 , the tearing effect is important: the psia function is different in these applications

The Re air,δ is based on the boundary layer thickness; the Rewater,a and the W ea on the liquid sheet thickness.