Alexandre Chiapolino
email: alexandre.chiapolino@rs2n.eu

Richard Saurel
email: richard.saurel@univ-amu.fr

Extended Noble-Abel-Stiffened-Gas equation of state for sub-and-supercritical liquid-gas systems far from the critical point

Keywords: Noble-Abel, Stiffened-Gas, convexity, two-phase flows, sub-supercritical flows, hyperbolic systems, phase change

) is extended to variable attractive and repulsive effects to improve the liquid phase accuracy when large temperature and pressure variation ranges are under consideration. The transition from pure phase to supercritical state is of interest as well. The gas phase is considered through the ideal gas assumption with variable specific heat rendering the formulation valid for high temperatures. The liquid equation-of-state constants are determined through the saturation curves making the formulation suitable for two-phase mixtures at thermodynamic equilibrium. The overall formulation is compared to experimental characteristic curves of the phase diagram showing good agreement for various fluids (water, oxygen). Compared to existing cubic equations of state the present one is convex, a key feature for computations with hyperbolic flow models.

Introduction

Modeling liquid-gas systems with or without phase transition is an old research topic in the physics community but still challenging at both theoretical and computational levels. The most common thermodynamical approach relies on cubic equations of state (EOS), the van der Waals one being the basic prototype. Indeed, this EOS involves all relevant molecular effects present in matter, i.e., thermal agitation, short distance repulsive forces and long range attractive ones. It is thus able to deal, at least qualitatively, with pure liquid, pure gas and two-phase mixture. This EOS, as all cubic ones, is aimed to close flow models based on balance equations of mass, momentum and energy for the mixture. The Euler equations are one of the relevant possible options, as well as more sophisticated ones aimed to model capillary effects, such as the [START_REF] Cahn | Free energy of a nonuniform system. I. Interfacial free energy[END_REF] [1] model for example. In this context, the thermodynamical state is determined from two internal variables only, the density and the internal energy of the mixture, or alternatively the density and the temperature, depending on the formulation of the equations. This approach consequently seems simple, but involves serious difficulties and limitations:

-The first and certainly the most obvious and limiting is related to its inability to deal with liquid and non-condensable gas separated by well-defined interfaces, such as for example interfacial flows of liquid water and air. The thermodynamics of these two media being considered as discontinuous, specific theoretical and numerical treatments have been addressed. In this context, Arbitrary Lagrangian Eulerian [START_REF] Hirt | An arbitrary Lagrangian-Eulerian computing method for all flow speeds[END_REF] [START_REF] Hirt | An arbitrary Lagrangian-Eulerian computing method for all flow speeds[END_REF]), Interface Reconstruction [START_REF] Youngs | An interface tracking method for a 3D Eulerian hydrodynamics code[END_REF] [START_REF] Youngs | An interface tracking method for a 3D Eulerian hydrodynamics code[END_REF]), Front Tracking [START_REF] Glimm | Three-dimensional front tracking[END_REF] [START_REF] Glimm | Three-dimensional front tracking[END_REF]), Level-Set [START_REF] Fedkiw | A non-oscillatory Eulerian approach to interfaces in multimaterial flows (the Ghost Fluid Method)[END_REF] [5]), anti-diffusion (Kokh and Lagoutiere (2010) [START_REF] Kokh | An anti-diffusive numerical scheme for the simulation of interfaces between compressible fluids by means of a five-equation model[END_REF]) methods are possible options. Another approach relies on continuous models with extra internal variables, such as volume and mass fractions and extended equation of state. Examples of such models are the Kapila et al. ( 2001) [START_REF] Kapila | Two-phase modeling of deflagration-to-detonation transition in granular materials: Reduced equations[END_REF] one and its extension with phase transition [START_REF] Saurel | Modelling phase transition in metastable liquids: application to cavitating and flashing flows[END_REF] [START_REF] Saurel | Modelling phase transition in metastable liquids: application to cavitating and flashing flows[END_REF]) to cite a few. With these formulations the same equations are solved everywhere routinely, in pure liquid, pure gas and interface which becomes a diffuse zone. These models are indeed often named "diffuse interface methods" (Saurel and Pantano (2018) [START_REF] Saurel | Diffuse Interfaces and Capturing Methods in Compressible Two-Phase Flows[END_REF]). In this approach, hyperbolic models with relaxation are considered and each phase evolves in its own volume, with its own thermodynamics. In particular there is no need to address cubic formulations. When phase transition is addressed, it occurs through mass transfer terms that can be considered finite rate [START_REF] Saurel | Modelling phase transition in metastable liquids: application to cavitating and flashing flows[END_REF] [8], Furfaro and Saurel (2016) [START_REF] Furfaro | Modeling droplet phase change in the presence of a multi-component gas mixture[END_REF]) or assumed stiff when the physical knowledge of the phase change kinetics is not enough documented (Le [START_REF] Métayer | Dynamic relaxation processes in compressible multiphase flows[END_REF] [START_REF] Métayer | Dynamic relaxation processes in compressible multiphase flows[END_REF], Chiapolino et al. (2017) [START_REF] Chiapolino | A simple phase transition relaxation solver for liquid-vapor flows[END_REF], [START_REF] Chiapolino | A simple and fast phase transition relaxation solver for compressible multicomponent two-phase flows[END_REF]) or unnecessary. -The second limitation is related to the lack on convexity of cubic EOSs, having dramatic consequences on sound propagation during phase transition. The square sound speed becomes negative in the spinodal decomposition zone, such behavior being unphysical. -The third limitation is related to the description of phase transition with such EOSs. Cubic equations of state consider phase transition as a thermodynamic process and not a kinetic one.

It is unclear at this level whether cubic EOSs are limited to the description of global two-phase mixtures with many interfaces and not local ones, at the scale of a single interface. -The fourth, but possibly not the last, is related to the numerical treatment of boundary conditions (BC) in practical compressible flow computations. Subsonic inflow and outflow BCs rely on stagnation enthalpy and entropy invariance coupled to Riemann invariants that can be defined and computed correctly only if the equation of state is well-posed. The second issue related to EOS convexity consequently reemerges at this level. Moreover, the practical expression of Riemann invariants may be inextricable with these EOSs.

This list of arguments gives motivations to the present work where an extended version of the Noble-Abel-Stiffened-Gas (NASG, Le Métayer and Saurel (2016) [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF]) EOS is examined to:

-Represent the thermodynamics of pure liquid, pure vapor and supercritical fluid. Combination of the pure liquid and pure vapor EOSs must be able to represent as accurately as possible the two-phase region. -Each phase EOS must be convex in its respective domain.

-The EOS must be as simple as possible, while remaining accurate, to simplify practical computations and building of mixture EOS in hyperbolic multiphase flow models.

Hyperbolic multiphase flow models have demonstrated their ability to solve a wide range of complex flow situations in severe conditions. Material interface problems [START_REF] Chiapolino | Sharpening diffuse interfaces with compressible fluids on unstructured meshes[END_REF], chemical reactions [START_REF] Petitpas | Modelling detonation waves in condensed energetic materials: Multiphase CJ conditions and multidimensional computations[END_REF], phase change [START_REF] Saurel | A general formulation for cavitating, boiling and evaporating flows[END_REF], surface tension [START_REF] Perigaud | A compressible flow model with capillary effects[END_REF], solid-fluid [START_REF] Favrie | Solid-fluid diffuse interface model in cases of extreme deformations[END_REF], plastic transformation [START_REF] Ndanou | Multi-solid and multi-fluid diffuse interface model: Applications to dynamic fracture and fragmentation[END_REF], dense and dilute flows [START_REF] Saurel | Modelling compressible dense and dilute two-phase flows[END_REF], shallow water flows [START_REF] Chiapolino | Models and methods for two-layer shallow water flows[END_REF] can be cited for instance. In these flow models, compressibility of each phase is responsible for the hyperbolic character of the equations and an appropriate and convex EOS is required for each fluid.

The NASG EOS combines relevant physics and simplicity. Its predictions are in a good agreement with experimental data but in restricted temperature range, [300 -500] K for example with liquid water at saturation. This limitation is linked to constant attractive and repulsive effects. Indeed, this assumption no longer holds when larger pressure and temperature ranges are addressed.

The present contribution aims at extending the liquid NASG EOS to variable attractive and repulsive effects to improve its range of validity, a necessary improvement in view of future engineering applications.

The thermodynamics of the liquid and vapor phases must be combined correctly to reproduce the phase diagram and relevant properties, such as the latent heat of phase change and saturation pressure for example. Building of appropriate EOSs in this direction was done by Le [START_REF] Métayer | Elaborating equations of state of a liquid and its vapor for two-phase flow models; Élaboration des lois d'état d'un liquide et de sa vapeur pour les modèles d'écoulements diphasiques[END_REF] [START_REF] Métayer | Elaborating equations of state of a liquid and its vapor for two-phase flow models; Élaboration des lois d'état d'un liquide et de sa vapeur pour les modèles d'écoulements diphasiques[END_REF] and Le [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF] [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF] with extended "Stiffened-Gas" (SG) formulations.

In the present contribution, the proposed EOS is also meant to describe transitions from pure fluids to supercritical state.

Such transition is indeed essential in some industrial applications such as flows evolving in combustion chambers of cryotechnic rocket engines as well as combustion systems of modern automotive engines for the sake of reduced pollutant emissions and fuel consumption.

The determination of the corresponding EOS parameters is of interest as well. For gases, the new formulation reduces to the ideal gas description. Variable heat capacities can easily be considered with such formulation making the equation of state able to deal with supercritical fluids at high temperatures.

This paper is organized as follows. The determination of the novel liquid EOS is described in Section 2. Among the different relations, the Gibbs free energies of the liquid-vapor couple provide the saturation conditions. Those latter ones are developed in Section 3 and the overall formulation is summarized in Section 4. Experimental and theoretical curves are compared in Section 5, considering water and oxygen at saturation. The abilities of the proposed EOS are illustrated in Section 6 with transitions from sub to supercritical state. Practical applications are illustrated in Section 7.

Extended NASG EOS

The following Extended NASG EOS (ENASG) is considered as a postulate,

p(e, v) = (γ -1)(e -q) v -b(v) -γp ∞ (e, v), (1)
with p ∞ (e, v) = p ∞ (T) = p ∞,1 T + p ∞,0 and b(v) = b 1 v + b 0 , (2)
where p, T, v, e and q represent respectively, the pressure, temperature, specific volume, specific internal energy and reference energy of a corresponding single phase fluid. q and γ are parameters considered as constant coefficients and are meant to be characteristics of the thermodynamic properties of the fluid. b(v) represents the covolume, modeling short range repulsive effects. The term γp ∞ (T) represents attractive ones, present in condensed matter only.

In the NASG formulation, the parameters b and p ∞ are considered constant, yielding simplicity while ensuring presence of the main molecular forces present in a fluid. Besides, as the formulation is close to the ideal gas expression, it facilitates the resolution of the Riemann problem [START_REF] Plohr | Shockless acceleration of thin plates modeled by a tracked random choice method[END_REF] [24], [START_REF] Menikoff | The Riemann problem for fluid flow of real materials[END_REF] [START_REF] Menikoff | The Riemann problem for fluid flow of real materials[END_REF], [START_REF] Cocchi | A Riemann problem based method for the resolution of compressible multimaterial flows[END_REF] [START_REF] Cocchi | A Riemann problem based method for the resolution of compressible multimaterial flows[END_REF]). The Riemann problem is indeed the cornerstone of numerical methods used to solve hydrodynamic problem, see [START_REF] Toro | Riemann solvers and numerical methods for fluid dynamics: A practical introduction[END_REF] [START_REF] Toro | Riemann solvers and numerical methods for fluid dynamics: A practical introduction[END_REF] for example.

The simplicity of the corresponding formulas is beneficial to the theoretical analysis and computational efficiency. This section aims at extending the liquid NASG EOS to deal with large pressure and temperature variations while remaining simple and convex.

In this work, simple linear dependencies on the specific volume and temperature have been added to the NASG EOS regarding respectively the covolume b(v) and attractive pressure p ∞ (T). The formulation is meant to be convex and thermodynamically consistent. These two points are addressed hereafter. Inverting Eq. ( 1), the internal energy reads,

e(v, T) = p(v, T) + γp ∞ (T) γ -1 v -b(v) + q. ( 3
)
From postulate (1) or its alternative (3) form, the aim is now to derive the thermal EOS p = p(v, T) and the caloric one e = e(v, T).

Thermal and caloric EOSs

Thermal and caloric EOSs must fulfill the compatibility condition,

∂ ∂v ∂ f ∂T v T = ∂ ∂T ∂ f ∂v T v , (4)
that is precisely the first Maxwell's relation, where f represents the Helmholtz free energy defined by,

f = e -Ts,
where s denotes the specific entropy. With the help of the thermodynamic definition of pressure and

entropy, p = -∂ f ∂v T and s = -∂ f ∂T v
, identity (4) transforms to a more convenient expression linking the thermal EOS p(v, T) and the caloric one e(v, T),

∂e ∂v T = T ∂p ∂T v -p. (5)
From postulate (3), the following partial derivatives arise,

∂e ∂T v = v -b(v) γ -1 ∂p ∂T v + γp ∞,1 , (6)
∂e ∂v T = v -b(v) γ -1 ∂p ∂v T + 1 -b 1 γ -1 p + γp ∞ (T) . ( 7
)
As the thermal capacity at constant volume is defined as,

C v = ∂e ∂T v , (8)
combining Eqs. ( 6) and ( 8) results in,

∂p ∂T v = (γ -1)C v v -b(v) -γp ∞,1 . (9)
The preceding relation ( 9) is now integrated over the temperature T leading to,

p(v, T) = (γ -1)C v T v -b(v) -γp ∞,1 T + K(v), (10)
where K(v) is a function depending on the specific volume v. Expression ( 10) is differentiated over v and at constant temperature T yielding,

∂p ∂v T = - (1 -b 1 )(γ -1)C v T v -b(v) 2 + dK(v) dv . ( 11
)
Afterwards, relation ( 7) is inserted into Maxwell's relation [START_REF] Fedkiw | A non-oscillatory Eulerian approach to interfaces in multimaterial flows (the Ghost Fluid Method)[END_REF] resulting in,

∂p ∂v T = (γ -1)T v -b(v) ∂p ∂T v - γ p + p ∞ (T) v -b(v) + b 1 p + p ∞ (T) v -b(v) . ( 12
)
Expressions ( 9) and ( 10) are now introduced into relation (12) leading to,

∂p ∂v T = - (γ -1)C v T(1 -b 1 ) v -b(v) 2 - K(v)(γ -b 1 ) + γp ∞ (T)(1 -b 1 ) v -b(v) + p ∞,1 γT(1 -b 1 ) v -b(v) . ( 13
)
The equality between Eqs. ( 11) and ( 13) yields a first-order ordinary differential equation,

dK(v) dv + K(v)(γ -b 1 ) + γp ∞,0 (1 -b 1 ) v -b(v) = 0. ( 14
)
The solution of Eq. ( 14) is given by,

K(v) = cst (γ -b 1 ) v -b(v) γ-b 1 1-b 1 - γp ∞,0 (1 -b 1 ) γ -b 1 . ( 15
)
Inserting Eq. ( 15) into [START_REF] Furfaro | Modeling droplet phase change in the presence of a multi-component gas mixture[END_REF], the thermal equation of state reads,

p(v, T) = (γ -1)C v T v -b(v) -p ′ ∞ (T) - d v -b(v) γ-b 1 1-b 1 , (16)
where d = -cst/(γb 1 ) is a constant to be determined and the "attractive" pressure p ′ ∞ (T) is defined as,

p ′ ∞ (T) = γp ∞,1 T + γp ∞,0 (1 -b 1 ) γ -b 1 . ( 17
)
In the present approach, as in Le Métayer et al. ( 2004) [START_REF] Métayer | Elaborating equations of state of a liquid and its vapor for two-phase flow models; Élaboration des lois d'état d'un liquide et de sa vapeur pour les modèles d'écoulements diphasiques[END_REF] and Le Métayer and Saurel (2016) [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF], each fluid, liquid and gas, is governed by the same EOS, here Eq. ( 16) but with different parameters unlike cubic EOSs. The term p ′ ∞ (T) is important for the liquid state whereas the second attractive term

d/ v -b(v) γ-b 1 1-b 1 ,
reminiscent of cubic EOSs, is important for dense gases. However, this coefficient yields conditional convexity (see Appendix E). The same observation holds for cubic EOSs. As this section aims to build an unambiguously convex EOS, the parameter d is set to zero. Note that d = 0 is a particular solution of Eq. ( 14). The corresponding ENASG EOSs then read,

p(v, T) = (γ -1)C v T v -b(v) -p ′ ∞ (T), (18)
e(v, T)

= C v T + γp ∞,0 v -b(v) γ -b 1 + q. ( 19
)
With the help of the caloric EOS [START_REF] Favrie | Solid-fluid diffuse interface model in cases of extreme deformations[END_REF], the temperature expresses as,

T(e, v) = e -q C v - γp ∞,0 v -b(v) C v (γ -b 1 ) , (20)
and yields,

p(e, v) = (γ -1)(e -q) v -b(v) -γp ∞ T(e, v) . ( 21
)
Note that Eqs. ( 18), (19), (20) and ( 21) reduce to the NASG expressions if p ∞,1 = 0 and b 1 = 0. Also, inserting Eq. ( 18) into Eq. ( 21), the internal energy expresses,

e(p, T) = p + γp ∞ (T) p + p ′ ∞ (T) C v T + q. ( 22
)
The expressions of the thermal and caloric EOSs being now available, the other thermodynamic variables may be obtained from the knowledge of the two independent variables p and T. This task is addressed hereafter.

Expression of the entropy

Expression of the specific entropy is mandatory to express the Gibbs free energy, a key function to address phase transition. The entropy formulation must fulfill the compatibility relation,

∂ ∂p ∂g ∂T p T = ∂ ∂T ∂g ∂p T p , (23)
that is precisely the second Maxwell's relation, where g represents the Gibbs free energy defined by, , identity [START_REF] Métayer | Elaborating equations of state of a liquid and its vapor for two-phase flow models; Élaboration des lois d'état d'un liquide et de sa vapeur pour les modèles d'écoulements diphasiques[END_REF] transforms to a more convenient expression,

g = h -Ts,
∂s ∂p T = - ∂v ∂T p . ( 24
)
With the help of Eq. ( 18), the partial derivative expresses,

∂v ∂T p = - 1 1 -b 1 -(γ -1)C v p + p ′ ∞ (T) + γp ∞,1 (γ -1)C v T p + p ′ ∞ (T) 2 . ( 25
)
Also, by use of Maxwell's rule [START_REF] Plohr | Shockless acceleration of thin plates modeled by a tracked random choice method[END_REF], the next equations arise,

∂s ∂p T = 1 1 -b 1 -(γ -1)C v p + p ′ ∞ (T) + γp ∞,1 (γ -1)C v T p + p ′ ∞ (T) 2 , (26)
s(p, T) = 1 1 -b 1 -(γ -1)C v ln p + p ′ ∞ (T) - γp ∞,1 (γ -1)C v T p + p ′ ∞ (T) + K(T). (27)
Equation ( 18) is now inserted into Eq. ( 27) and yields,

s(v, T) = 1 1 -b 1 -(γ -1)C v ln (γ -1)C v T v -b(v) -γp ∞,1 v -b(v) + K(T). (28)
This last equation admits the partial derivative,

∂s ∂T v = - (γ -1)C v T(1 -b 1 ) + dK(T) dT . ( 29
)
The definition of the thermal capacity at constant volume may be used under the following form,

∂s ∂T v = C v T . ( 30
)
Consequently, the next equation arises,

dK(T) = C v dT T + (γ -1)C v (1 -b 1 ) dT T ,
and is directly integrated yielding a temperature-dependent function K(T),

K(T) = C v ln(T) + (γ -1)C v (1 -b 1 ) ln(T) + q ′ , ( 31
)
where q ′ is defined as a constant (reference entropy). Equation ( 31) is now embedded in Eq. [START_REF] Mcbride | Coefficients for Calculating Thermodynamic and Transport Properties of Individual Species[END_REF].

After some algebraic manipulations, the resulting equation reads,

s(v, T) = C v ln(T) - γp ∞,1 v -b(v) 1 -b 1 + (γ -1)C v (1 -b 1 ) ln v -b(v) - (γ -1)C v (1 -b 1 ) ln (γ -1)C v + q ′ .
As the last term of this equation is constant, it is convenient to define,

q ′′ = - (γ -1)C v (1 -b 1 ) ln (γ -1)C v + q ′ . ( 32
)
The entropy equation consequently reads,

s(v, T) = C v ln(T) + (γ -1)C v (1 -b 1 ) ln v -b(v) - γp ∞,1 v -b(v) 1 -b 1 + q ′′ . ( 33
)
Obviously, definition (30) is satisfied. Equation ( 18) is now inserted into (33), yielding,

s(p, T) = C v ln(T) + (γ -1)C v (1 -b 1 ) ln (γ -1)C v T p + p ′ ∞ (T) - γp ∞,1 (γ -1)C v T 1 -b 1 p + p ′ ∞ (T) + q ′′ . ( 34
)
The relation s(p, T) being now available, it is worth analyzing the expression of heat capacity at constant pressure. Equation ( 34) admits as partial derivative,

∂s ∂T p = C v T + (γ -1)C v (1 -b 1 )T - γp ∞,1 (γ -1)C v 1 -b 1 p + p ′ ∞ (T) 1 + p + p ′ ∞ (T) -γp ∞,1 T p + p ′ ∞ (T) . ( 35
)
As, the heat capacity is defined as,

∂s ∂T p = C p T , (36)
the ENASG thermal capacity at constant pressure consequently reads,

C p = C v + (γ -1)C v (1 -b 1 ) - γp ∞,1 (γ -1)C v T 1 -b 1 p + p ′ ∞ (T) 1 + p + p ′ ∞ (T) -γp ∞,1 T p + p ′ ∞ (T) . ( 37
)
It then appears that C p = γC v . However, if p ∞,1 = 0 and b 1 = 0, then the preceding relation reduces to,

C p = C v + (γ -1)C v = γC v ,
and the NASG thermal capacity at constant pressure is recovered. Note also that this feature is valid for the Stiffened-Gas (SG) (b 0 = 0) and ideal gas (p ∞,0 = 0) equations of state as well. Equation (34) can be manipulated to obtain an entropy relation closer to the NASG one. Indeed, after some algebraic manipulations, Eq. ( 34) can be written as,

s(p, T) = C v ln(T) + ln T p + p ′ ∞ (T) γ-1 1-b 1 - γp ∞,1 (γ -1)C v T 1 -b 1 p + p ′ ∞ (T) + (γ -1)C v 1 -b 1 ln (γ -1)C v + q ′′ .
Using Eq. ( 32), the last term of this equation reduces to,

(γ -1)C v 1 -b 1 ln (γ -1)C v + q ′′ = q ′ .
Equation [START_REF] Callen | an introduction to the physical theories of equilibrium thermostatics and irreversible thermodynamics[END_REF] consequently transforms to,

s(p, T) = C v ln   T γ-b 1 1-b 1 p + p ′ ∞ (T) γ-1 1-b 1   - γp ∞,1 (γ -1)C v T 1 -b 1 p + p ′ ∞ (T) + q ′ . ( 38
)
Under form (38), it is straightforward to see that the relation does reduce to the NASG equation if p ∞,1 = 0 and b 1 = 0. At this point, the caloric, thermal and entropy equations of state are determined. The next step is to check positivity of the sound speed as it is a key feature in fluid dynamics.

Speed of sound

The sound speed is defined as,

c 2 = -v 2 ∂p ∂v s . ( 39
)
The pressure is expressed as a function of the specific volume and the specific entropy by combining relations ( 18) and (38),

p(v, s) = exp s-q ′′ C v exp γp ∞,1 v-b(v) C v (1-b 1 ) (γ -1)C v -γp ∞,1 v -b(v) v -b(v) γ-b 1 1-b 1 - γp ∞,0 (1 -b 1 ) γ -b 1 , (40)
with

q ′′ = - (γ -1)C v (1 -b 1 ) ln (γ -1)C v + q ′ . ( 41
)
The ENASG speed of sound consequently expresses after some algebraic manipulations as,

c 2 (p, v) = -v 2 γp ∞,1 p + γp ∞,0 (1 -b 1 ) γ -b 1 γ -1 (γ -1)C v -γp ∞,1 v -b(v) + 1 C v +   p + γp ∞,0 (1-b 1 ) γ-b 1 (γ -1)C v -γp ∞,1 v -b(v)   v 2 (γ -b 1 )(γ -1)C v v -b(v) . ( 42
)
It is worth mentioning that Eq. (42) reduces to,

c 2 (p, v) = v 2 (p + p ∞,0 )γ v -b 0 , (43)
if p ∞,1 = 0 and b 1 = 0, that corresponds to the NASG speed of sound. It also appears that

c 2 > 0 unambiguously if p ∞,1 ≤ 0, p ∞,0 ≥ 0 and b 1 < 1. Obviously v -b(v) must be positive as well.
For a liquid state, the attractive effects summarized by the terms p ∞ (T) and p ′ ∞ (T) are expected to decrease when the temperature rises. The conditions p ∞,1 ≤ 0 and p ∞,0 ≥ 0 are consequently in agreement with the physics to represent, at least qualitatively. The same observation holds for the repulsive effects summarized by the covolume b(v). When the density decreases, those short distance effects are expected to vanish as the liquid tends to become a dense gas. Condition b 1 < 1 is then not restrictive for liquids. Also, it is worth noticing from Eq. (42), that

p ∞,1 = γ -1 v -b(v) C v γ ,
must be satisfied for the ENASG sound speed function to be defined. However, as the right-hand side of this relation is necessarily positive, considering p ∞,1 ≤ 0 satisfies unambiguously this condition. The two linear dependencies p ∞ (T) and b(v) are thus in agreement with the description of a liquid state. As it will be seen further, those simple functions result in predictions in very good agreement with experimental data. This set of liquid EOS relationships results in a convex formulation, this feature being essential for both theoretical and numerical considerations. Indeed, sufficient conditions to ensure convexity are summarized by,

p ∞,1 ≤ 0, p ∞,0 ≥ 0 and b 1 < 1, (44)
related calculations being given in Appendix A. For the gas phase, the attractive effects are expected to increase with the temperature but this evolution is not in agreement with the convexity condition (44). Gas attractive effects are thus removed by setting p ∞,0 = 0 and p ∞,1 = 0, reducing the formulation to the Noble-Abel (NA) EOS with variable covolume b(v).

However, covolume effects alone are not enough to describe dense gases near the critical point. Attractive effects are needed in addition (see Appendix E) but yield conditional convexity. As it will be seen further, the ideal gas EOS is well-suited for fluids evolving away from the critical point, either at low temperatures where thermal capacities can be considered constant or at much higher ones where heat capacities are meant to evolve with the temperature.

Thereby, for the sake of simplicity, covolume effects will be removed as well reducing the formulation to the ideal gas description. Consequently, the saturated vapor phase lacks of accuracy near the critical point since attractive effects are absent but the overall formulation remains convex, a key feature for computational fluid dynamics.

For ENASG formulation completion, the expression of the saturation conditions of the liquid-vapor couple must be determined. This task is addressed in the next section.

Saturation condition of the liquid-vapor couple

Thermodynamic equilibrium is considered when the fluids are in pressure, temperature and Gibbs free energies (g l = g v ) equilibrium. The saturation condition results from these equilibria. As,

g(p, T) = h(p, T) -Ts(p, T), (45)
the enthalpies of pure constituents must be determined first. The enthalpy is defined as,

h(p, T) = e(p, T) + pv(p, T). (46)
Combining Eqs. ( 18) and ( 21), the enthalpy of the ENASG EOS reads,

h(p, T) = C v T p + p ′ ∞ (T) γ p + p ∞ (T) -pb 1 -γb 1 p ∞ (T) 1 -b 1 + pb 0 1 -b 1 + q. ( 47
)
Note that the partial derivative of Eq. ( 47) gives after some algebraic manipulations,

∂h ∂T p = C v + (γ -1)C v (1 -b 1 ) - γp ∞,1 (γ -1)C v T 1 -b 1 p + p ′ ∞ (T) 1 + p + p ′ ∞ (T) -γp ∞,1 T p + p ′ ∞ (T) . ( 48
)
Equation ( 37) is then recovered and the definition of the thermal capacity at constant pressure is satisfied, providing extra verification,

C p = T ∂s ∂T p = ∂h ∂T p . ( 49
)
Note that if p ∞,1 ≤ 0, b 1 < 1 and γ > 1, then the thermodynamic condition C p > C v is ensured. These conditions are the same that preserve convexity of the formulation. Thanks to Eqs. ( 38), ( 45) and (47), the Gibbs free energy of a pure constituent expresses as,

g(p, T) = C v 1 -b 1 γ p + p ∞ (T) -pb 1 -γb 1 p ∞ (T) p + p ′ ∞ (T) -q ′ T -C v T ln   T γ-b 1 1-b 1 p + p ′ ∞ (T) γ-1 1-b 1   + pb 0 1 -b 1 + q + γp ∞,1 (γ -1)C v T 2 1 -b 1 p + p ′ ∞ (T)
.

(50) Solution of the equation,

g l (p, T) = g v (p, T), (51)
provides the saturation pressure as a function of temperature p sat (T). Subscripts l and v denote respectively the liquid and vapor states.

The equality of Gibbs free energies of both phases corresponds to phase equilibrium and leads to the following expression linking the pressure and temperature,

ln p + p ′ ∞,v (T) = C v,l (1 -b 1,v ) (1 -b 1,l )(γ v -1)C v,v γ l p + p ∞,l (T) -pb 1,l -γ l b 1,l p ∞,l (T) p + p ′ ∞,l (T) - 1 γ v -1 γ v p + p ∞,v (T) -pb 1,v -γ v b 1,v p ∞,v (T) p + p ′ ∞,v (T) + (q ′ v -q ′ l )(1 -b 1,v ) (γ v -1)C v,v B + Ep T + C ln(T) + D ln p + p ′ ∞,l (T) -T γ v p ∞,1,v p + p ′ ∞,v (T) - γ l p ∞,1,l D p + p ′ ∞,l (T) , (52) with
B = (q l -q v )(1 -b 1,v ) (γ v -1)C v,v , C = (γ v -b 1,v )C v,v 1 -b 1,v - (γ l -b 1,l )C v,l 1 -b 1,l 1 -b 1,v (γ v -1)C v,v , D = (γ l -1)C v,l (1 -b 1,v ) (γ v -1)C v,v (1 -b 1,l ) , E = b 0,l 1 -b 1,l - b 0,v 1 -b 1,v 1 -b 1,v (γ v -1)C v,v . (53)
Relation (52) provides a unique value of the pressure for a given temperature and implicitly represents the theoretical saturated pressure as a function depending on the temperature. Numerical resolution is needed to compare the predictions with experiments, as will be examined later. When p ∞,1,k and b 1,k are set to zero, the preceding relation reduces to,

ln p + p ∞,0,v = A + B + Ep T + C ln(T) + D ln p + p ∞,0,l , (54)
with

A = γ l C v,l -γ v C v,v + q ′ v -q ′ l γ v C v,v -C v,v , B = q l -q v γ v C v,v -C v,v , C = γ v C v,v -γ l C v,l γ v C v,v -C v,v , D = γ l C v,l -C v,l γ v C v,v -C v,v , E = b 0,l -b 0,v γ v C v,v -C v,v . (55)
The NASG relation is then recovered. The whole ENASG formulation is summarized in the next section and is compared with experimental data in the following ones.

Summary of the Extended NASG state functions

The different liquid ENASG functions of common use are,

                                                                                         p(e, v) = (γ -1)(e -q) v -b(v) -γ p ∞ T(e, v) , v(p, T) = (γ -1)C v T (1 -b 1 )(p + p ′ ∞ (T)) + b 0 (1 -b 1 ) , e(p, T) = p + γ p ∞ (T) p + p ′ ∞ (T) C v T + q, h(p, T) = C v T p + p ′ ∞ (T) γ p + p ∞ (T) -pb 1 -γb 1 p ∞ (T) 1 -b 1 + pb 0 1 -b 1 + q, s(p, T) = C v ln   T γ-b 1 1-b 1 p + p ′ ∞ (T) γ-1 1-b 1   - γp ∞,1 (γ -1)C v T 1 -b 1 p + p ′ ∞ (T) + q ′ , g(p, T) = C v 1 -b 1 γ p + p ∞ (T) -pb 1 -γb 1 p ∞ (T) p + p ′ ∞ (T) -q ′ T -C v T ln   T γ-b 1 1-b 1 p + p ′ ∞ (T) γ-1 1-b 1   + pb 0 1 -b 1 + q + γp ∞,1 (γ -1)C v T 2 1 -b 1 p + p ′ ∞ (T) , c 2 (p, v) = -v 2 γ p ∞,1 p + γ p ∞,0 (1 -b 1 ) γ -b 1 γ -1 (γ -1)C v -γ p ∞,1 v -b(v) + 1 C v +   p + γ p ∞,0 (1-b 1 ) γ-b 1 (γ -1)C v -γ p ∞,1 v -b(v)   v 2 (γ -b 1 )(γ -1)C v v -b(v) , (56) with
T(e, v) = e -q C v - γ p ∞,0 v -b(v) C v (γ -b 1 ) , p ∞ (e, v) = p ∞ (T) = p ∞,1 T + p ∞,0 , p ′ ∞ (T) = γp ∞,1 T + γp ∞,0 (1 -b 1 ) γ -b 1 and b(v) = b 1 v + b 0 .
Those different functions are in agreement with the fundamental relations of Maxwell analyzed in Appendix B and are thermodynamically consistent and convex under conditions p ∞,1 ≤ 0, p ∞,0 ≥ 0 and b 1 < 1. The whole formulation reduces to the NASG EOS if p ∞,1 = 0 and b 1 = 0. In addition, the ideal gas description is recovered if p ∞,0 = 0, b 0 = 0 and is used for the sake of convexity and simplicity of the gas-phase formulation. The corresponding formulas thus read,

                                         p(e, v) = (γ -1)(e -q) v , v(p, T) = (γ -1)C v T p , e(T) = C v T + q, h(T) = γC v T + q, s(p, T) = C v ln T γ p γ-1 + q ′ , g(p, T) = γC v -q ′ T -C v T ln T γ p γ-1 + q, c 2 (p, v) = γpv. (57)
In these formulations (ENASG and ideal gas), heat capacities are considered constant. This assumption is fair for the liquid phase. It is also valid for the gas phase evolving at low temperatures. However, this assumption fails at high temperatures motivating consideration of variable heat capacities as introduced in Section 6. This situation is typical of supercritical fluids at high temperatures.

Extended NASG parameters

The method used in this work to determine the different EOS parameters is summarized in Appendix C. The liquid parameters are computed with the experimental saturation curve as in Le Métayer and Saurel (2016) [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF], but unlike this last reference the gas parameters are chosen regardless of the saturation conditions. The present method is directly applied to water and oxygen liquid-gas couples as countless engineering applications involve those two fluids. Safety studies of thermohydraulic systems of power plants and flows in cryotechnic rocket engines can be cited for instance. This latter example involves specific situations where transitions from pure fluid into two-phase mixture are present as well as transition to supercritical state. In the same context, combustion systems of modern automotive engines also involve transitions from pure phase to both two-phase mixture and supercritical fluid. Tables 1 and2 provide the associated parameters of the ENASG EOS (56), (57).

Table 1. Extended NASG (ENASG) coefficients for water. The NASG parameters are also given and determined with the method given in Le Métayer and Saurel (2016) [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF] except for the liquid reference entropy q ′ that is computed with the NASG reduction of Eq. (A31) (see Appendix C). The NASG water parameters are determined with N = 201 experimental saturation points in the temperature range T ex p ∈ [300 -500] K.

Coefficients

T (K) T (K) p sat (bar) L v (kJ/kg) h l (kJ/kg) h v (kJ/kg) ρ l (kg/m 3 ) ρ v (kg/m 3 ) Figure 1.
Comparison between experimental and theoretical saturation curves for liquid l and vapor v water. Symbols represent experimental data. The thick lines represent the theoretical saturation curves obtained with the Extended NASG EOS (ENASG) reducing to the ideal gas description for the vapor phase, Eqs. ( 56) and (57). The thin lines represent results obtained with the original NASG EOS also reducing to the ideal gas formulation for the vapor phase. p sat denotes the saturation pressure, L v the latent heat, h the specific enthalpy and ρ the density.

The ENASG EOS (56) presents good agreement with liquid experimental data at saturation. The saturated pressure resulting from the equality of the liquid and vapor Gibbs energies is rather good as well. Away from the critical point, the vapor phase, described by the ideal gas expressions (57), is also in good agreement with experimental data. However, the saturated vapor phase necessarily lacks of accuracy near the critical point as the attractive effects have been removed in order to keep an unambiguously convex formulation.

The results of the original NASG EOS (Le Métayer and Saurel (2016) [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF]) are plotted as well in Figs. 1 and2 for comparison. The corresponding parameters are given in Tables 1 and2. As the attractive pressure is constant in such formulation, liquid density necessarily lacks of accuracy away from its reference temperature range. However, as the ideal gas parameters have been determined thanks the saturation curve in [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF] (unlike the present work, see Appendix C), the vapor enthalpy and latent heat are in slightly better agreement than the present ENASG EOS (56), (57).

The present paper aims at building an overall EOS able to deal with pure liquid, pure vapor and supercritical phases, while being as accurate as possible at saturation. When thermodynamic conditions remain close to the saturation ones and away from the critical point, the original NASG EOS with its associated parameters (Le Métayer and Saurel (2016) [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF]) is preferred as the formulation is simpler than the ENASG one and yields excellent results as seen in Figs. 1 and2. However, as the original method [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF] uses the saturation curves for both liquid and vapor phases when determining the corresponding parameters, the NASG EOS lacks of accuracy away from the saturation conditions. This will be illustrated in the following section. As it will be seen later, the overall ENASG EOS presents good agreement with experimental data away from the saturated conditions while being rather satisfying at saturation (except for the vapor phase near the critical point as discussed earlier).

In the next section, the theoretical behavior of the ENASG EOS is analyzed with thermodynamic conditions corresponding to the transition from single phase to supercritical state.

Transition to supercritical fluids

This section deals with fluids transitioning from a pure phase to supercritical state as schematically illustrated in Fig. 3.

Liquid-vapor mixture

Supercritical state

Figure 3. The saturation curve is composed of the boiling and the dew curves separating the two-phase mixture zone and the pure phase zones. Beyond the critical isotherm, there is no transition between the liquid and the gaseous state. The fluid is neither liquid nor gas, it is said to be supercritical. Phase transition can happen either through the saturation dome corresponding to liquid-vapor phase change, or through the critical isotherm corresponding to a pure-phase-to-supercritical-state transition.

At high temperatures, the assumption of constant heat capacities no longer holds for the supercritical phase. The ideal gas description can still be fairly assumed in the supercritical state but variable heat capacities are required when reaching a certain temperature. Consequently, the following definition of heat capacity (at constant volume or pressure) is considered:

C v = C v,0 if T ≤ T 0 , C v = C v (T) otherwise, (58)
where C v,0 denotes the constant heat capacity given in Tables 1 and2. T 0 is the temperature at which the assumption of constant heat capacities starts to fail. At such temperature, the fluid is necessarily supercritical. These temperatures are reported in Table 3 for water and oxygen.

For ideal gases, Mayer's relation C p (T) -C v (T) = R holds and the ratio of heat capacities reads

γ(T) = C p (T)
C v (T) . In these relations, R = R/W with R denoting the universal gas constant and W the molar mass. In the present work, C p (T) is estimated via the NASA polynomial expression [START_REF] Mcbride | Coefficients for Calculating Thermodynamic and Transport Properties of Individual Species[END_REF],

C p (T) = R a 1 + a 2 T + a 3 T 2 + a 4 T 3 + a 5 T 4 , (59)
with corresponding parameters reported in Table 3.

Table 3. Parameters of the NASA polynomial expression [START_REF] Mcbride | Coefficients for Calculating Thermodynamic and Transport Properties of Individual Species[END_REF] for the heat capacity at constant pressure, Eq. ( 59). Proceeding similar derivations as in Section 2, the ideal gas formulation yields,

Fluid T 0 (K) a 1 a 2 (K -1 ) a 3 (K -2 ) a 4 (K -3 ) a 5 (K -
                                   p(v, T) = RT v , e(T) = C v (T)dT + q, h(T) = C p (T)dT + q, s(p, T) = C p (T) dT T -R ln(p) + q ′ , c 2 (T) = γ(T)RT = C p (T) C v (T) RT, C p (T) -C v (T) = R. (60)
As explained in Appendix C, the γ parameter of the gas phase is determined thanks to Mayer's relation as to ensure γ(T) -1 C v (T) = R. Note that the ideal gas reduction of the NASG EOS (Le Métayer and Saurel (2016) [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF]), with its associated original method to determine the different coefficients, does not ensure the preceding Mayer's relation as the gas parameters are determined with the saturation curves.

Liquid-to-supercritical-state transition

The liquid phase is described with the ENASG EOS as its particular interest resides in variable repulsive and attractive effects. In the transcritical zone, attractive and repulsive molecular forces are the dominant effects of the fluid. Thereby, the ENASG EOS is also used to describe liquids transitioning to their supercritical state. The transition is then continuous.

At much higher temperatures, thermal agitation becomes the dominant effect. The ideal gas description is then to be used in this context. However, as two different EOSs are used through different parameters, the continuity between the ENASG EOS and its ideal gas reduction is not trivial.

Indeed, the two EOSs must be connected in order to make a continuous formulation. For a given pressure, there exists a connection temperature where the two EOSs are equal. Nevertheless, those connection temperatures are a priori dependent on the variable of interest.

Equations ( 56) and (57) provide expressions of the different variables for the ENASG and ideal gas EOSs. Equality of both expressions provides the connection temperature that is the positive solution of a quadratic equation,

aT 2 + bT + c = 0. ( 61
)
Note that the admissible range of the sought-after temperature is known as this latter is necessarily higher than the critical one (T c ) and must ensure p + p ∞,l (T) > 0 and p + p ′ ∞,l (T) > 0. Note also that Eq. ( 61) is available only if the heat capacities are constant. An iterative method is required otherwise, but this situation (high temperatures, T ≥ T 0 ) is not to be encountered in practice as the fluid is necessarily supercritical and the ENASG EOS is not to be used (see section 7). The different parameters of the quadratic equation (61) are provided in Appendix D.

In the following, two isobars are considered for both fluids (water and oxygen). The first one is rather close (230 bars for water and 60 bars for oxygen) to the critical pressure (220 bars for water and 50 bars for oxygen) and the second is much higher (500 bars for water and 200 bars for oxygen). Figures 4,5, 6 and 7 show results corresponding to the transition from pure liquid to the supercritical state as schematically represented in Fig. 3.

As seen in Figs. 4, 5, 6 and 7, the liquid ENASG EOS (56) presents excellent agreement with experimental data. Passed the critical temperature T c , the ENASG EOS is also used until connection with the ideal gas EOS and yields good agreement as well.

At the temperature of connection, the ideal gas EOS is considered with constant heat capacities until the temperature T 0 is reached. From this temperature, variable thermal capacities are used. Consequently, the ideal gas formulation is rather well-suited in the supercritical state. At such high temperatures, thermal agitation is indeed expected to be the dominant effect determining the properties of the fluid. The fundamental assumption of the ideal gas description (molecules free to evolve regardless of the others) then reappears in such thermodynamic conditions. The results computed with the ideal gas expression in the supercritical state are in excellent agreement with experimental data with the exception of the specific volume of supercritical water that presents a lesser agreement.

The results of the overall ENASG formulation are in good agreement with experimental data at both pressures close to the critical one and much higher. They also show the good behavior of the ENASG EOS when dealing with conditions away from the saturation ones.

As seen in Figs. 4, 5, 6 and 7, the extension of the liquid ENASG EOS results in good agreement with experimental data and provides a continuous formulation in the transcritical zone.

The only discontinuous thermodynamic variable is the speed of sound. This is clearly seen in Fig. 7 for example. Regarding the sound speed, the liquid ENASG EOS is not extended beyond the critical temperature as the formulation may not connect to the ideal gas expression. The discontinuous speed of sound at the critical temperature appears to be in practice similar to situations involving large sound-speed variations such as the transition from a pure fluid into a two-phase mixture at equilibrium, which does not cause practical difficulties. Obviously, this flaw is not encountered with cubic EOSs.

The ENASG EOS is also compared to cubic ones in the preceding figures. The van der Waals (VdW) [START_REF] Van Der Waals | On the Continuity of the Gaseous and Liquid States[END_REF] and Soave-Redlich-Kwong (SRK) [START_REF] Soave | Equilibrium constants from a modified Redlich-Kwong equation of state[END_REF] are used in this work. Detailed reviews of cubic equations of state can be found in Wei and Sadus (2000) [START_REF] Wei | Equations of state for the calculation of fluid-phase equilibria[END_REF] for instance.

As shown in Figs. 4 and5, the cubic EOSs present poor accuracy regarding liquid water. However, the supercritical phase is well-described and the transition from liquid to supercritical state is naturally continuous since a unique formulation is used for both phases.

When oxygen is considered (Figs. 6 and7), the VdW EOS shows very good results regarding the supercritical phase but is unable to represent properly the liquid state. Nevertheless, the SRK EOS presents excellent agreement with experimental data for both liquid and supercritical phases.

Vapor-to-supercritical-state transition

Much lower pressures (30 bars for water and 10 bars for oxygen) are considered in Figs. 8 and9. Consequently, the three states of the corresponding fluids are involved (liquid, vapor, supercritical) and the vapor phase transforms to supercritical fluid beyond the critical temperature (Fig. 3). 56), ( 57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation for the supercritical phase. The dashed lines represent the van der Waals (VdW) theoretical predictions and the dashed-dotted lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. From this temperature, the liquid ENASG EOS is extended and joins the ideal gas EOS (except for the sound speed). The temperature T 0 = 1000 K at which variable heat capacities are considered is indicated in dotted lines as well. 56), ( 57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation for the supercritical phase. The dashed lines represent the van der Waals (VdW) theoretical predictions and the dashed-dotted lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. From this temperature, the liquid ENASG EOS is extended and joins the ideal gas EOS (except for the sound speed). The temperature T 0 = 1000 K at which variable heat capacities are considered is indicated in dotted lines as well. 56), ( 57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation for the supercritical phase. The dashed lines represent the van der Waals (VdW) theoretical predictions and the dashed-dotted lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. From this temperature, the liquid ENASG EOS is extended and joins the ideal gas EOS (except for the sound speed). The temperature T 0 = 400 K at which variable heat capacities are considered is indicated in dotted lines as well. 56), ( 57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation for the supercritical phase. The dashed lines represent the van der Waals (VdW) theoretical predictions and the dashed-dotted lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. From this temperature, the liquid ENASG EOS is extended and joins the ideal gas EOS (except for the sound speed). The temperature T 0 = 400 K at which variable heat capacities are considered is indicated in dotted lines as well. As seen in Fig. 8, the ENASG EOS is able to represent correctly the liquid water unlike the VdW and SRK ones. The vapor and supercritical phases are rather well-described with all EOSs (ideal gas, VdW, SRK) and with continuous formulations (except for the sound speed with the ideal gas EOS). 57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation. The dashed-dotted lines represent the van der Waals (VdW) theoretical predictions and the dashed lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. The temperature T 0 = 400 K at which variable heat capacities are considered is indicated in dotted lines as well.

The ENASG EOS is also well-suited for oxygen as seen in Fig. 9. The VdW EOS is again unable to represent properly the liquid state but the SRK one presents excellent results.

Concluding remarks

The results of the present section illustrate the good behavior of the ENASG EOS (56) and its reduction to the ideal gas expression (57), (60) in situations away from the saturation thermodynamic conditions. Figures 4,5, 6 and 7 show that the proposed EOS is able to deal with pure liquids and supercritical states in pressure conditions both close and much higher than the critical pressure. Figures 8 and9 illustrate the good behavior of the overall formulation with much lower pressures as well. The liquid expression seems able to represent the whole liquid phase diagram including the saturation conditions and the transcritical zone with a unique set of parameters. The corresponding vapor and supercritical phases are described accurately as well with the ideal gas expressions, except near the critical point as attractive terms have been removed in order to remain unambiguously convex, as discussed earlier.

The results of the original NASG EOS (Le Métayer and Saurel (2016) [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF]), with associated parameters given in Tables 1 and2, are also plotted in Figs. 4, 5, 6, 7, 8 and 9 for comparison. As already discussed, the supercritical phase is inaccurate as the different parameters have been determined with the help of the saturation curve.

Nevertheless, for flows evolving away from saturation and at pressures much lower than the critical one, the NASG EOS is preferred to the ENASG one for the sake of simplicity. However, its reduction to the ideal gas formulation should use parameters determined away from the saturation curves since the ones determined at the thermodynamic equilibrium [START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF] induce inaccuracy as seen in Figs. 4, 5, 6, 7, 8 and 9. The results at saturation may be slightly degraded but the overall formulation is expected to deal with multiple thermodynamic conditions.

Two-phase flow illustrations

In the following, two-phase flows subject to phase change are of interest. When evaporation or condensation phenomena appear, instantaneous phase transition is considered through the stiff thermochemical relaxation solver of Chiapolino et al. (2017) [START_REF] Chiapolino | A simple phase transition relaxation solver for liquid-vapor flows[END_REF], [START_REF] Chiapolino | A simple and fast phase transition relaxation solver for compressible multicomponent two-phase flows[END_REF]. For the sake of simplicity, the Homogeneous Relaxation Model (HRM) [START_REF] Downar-Zapolski | The non-equilibrium relaxation model for one-dimensional flashing liquid flow[END_REF] is considered and is reminiscent of the reactive (or multicomponent) Euler equations widely used in chemically reacting flows. The corresponding system reads,

                                       ∂ρ ∂t + div(ρu) = 0, ∂(ρu) ∂t + div ρu ⊗ u + pI = 0, ∂(ρE) ∂t + div (ρE + p)u = 0, ∂(ρY 1 ) ∂t + div(ρY 1 u) = ρν(g 2 -g 1 ), ∂(ρY 2 ) ∂t + div(ρY 2 u) = -ρν(g 2 -g 1 ), ∂(ρY k ) ∂t + div(ρY k u) = 0, (62)
with

E = e + 1 2 u 2 , e = N ∑ k=1 Y k e k .
System (62) considers implicitly mechanical and thermal equilibrium. The thermodynamic equilibrium is reached through the instantaneous relaxation (ν → ∞) of Gibbs free energies g 1 = g 2 where the indexes 1 and 2 denotes respectively the liquid and vapor phases (see [START_REF] Chiapolino | A simple phase transition relaxation solver for liquid-vapor flows[END_REF], [START_REF] Chiapolino | A simple and fast phase transition relaxation solver for compressible multicomponent two-phase flows[END_REF]). The other constituents of the flow (N = 3 → N) are considered as non-condensable gases. u represents the mixture centre of mass velocity and E the mixture total energy. System (62) is closed by a mixture equation of state made from the mechanical and thermal equilibrium. In the first place, let us considered gaseous flows transitioning to a supercritical state. When the critical temperature is reached, liquid is no longer present and the ENASG EOS is not to be used. Following the strategy of Chiapolino et al. (2017) [START_REF] Chiapolino | A simple phase transition relaxation solver for liquid-vapor flows[END_REF], [START_REF] Chiapolino | A simple and fast phase transition relaxation solver for compressible multicomponent two-phase flows[END_REF], two expressions of the mixture temperature can be found according to the definitions of the mixture mass and mixture energy,

T = T k ∀k, p = p k ∀k, v = ∑ (Y k v k ), e = ∑ (Y k e k ), (63)
with Y k denoting the mass fraction of the chemical species k.

In practical computations, gaseous mixture can be considered if Y 1 < ǫ with ǫ ≃ 10 -8 . In that case, k = 2 → N and the combination of Eqs. ( 63) and (56) leads to,

T v = pv ∑ N k=2 Y k (γ k -1)C v,k
and

T e = e -∑ N k=2 Y k q k ∑ N k=2 Y k C v,k . (64)
Equality of T v and T e provides the mixture gas pressure,

p(e, v) = e -∑ N k=2 Y k q k ∑ N k=2 Y k (γ k -1)C v,k v ∑ N k=2 Y k C v,k . ( 65
)
Figure 10 shows such a situation where water vapor transforms into supercritical state through compression effects of a shock wave. A shock tube is indeed considered with liquid water, vapor water and air. In the high pressure chamber, air is initially in major proportion, Y 3 → 1 with thermodynamic conditions p = 30 bars and T = 800 K. In the second chamber, water vapor is in major proportion Y 2 → 1 with p = 1 bar and T = 600 K. The mixture is initially at thermodynamic equilibrium according to the conditions detailed in [START_REF] Chiapolino | A simple and fast phase transition relaxation solver for compressible multicomponent two-phase flows[END_REF]. The ideal gas reduction of the ENASG EOS is used with parameters for water given in Table 1 while the only coefficients needed for air are: C v,3 = 719 J/kg/K and γ 3 = 1.4. Liquid is present but in negligible proportions so mixture equation of state (65) is used in practice.

As seen in Fig. 10, the transition from "pure" vapor to supercritical state is naturally continuous when the temperature becomes higher than the critical one T c = 646 K.

Let us now consider a situation where liquid is in major proportion. In that case, the combination of the mixture definitions (63) and the ENASG relations (56) leads to two quadratic expressions for the mixture temperature (note that only one liquid is considered in this work),

aT 2 + bT + c = 0, (66)
with the corresponding coefficients,

                     a v = γ 1 p ∞,1,1 p N ∑ k=2 Y k (γ k -1)C v,k , b v = Y 1 (γ 1 -1)C v,1 1 -b 1,1 + N ∑ k=2 Y k (γ k -1)C v,k 1 + γ 1 p ∞,0,1 (1 -b 1,1 ) (γ 1 -b 1,1 )p + b0 -v γ 1 p ∞,1,1 , c v = b0 -v p + γ 1 p ∞,0,1 (1 -b 1,1 ) γ 1 -b 1,1 , (67)
                     a e = Y 1 γ 1 p ∞,1,1 C v,1 + N ∑ k=2 Y k C v,k γ 1 p ∞,1,1 , b e = Y 1 (p + γ 1 p ∞,0,1 ) C v,1 + N ∑ k=2 Y k C v,k p + γ 1 p ∞,0,1 (1 -b 1,1 ) γ 1 -b 1,1 + ( q -e) γ 1 p ∞,1,1 , c e = ( q -e) p + γ 1 p ∞,0,1 (1 -b 1,1 ) γ 1 -b 1,1 , (68)
where mixture quantities are introduced: ≃ 0.99999989. The test was carried out with Godunov time integration method and HLLC Riemann solver extended to second order: MUSCL scheme with Minmod flux limiter. The solution is given at t ≈ 0.3 ms on a 1000-cell mesh using CFL=0.8.

q = N ∑ k=1 Y k q k , b = N ∑ k=1 Y k b k . ( 69
)
Equality of the two positive solutions provides the mixture pressure p(e, v). An iterative method is required nonetheless. However, Eq. ( 65) is to be used where Y 1 ≤ ǫ → 0 corresponds to a gaseous mixture.

The transition from supercritical state to pure liquid is now considered through a double expansion test. In Fig. 11, vapor water and air are present in negligible proportions and supercritical water undergoes expansion waves. Those induce pressure drop from 350 bars to about 226 bars. The final pressure then remains slightly above the critical one, p c = 220 bars. They also induce temperature drop from 655 K to about 641 K. The final temperature is consequently inferior to the critical one T c = 646 K resulting in transition from supercritical to liquid phase that is computed continuously. A configuration where liquid-gas interfaces are present is now considered. Phase change is illustrated on the evaporating liquid jet configuration detailed in Chiapolino et al. (2017) [START_REF] Chiapolino | A simple and fast phase transition relaxation solver for compressible multicomponent two-phase flows[END_REF]. In this last reference, the mixture EOS is made from the NASG EOS for each fluid (reduced to SG in this test) and is reconsidered in the following in the frame of the ENASG EOS.

The conditions are typical of cryotechnic rocket engines during the ignition phase (for which the engine has not yet reached supercritical conditions). The flow consists of a coaxial liquid oxygen jet surrounded by a high-speed hydrogen flow (non-condensable gas), injected in conditions above the saturation point of the inner oxygen core, which then evaporates whilst being destabilized. Such a case is very challenging, because there is initially no vapor oxygen, and mass transfer is the only possible term for vapor production. The ENASG EOS is used with parameters for oxygen given in Table 2 while the only coefficients needed for hydrogen are: C v,3 = 10183 J/kg/K, γ 3 = 1.4 and q 3 = -1.2 • 10 6 J/kg. Mass transfer is treated with the thermochemical relaxation solver detailed in [START_REF] Chiapolino | A simple and fast phase transition relaxation solver for compressible multicomponent two-phase flows[END_REF]. Figure 12 As expected, the filaments separating the main liquid core and the gas gradually vanish as a consequence of evaporation and the created vapor mass fraction is of utmost importance for future works, which shall include the gaseous combustion between vapor oxygen and hydrogen. In that case, only the ideal gas reduction of the ENASG EOS is to be used with variable heat capacities.

Conclusion

The Noble-Abel-Stiffened-Gas (NASG) equation of state has been extended to variable attractive and repulsive effects to deal with liquids when large temperature and pressure ranges are under consideration. The liquid phase is well-described at thermodynamic conditions both near and away from the saturation ones with a convex formulation. The overall ENASG EOS reduces to the ideal gas description for both vapor and supercritical phases for the sake of convexity.

The transition from pure fluid to supercritical state is of interest as well, including at high pressures where the liquid directly transforms to supercritical fluid. The ENASG EOS proposes a solution in the direction of such transition while remaining convex, an essential property in computational fluid dynamics.

Two different liquid-gas couples have been addressed, water and oxygen, presenting respectively triatomic and diatomic molecular fluids. The overall formulation presents good agreement with experimental data. However, the saturated vapor phase necessarily lacks of accuracy near the critical point as attractive effects are absent.

Those latter ones are nonetheless responsible for conditional convexity, a feature reminiscent of cubic equations of state.

The Extended NASG (ENASG) equation of state recovers the NASG one when the new introduced coefficients are set to zero. Its formulation remains quite simple, convex and is beneficial to the introduction of phase transition solvers such as the ones promoted in Chiapolino et al. (2017) [START_REF] Chiapolino | A simple phase transition relaxation solver for liquid-vapor flows[END_REF], [START_REF] Chiapolino | A simple and fast phase transition relaxation solver for compressible multicomponent two-phase flows[END_REF].

Convexity of the equation of state requires fulfillment of five different conditions [START_REF] Godunov | Résolution numérique des problèmes multidimensionnels de la dynamique des gaz[END_REF] [START_REF] Godunov | Résolution numérique des problèmes multidimensionnels de la dynamique des gaz[END_REF], [START_REF] Menikoff | The Riemann problem for fluid flow of real materials[END_REF] [START_REF] Menikoff | The Riemann problem for fluid flow of real materials[END_REF] that are analyzed hereafter,

a) ∂ 2 e ∂v 2 s > 0, b) ∂ 2 e ∂s 2 v > 0, c) ∂ ∂s ∂e ∂v s v < 0, d) ∂ 2 e ∂s 2 v ∂ 2 e ∂v 2 s - ∂ ∂s ∂e ∂v s v 2 > 0, e) ∂ 3 e ∂v 3 s < 0. (A1)
Combining Eqs. ( 20), ( 21) and (40), the internal energy expresses as,

e(v, s) =q + (γ -1)C v γp ∞,0 v -b(v) (γ -b 1 ) (γ -1)C v -γp ∞,1 v -b(v) - γ 2 p ∞,1 p ∞,0 v -b(v) 2 (γ -b 1 ) (γ -1)C v -γp ∞,1 v -b(v) + C v exp s-q ′′ C v exp γp ∞,1 v-b(v) C v (1-b 1 ) v -b(v) γ-b 1 1-b 1 .
(A2) After some algebraic manipulations, the first partial derivative reads,

∂e ∂v s = exp s-q ′′ C v exp γp ∞,1 v-b(v) C v (1-b 1 ) -(γ -1)C v + γp ∞,1 v -b(v) v -b(v) γ-b 1 1-b 1 + γp ∞,0 (1 -b 1 ) γ -b 1 .
(A3) Analyzing Eqs. (40) and (A3), it appears that the thermodynamic definition of the pressure is satisfied, p = -∂e ∂v s . Continuing the calculations, the second derivative reads,

∂ 2 e ∂v 2 s exp s-q ′′ C v =(γ -b 1 )(γ -1)C v v -b(v) -γ-1+2b 1 1-b 1 exp γp ∞,1 v -b(v) C v (1 -b 1 ) + (γp ∞,1 ) 2 C v exp γp ∞,1 v -b(v) C v (1 -b 1 ) v -b(v) 1-γ 1-b 1 -2γp ∞,1 (γ -1) v -b(v) -γ+b 1 1-b 1 exp γp ∞,1 v -b(v) C v (1 -b 1 ) . (A4)
Analyzing Eq. (A4), it appears that condition (A1) (a) is satisfied unambiguously if p ∞,1 ≤ 0 and b 1 < γ. Equation (A2) is now used and yields the following partial derivative,

∂e ∂s v = exp s-q ′′ C v exp γp ∞,1 v-b(v) C v (1-b 1 ) v -b(v) γ-1 1-b 1 . (A5)
Furthermore, manipulating Eq. ( 33), the liquid temperature can be expressed as, . With the help of (A5), the second partial derivative is expressed as,

T(v, s) = exp s-q ′′ C v exp γp ∞,1 v-b(v) C v (1-b 1 ) v -b(v) γ-1 1-b 1 . ( A6
∂ 2 e ∂s 2 v = exp s-q ′′ C v exp γp ∞,1 v-b(v) C v (1-b 1 ) C v v -b(v) γ-1 1-b 1 . ( A7
∂v s v = exp s-q ′′ C v exp γp ∞,1 v-b(v) C v (1-b 1 ) γp ∞,1 v -b(v) -(γ -1)C v C v v -b(v) γ-b 1 1-b 1 . (A8)
Condition (A1) (c) is then unambiguously satisfied and defined if p ∞,1 ≤ 0 and b 1 = 1. Besides, combining Eqs. (A4), (A7) and (A8) leads to the next relation,

∂ 2 e ∂s 2 v ∂ 2 e ∂v 2 s - ∂ ∂s ∂e ∂v s v 2 = exp s -q ′′ C v 2 exp γp ∞,1 v -b(v) C v (1 -b 1 ) 2 (1 -b 1 )(γ -1) v -b(v) -2γ+2b 1 1-b 1 . (A9)
Condition (A1) (d) is then satisfied as well if b 1 < 1. Finally, from Eq. (A4), relation (A1) (e) reads,

∂ 3 e ∂v 3 s exp s-q ′′ C v = γp ∞,1 C v exp γp ∞,1 v -b(v) C v (1 -b 1 ) (γ -b 1 )(γ -1)C v v -b(v) -γ-1+2b 1 1-b 1 + (γp ∞,1 ) 2 C v v -b(v) 1-γ 1-b 1 -2γp ∞,1 (γ -1) v -b(v) -γ+b 1 1-b 1 + exp γp ∞,1 v -b(v) C v (1 -b 1 ) (γ -b 1 )(γ -1)C v (-γ -1 + 2b 1 ) v -b(v) -γ-2+3b 1 1-b 1 -(γ -1) (γp ∞,1 ) 2 C v v -b(v) -γ+b 1 1-b 1 + 2γp ∞,1 (γ -1)(γ -b 1 ) v -b(v) -γ-1+2b 1 1-b 1 . (A10) Analyzing Eq. (A10), condition (A1) (e) is satisfied unambiguously if p ∞,1 ≤ 0, b 1 < γ and b 1 < 1 2 + γ 2 .
As γ > 1, the most restrictive condition regarding the covolume remains b 1 < 1.

The present formulation is then unambiguously convex if,

p ∞,1 ≤ 0, p ∞,0 ≥ 0 and b 1 < 1. (A11)

Appendix B Maxwell's relations

Maxwell's relations arise from the equality of the mixed partial derivatives of the fundamental thermodynamic relations [START_REF] Callen | an introduction to the physical theories of equilibrium thermostatics and irreversible thermodynamics[END_REF]. The different functions of common use read, ) This section aims at verifying that those fundamental relations are satisfied with the ENASG formulation. Using Eqs. ( 18) and (38), the next relation directly arises,

a) ∂s ∂p T = - ∂v ∂T p , b) ∂p ∂T v = ∂s ∂v T , c) ∂T ∂v s = - ∂p ∂s v , d) ∂T ∂p s = ∂v ∂s p . ( A12
∂s ∂p T = - (γ -1)C v (1 -b 1 ) p + p ′ ∞ (T) + γp ∞,1 (γ -1)C v T (1 -b 1 ) p + p ′ ∞ (T) 2 = - ∂v ∂T p . ( A13
)
It is then clear that Maxwell's relation (A12) (a) is satisfied. Equations ( 18) and ( 33) are now used and lead to,

∂p ∂T v = (γ -1)C v v -b(v) -γp ∞,1 = ∂s ∂v T . ( A14
)
Consequently, Maxwell's relation (A12) (b) is also satisfied. This is not surprising as these two equations, (A12) (a)-(b), are precisely Eqs. ( 23) and ( 4), which are the basis of the theoretical derivations. Equations (A6) and (40) being now considered, the following partial derivatives are obtained, Let us then analyze the fourth relation. Considering Eq. (40), it is noted that v(s, p) can not be directly formulated unless p ∞,1 = 0, that reduces the formulation to the NASG equation of state with variable covolume. However, the partial derivatives can be directly formulated. Indeed, as the left-hand side of Eq. (A12) (d)

∂T ∂v s = exp s -q ′′ C v exp γp ∞,1 v -b(v) C v (1 -b 1 ) -(γ -1) v -b(v) -γ+b 1 1-b 1 + γp ∞,1 C v v -b(v) -γ+1 1
= - -(γ -1)C v p + p ′ ∞ (T) + γp ∞,1 (γ -1)C v T C v T (γ -b 1 ) p + p ′ ∞ (T) 2 -γp ∞,1 (γ -1)C v 2p + γp ∞,1 T + 2γp ∞,0 (1-b 1 ) γ-b 1 . ( A18
)
∂p ∂s v has also been determined previously, Eq. (A15) and thanks to relation (40), the next derivative arises after some algebraic manipulations,

∂p ∂v s = exp s -q ′′ C v exp γp ∞,1 v -b(v) C v (1 -b 1 ) -(γ -b 1 )(γ -1)C v v -b(v) -γ-1+2b 1 1-b 1 + 2(γ -1)γp ∞,1 v -b(v) -γ+b 1 1-b 1 - (γp ∞,1 ) 2 C v v -b(v) 1-γ 1-b 1 .
(A19) Equation (A17) then reads after calculations,

∂v ∂s p = - (γ -1)C v -γp ∞,1 v -b(v) -C 2 v (γ -1)(γ -b 1 ) v -b(v) -1 + 2γp ∞,1 (γ -1)C v -(γp ∞,1 ) 2 v -b(v) . ( A20
)
Inserting Eq. ( 18) into Eq. (A20), the following result is obtained after some algebraic manipulations, + γ l p ∞,1,l (γ l -1)C v,l T ex p,l,i 1b 1,l p ex p,l,i + p ′ ∞ (T ex p,l,i )

∂v ∂s p = - -(γ -1)C v p + p ′ ∞ (T) + γp ∞,1 (γ -1)C v T C v T (γ -b 1 ) p + p ′ ∞ (T) 2 -γp ∞,1 (γ -1)C v 2p + γp ∞,1 T + 2γp ∞,0 (
.

(A31) The different reference state values used for the calculation of the liquid ENASG coefficients are summarized in Tables A1 andA2. In this work, all experimental data come from NIST1 website. For the liquid phase, the saturation (boiling) curve is considered.

Gas phase

In the present formulation (57), the gas phase is considered as ideal and the different parameters are determined regardless of the saturation conditions. Only four parameters are required for the gas phase, C v , γ, q and q ′ . The atmospheric conditions are used in this work via the experimental isobar p 0 = 1 bar. According to the experimental data of water and oxygen at such low pressure, there exists a significant temperature range where the heat capacity (C v ) is quite constant. The parameter C v is thereby chosen as a constant, representative of the present thermodynamic conditions. The parameter γ is then determined as γ = C p C v with C p -C v = R/W according to Mayer's relation. In the previous relation, R denotes the universal gas constant and W the molar mass. The values reported in Tables 1 and2 are consequently close to the expected triatomic (H 2 O) and diatomic (O 2 ) predictions for ideal gases (γ = 1.3079 ≃ 9/7 and C v = 1500 J/kg/K ≃ (7/2)R for water and γ = 1.3985 ≃ 1.4 and C v = 652 J/kg/K ≃ (5/2)R for oxygen).

A reference point on the present isobar p 0 = 1 bar is used to determine the coefficient q, q = e 0 -C v T 0 . (A32)

In this work T 0 = 393.38 K and e 0 = 2537.7 kJ/kg are used for water and T 0 = 100.07 K and e 0 = 63.657 kJ/kg for oxygen. The last coefficient q ′ g is finally determined with the least squares method, corresponding to the ideal gas reduction of Eq. (A31). The isobar p 0 = 1 bar is used one more time with N = 542 points for water corresponding to the temperature range T ex p ∈ [372.76 -1275] K and N ex p = 579 points, T ex p ∈ [90.062 -1000] K for oxygen.

where h represents the

  specific enthalpy. As the thermodynamic definition of entropy and specific volume implies s = -∂g ∂T p and v = ∂g ∂p T

Figure 2 .

 2 Figure 2. Comparison between experimental and theoretical saturation curves for liquid l and vapor v oxygen. Symbols represent experimental data. The thick lines represent the theoretical saturation curves obtained with the Extended NASG EOS (ENASG) reducing to the ideal gas description for the vapor phase, Eqs. (56) and (57). The thin lines represent results obtained with the original NASG EOS also reducing to the ideal gas formulation for the vapor phase. p sat denotes the saturation pressure, L v the latent heat, h the specific enthalpy and ρ the density.

  -supercritical-state transition Vapor-to-supercritical-state transition

h

  (kJ/kg) e (kJ/kg)ρ (kg/m 3 ) v (m 3 /kg) s (J/kg/K) c (m/s)

Figure 4 .

 4 Figure 4. Comparison between experimental and theoretical isobar for water. The symbols represent the experimental isobar p = 230 bars. Beyond the critical temperature T c = 646 K, the liquid transforms to supercritical state. The thick solid lines represent the Extended NASG (ENASG) EOS, reducing to the ideal gas description for the supercritical phase, Eqs. (56), (57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation for the supercritical phase. The dashed lines represent the van der Waals (VdW) theoretical predictions and the dashed-dotted lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. From this temperature, the liquid ENASG EOS is extended and joins the ideal gas EOS (except for the sound speed). The temperature T 0 = 1000 K at which variable heat capacities are considered is indicated in dotted lines as well.

hFigure 5 .

 5 Figure 5. Comparison between experimental and theoretical isobar for water. The symbols represent the experimental isobar p = 500 bars. Beyond the critical temperature T c = 646 K, the liquid transforms to supercritical state. The thick solid lines represent the Extended NASG (ENASG) EOS, reducing to the ideal gas description for the supercritical phase, Eqs. (56), (57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation for the supercritical phase. The dashed lines represent the van der Waals (VdW) theoretical predictions and the dashed-dotted lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. From this temperature, the liquid ENASG EOS is extended and joins the ideal gas EOS (except for the sound speed). The temperature T 0 = 1000 K at which variable heat capacities are considered is indicated in dotted lines as well.

hFigure 6 .

 6 Figure 6. Comparison between experimental and theoretical isobar for oxygen. The symbols represent the experimental isobar p = 60 bars. Beyond the critical temperature T c = 154 K, the liquid transforms to supercritical state. The thick solid lines represent the Extended NASG (ENASG) EOS, reducing to the ideal gas description for the supercritical phase, Eqs. (56), (57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation for the supercritical phase. The dashed lines represent the van der Waals (VdW) theoretical predictions and the dashed-dotted lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. From this temperature, the liquid ENASG EOS is extended and joins the ideal gas EOS (except for the sound speed). The temperature T 0 = 400 K at which variable heat capacities are considered is indicated in dotted lines as well.

hFigure 7 .

 7 Figure 7. Comparison between experimental and theoretical isobar for oxygen. The symbols represent the experimental isobar p = 200 bars. Beyond the critical temperature T c = 154 K, the liquid transforms to supercritical state. The thick solid lines represent the Extended NASG (ENASG) EOS, reducing to the ideal gas description for the supercritical phase Eqs. (56), (57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation for the supercritical phase. The dashed lines represent the van der Waals (VdW) theoretical predictions and the dashed-dotted lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. From this temperature, the liquid ENASG EOS is extended and joins the ideal gas EOS (except for the sound speed). The temperature T 0 = 400 K at which variable heat capacities are considered is indicated in dotted lines as well.

hFigure 8 .

 8 Figure 8. Comparison between experimental and theoretical isobar for water. The symbols represent the experimental isobar p = 30 bars. Beyond the critical temperature T c = 646 K, the vapor transforms to supercritical state. The thick solid lines represent the Extended NASG (ENASG) EOS, reducing to the ideal gas description for vapor and supercritical phases, Eqs. (56), (57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation. The dashed-dotted lines represent the van der Waals (VdW) theoretical predictions and the dashed lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. The temperature T 0 = 1000 K at which variable heat capacities are considered is indicated in dotted lines as well.

hFigure 9 .

 9 Figure 9. Comparison between experimental and theoretical isobar for oxygen. The symbols represent the experimental isobar p = 10 bars. Beyond the critical temperature T c = 154 K, the vapor transforms to supercritical state. The thick solid lines represent the Extended NASG (ENASG) EOS, reducing to the ideal gas description for vapor and supercritical phases, Eqs. (56), (57) and (60). The thin solid lines represent the original NASG EOS also reducing to the ideal gas formulation. The dashed-dotted lines represent the van der Waals (VdW) theoretical predictions and the dashed lines represent the Soave-Redlich-Kwong (SRK) ones. The critical temperature is indicated with the dotted lines. The temperature T 0 = 400 K at which variable heat capacities are considered is indicated in dotted lines as well.

Figure 10 .f t 3 = 1 - 2 •f t 1 ≃ 10 -8 and Y le f t 2 ≃ 1 . 9 • 3 = 1 ≃

 1031211021931 Figure10. Shock tube test illustrating the transition from "pure" water vapor to supercritical state. The critical temperature is indicated with the dotted line. The thick lines represent the solution obtained with the mixture ENASG EOS reducing to Eq. (65) in the present example as liquid mass fraction is non-zero but in negligible proportions. The dashed lines represent the initial conditions. In the left chamber, air is initially in major proportions with Y le f t 3 = 1 -2 • 10 -7 , p = 30 bars and

Figure 11 .

 11 Figure 11. Double expansion test illustrating the transition from supercritical state to "pure" liquid water. The critical pressure and temperature are indicated with the dotted lines. The thick lines represent the solution obtained with the mixture ENASG EOS. The dashed lines represent the initial conditions. Liquid water is initially in major proportions withY 1 = 1 -2 • 10 -6 , Y 2 = Y 3 = 10 -6, p = 350 bars, T = 655 K and u = ±45 m/s. The test was carried out with Godunov time integration method and HLLC Riemann solver extended to second order: MUSCL scheme with Minmod flux limiter. The solution is given at t ≈ 0.3 ms on a 1000-cell mesh using CFL=0.8.

3 )

 3 shows the density contours and the vapor mass fraction created. Density (kg.m -Vapor mass fraction

Figure 12 .

 12 Figure12. Density and vapor mass fraction profiles of a liquid oxygen jet surrounded by hydrogen at high speed entering a combustion chamber of a cryotechnic rocket engine. Shear effects induce jet fragmentation. The filaments separating the main liquid core and the gas gradually vanish as a consequence of evaporation. The computation was done with the MUSCL scheme with Superbee limiter and CFL=0.7. The solution is given at t ≈ 4.1 ms. The mesh is unstructured and made of about 360, 000 triangles.

  )Analyzing Eqs. (A5) and (A6), it appears that the thermodynamic definition of the temperature is satisfied, T = ∂e ∂s v

  ) Condition (A1) (b) is then satisfied ∀ p ∞,1 and b 1 = 1. Also, from Eq. (A3), relation (A1) (c) transforms to,

  showing that Maxwell's relation (A12) (c) is satisfied as well.

1 -b 1 ) γ-b 1 .

 111 Eqs. (A18) and (A21), the last Maxwell's relation (A12) (d) is satisfied. ,l,i + p ′ ∞ (T ex p,l,i )

Table 2 .

 2 Extended NASG (ENASG) coefficients for oxygen. The NASG parameters are also given and determined with the method given in Le Métayer and Saurel (2016)[START_REF] Métayer | The Noble-Abel Stiffened-Gas equation of state[END_REF] except for the liquid reference entropy q ′ that is computed with the NASG reduction of Eq. (A31) (see Appendix C). The NASG oxygen parameters are determined with N = 41 experimental saturation points in the temperature range T ex p ∈ [60 -100] K.

		ENASG Liq	ENASG gas	NASG Liq	NASG gas
	γ	1.0147	1.3079	1.1807	1.5377
	C v (J/kg/K) b 1 b 0 (m 3 /kg) p ∞,1 (Pa/K)	4014 -0.6050 1.5196 × 10 -3 -471025	1500 0 0 0	3630 0 6.8428 × 10 -4 0	856 0 0 0
	p ∞,0 (Pa)	307078403	0	664961465	0
	q (J/kg) q ′ (J/kg/K)	-1112426 -22049	1947630 1136	-1178154 -10742	2176064 4863
	Coefficients	ENASG Liq	ENASG gas	NASG Liq	NASG gas
	γ	1.0281	1.3985	1.6610	1.4730
	C v (J/kg/K) b 1 b 0 (m 3 /kg) p ∞,1 (Pa/K)	1535 -0.6721 1.3131 × 10 -3 -324997	652 0 0 0	1016 0 5.7003 × 10 -4 0	548 0 0 0
	p ∞,0 (Pa)	50890107	0	196815802	0
	q (J/kg) q ′ (J/kg/K)	-278134 -3691	-1589 4237	-285545 8171	6528 4650

  considers constant entropy, the following relation can be used,

		ds =	∂s ∂p T	dp +	∂s ∂T p	dT = 0.
	Consequently, the partial derivative is found as,
		dT dp	=	∂T ∂p s	= -	∂s ∂p T ∂T p ∂s	.	(A16)
	The same reasoning is repeated for ∂v ∂s p	yielding,
		dv ds	=	∂v ∂s p	= -	∂p ∂v s ∂s v ∂p	.	(A17)
	∂s ∂p T after calculations, and ∂s ∂T p	have been determined previously, Eqs. (A13) and (35). Equation (A16) then reads
	∂T
	∂p s

Table A1 .

 A1 Reference state values used for the determination of liquid ENASG coefficients. Fluid N T c (K) p c (bar) v c (m 3 /kg) p ′ ∞,c (Pa) b c (m 3 /kg) c 0 (m/s) p 0 (bar) v 0 (m 3 /kg)

	H 2 O O 2	374 646.16 101 154.36	221 50	0.0025101 0.0019522	0.01 0.01	10 -6 10 -6	1552.1 1065.7	1 1	0.0010182 0.00080871

Table A2 .

 A2 Reference state values used for the determination of liquid ENASG coefficients (continued).FluidT re f (K) p re f (Pa) v re f (m 3 /kg) e re f (kJ/kg) b re f (m 3 /kg)

	H 2 O O 2	300.16 70.631	3570.2 6684.7	0.0010035 0.00080952	113.23 -166.823	0.0009125 0.000769

Table A3 .

 A3 Coefficients for water and oxygen for the "alternative" ENASG EOS whose formulation is summarized in Eq. (A38). With such description, the gas attractive effects are taken into account via the parameter d but result in conditional convexity, Eq. (A39). The liquid ENASG EOS is unchanged, Eq. (56).CoefficientsENASG H 2 O, Liq ENASG H 2 O, vap ENASG O 2 , Liq ENASG O 2 , vap

	γ	1.0178	1.3189	1.033	1.3875
	C v (J/kg/K) b 1 b 0 (m 3 /kg) p ∞,1 (Pa/K)	3848 -0.5934 1.4905 × 10 -3 -607195	1719 0 3.3514 × 10 -4 0	1451 -0.6661 1.3013 × 10 -3 -405133	779 0 0 0
	p ∞,0 (Pa)	396642530	0	63642939	0
	q (J/kg) q ′ (J/kg/K) d (Pa m 3γ /kg γ )	-1065948 -20985 0	1975421 -3131 41200	-272675 -3277 0	-1597 2224 2950

http://webbook.nist.gov/chemistry/

Funding: This research received funding from ANR project SUBSUPERJET, ANR-14-CE22-0014.

Acknowledgments:

Part of this work has been carried out in the framework of the Labex MEC (ANR-10-LABX-0092) and of the A*MIDEX project (ANR-11-IDEX-0001-02), funded by the "Investissements d'Avenir" French Government program managed by the French National Research Agency (ANR).

Conflicts of Interest:

The authors declare no conflict of interest.

Appendix A Convexity of the ENASG formulation

Author Contributions:

Alexandre Chiapolino conceived the presented idea, developed the theory and performed the computations under the guidance and the supervision of Richard Saurel. Both authors discussed the results and contributed to the final manuscript.

Appendix C Methodology to determine the various Extended NASG (ENASG) parameters

This section details the procedure used in this work to determine the different ENASG parameters for liquid and gas phases. Depending on the studied application, the determination of the corresponding parameters can be different, but the use of the experimental curves is mandatory.

Liquid phase

In this section, let us introduce,

In the calculations that follow, these coefficients will be considered as known. Their numerical values will be addressed later. The liquid coefficients are determined with the help of an experimental saturation curve. In the following, the least squares method is used with the specific volume, Eq. [START_REF] Perigaud | A compressible flow model with capillary effects[END_REF]. Searching the optimum γ l coefficient, the next relation appears after some algebraic manipulations,

(A23)

In this section, let us introduce the following convention: ex p denotes the experimental values and N the number of experimental points considered. The experimental values of the internal energy are now used. Combining Eqs. ( 18), ( 20) and ( 21), the internal energy reads,

(A24)

A reference state re f is now used to express the liquid reference energy q l . Using (A24), the next relation arises,

Inserting Eq. (A25) into (A24), the internal energy transforms to,

) .

(A26) The least squares method is now applied to Eq. (A26). Searching the optimum C v,l coefficient, the following relation appears after some algebraic manipulations, -

,

,

(A27)

Equations ( A23) and (A27) then create a two-unknown equation system whose solution provides γ l and C v,l . An analytical solution is available and reads,

The expressions of γ l and C v,l are now available. Their numerical values will be determined with the help of the parameters A l , B l , C l and b 0,l , b 1,l . Those are addressed hereafter. From Eq. ( 17), the

. In this relation, the critical point is used via T c and p ′ ∞,c

that will be given arbitrary as p ′ ∞,c → 0. However, the coefficients b 0,l and b 1,l are included in B l and

The reference re f and critical c states are considered known either by an experimental point or arbitrary.

The coefficient B l is then known through b 0,l , b 1,l . The parameter A l depending only on C l , the coefficients γ l and C v,l become consequently, γ l = γ l (C l ) and C v,l = C v,l (C l ). The coefficient C l is then the only unknown at this point. To determine this latter, the speed of sound is used with another reference state denoted 0 (atmospheric conditions). With the help of Eq. ( 42), the next relation arises,

and can be solved with an iterative method. p ∞,0,l and p ∞,1,l are then determined via Eq. (A22), and the reference internal energy is computed with Eq. (A25). The reference entropy is the only unknown value at this point. The least squares method is used one more time with Eq. (38). Searching the optimum q ′ l coefficient, the following relation appears after some algebraic manipulations, Note that for the VdW and SRK EOSs, the specific internal energy, enthalpy and entropy require

v , q and q ′ as well for practical computations. In the present work

v,O 2 = 652 J/kg/K for both VdW and SRK EOSs. q H 2 O = 1799218 J/kg, q O 2 = 17918 J/kg, q ′ H 2 O = -3360 J/kg/K, q ′ O 2 = 2793 J/kg/K for the VdW EOS and q H 2 O = 1799885 J/kg, q O 2 = -682 J/kg, q ′ H 2 O = -3360 J/kg/K, q ′ O 2 = 2793 J/kg/K for the SRK EOS.

Appendix D Connection temperature between the ENASG EOS and ideal gas formulation

The different parameters of the quadratic equation ( 61) are provided hereafter for the specific volume (v), internal energy (e) and enthalpy (h). The solution of Eq. ( 61) provides the connection temperature between the ENASG EOS and ideal gas formulation during the liquid-to-supercritical state transition,

The subscripts l and g denote the liquid and gas phases respectively. Note that the determination of connection temperature for the entropy requires an iterative method because of the logarithmic function present in Eq. ( 56).

Appendix E Toward the critical point

Near the critical point, the vapor phase necessarily lacks of accuracy with the present ENASG EOS that is reduced to the ideal gas expression. The reason is linked to the absence of gas attractive effects. However, the introduction of those latter ones results in conditional convexity. They are thereby removed in this work. Nevertheless, they also result in much better agreement with experimental data as illustrated in the following.

Equation ( 16), recalled hereafter, does consider an attractive term via the parameter d,

Previously the coefficient d was set to zero for the sake of convexity and simplicity. This

1-b 1 extra term is reminiscent of cubic EOSs but seems nonetheless essential to describe dense gases near the critical point. With this parameter, the attractive pressures p ∞ (T) and p ′ ∞ (T) are no longer required for the gas phase and the thermal equation of state reads,

(A37) Equation (A37) does not provide an explicit formulation of the specific volume v(p, T). Cubic EOSs present the same flaw. Note that b(v) = b 0 = cst is accurate enough when dealing with gases. Following a mathematical procedure based on Maxwell's relations, similar to the derivations detailed in the previous sections, the present "alternative" formulation yields,

(A38) The v subscript denotes here the vapor phase. For the sake of space, the details of calculations are omitted. This "alternative" formulation respects Maxwell's relations and is thermodynamically consistent and convex under condition,

To represent correctly the physics of attractive terms, d v > 0 must be chosen and condition (A39) becomes restrictive. The same observation holds for cubic EOSs. It is worth mentioning that despite this conditional convexity, the speed of sound remains unambiguously positive unlike cubic EOSs and corresponds to the Noble-Abel (NA) sound speed. Table A3 provides the associated parameters of this "alternative" but conditionally convex formulation (A38) for water and oxygen. Figures A1 and A2 display the corresponding results at saturation. As seen in Figs. A1 and A2 the theoretical saturation pressure, liquid enthalpy, liquid and vapor specific densities are in very good agreement with experimental saturation data from the lowest temperature available to the critical one. Yet, the vapor enthalpy seems to present lesser agreement. It is interesting to note that vapor enthalpy is the only thermodynamic variable that presents a non-monotonic behavior. However, analyzing the range of variation, it seems that the theoretical order of magnitude is satisfied. The latent heat, very important during phase transition (L v = h vh l ) is also in very good agreement with experimental data. The introduction of the coefficient d involves a conditional convexity (Eq. (A39)) but does illustrate the significance of the gas attractive effects. 56) and its "alternative" but conditionally convex formulation for the vapor phase, Eq. (A38). p sat denotes the saturation pressure, L v the latent heat, h the specific enthalpy and ρ the density. 56) and its "alternative" but conditionally convex formulation for the vapor phase, Eq. (A38). p sat denotes the saturation pressure, L v the latent heat, h the specific enthalpy and ρ the density.