C Z Qiao

† H R Jiang

† S L Zhao
email: szhao@ecust.edu.cn

W Dong
email: wei.dong@ens-lyon.fr

A general approach based on morphological thermodynamics for a fluid confined in various porous media

We propose a general approach based on morphological thermodynamics for determining adsorption isotherms, i.e., the chemical potential of a confined fluid as a function of its density. The validity of this approach and its versatility are established by its remarkable accuracy compared to Monte-Carlo simulation results and its capability of accounting for a quite large variety of porous media, ranging from a simple slit pore to a random sponge matrix. It is also revealed that the contribution of the curvature terms to the chemical potential of the confined fluid is negligibly small when the interface curvature is not too large. This finding is of a particular importance for simplifying the treatment of experimental results of adsorption isotherms since no experimental technique is currently available for determining the curvatures of the pore surface inside a porous material.

Introduction

Porous materials are widely used in various domains, e.g., in chemical industry for molecular sieves and supported catalysts, in new clean energy technology for storing hydrogen, in some new therapy for long-lasting delivery of medicines by encapsulation. It is now well recognized that confinement can modify drastically some properties of adsorbed fluids. Accompanying the advent of many highperformance functionalized nanoporous materials, a large number of experimental and theoretical investigations have been made during the last decades. Nevertheless, a unifying picture highlighting the behavior of confined fluids emerges quite slowly due to the large diversity of structure and morphology for various porous materials. One salient characteristic of fluids confined in porous media is the large interface between fluid and pore wall, which has usually a complex morphology. Although one expects intuitively an important surface contribution to the thermodynamic potentials of such systems, it is not obvious whether it is still possible to define a meaningful surface tension when the characteristic pore size is of the order of a few molecular diameters of the confined fluid. Even if an approach based on thermodynamics is possible, one still needs to know which are the most relevant variables for characterizing the complex interface landscape of fluids adsorbed in porous media. The theoretical study of fluids confined in porous media is still carried out essentially in a case-by-case way. The morphological thermodynamics proposed and advocated by K. Mecke, R. Roth and co-workers [START_REF] König | Morphological Thermodynamics of Fluids: Shape Dependence of Free Energies[END_REF][START_REF] Mecke | Fluids in Porous Media: A Morphometric Approach[END_REF][START_REF] Roth | Morphometric Approach to the Solvation Free Energy of Complex Molecules[END_REF][START_REF] Oettel | Depletion Interaction of Two Spheres -Full Density Functional Theory vs. Morphometric Results[END_REF][START_REF] Kodama | Morphometric approach to thermodynamic quantities of solvation of complex molecules: Extension to multicomponent solvent[END_REF][START_REF] Jin | Shape Effect on Nanoparticle Solvation: A Comparison of Morphometric Thermodynamics and Microscopic Theories[END_REF][START_REF] Evans | Solvation of a sponge-like geometry[END_REF][START_REF] Reindl | Implications of Interface Conventions for Morphometric Thermodynamics[END_REF] provides a framework for a general thermodynamic description of complex interfacial systems. Starting from Hadwiger's theorem in integral geometry, morphological thermodynamics postulates that four geometrical measures are enough to characterize the thermodynamic potential of a complex interfacial system, i.e., volume, surface area, integrated mean curvature and integrated Gauss curvature [START_REF] König | Morphological Thermodynamics of Fluids: Shape Dependence of Free Energies[END_REF]. Moreover, Mecke et al assume that the thermodynamic variables conjugated to these geometrical measures, i.e. pressure, surface tension, as well as two surface bending rigidities, can be determined for a simple system then be used to describe systems with more complex morphology [START_REF] König | Morphological Thermodynamics of Fluids: Shape Dependence of Free Energies[END_REF][START_REF] Mecke | Fluids in Porous Media: A Morphometric Approach[END_REF][START_REF] Roth | Morphometric Approach to the Solvation Free Energy of Complex Molecules[END_REF]. Although morphological thermodynamics gives promising results for some systems [START_REF] König | Morphological Thermodynamics of Fluids: Shape Dependence of Free Energies[END_REF][START_REF] Mecke | Fluids in Porous Media: A Morphometric Approach[END_REF][START_REF] Roth | Morphometric Approach to the Solvation Free Energy of Complex Molecules[END_REF][START_REF] Oettel | Depletion Interaction of Two Spheres -Full Density Functional Theory vs. Morphometric Results[END_REF][START_REF] Kodama | Morphometric approach to thermodynamic quantities of solvation of complex molecules: Extension to multicomponent solvent[END_REF][START_REF] Jin | Shape Effect on Nanoparticle Solvation: A Comparison of Morphometric Thermodynamics and Microscopic Theories[END_REF][START_REF] Evans | Solvation of a sponge-like geometry[END_REF][START_REF] Reindl | Implications of Interface Conventions for Morphometric Thermodynamics[END_REF], the validity of one of its fundamental postulates has been questioned. Theoretical and simulation investigations have provided evidence for the existence of non-Hadwiger terms, i.e., high-order curvature contributions to surface tension, which are not taken into account in the morphological thermodynamics [START_REF] Laird | Interfacial Free Energy of a Hard-Sphere Fluid in Contact with Curved Hard Surfaces[END_REF][START_REF] Blokhuis | Existence of a Bending Rigidity for a Hard-Sphere Liquid near a Curved Hard Wall: Validity of the Hadwiger Theorem[END_REF][START_REF] Urrutia | Bending Rigidity and Higher-Order Curvature Terms for the Hard-Sphere Fluid near a Curved Wall[END_REF][START_REF] Hansen-Goos | Communication: Non-Hadwiger Terms in Morphological Thermodynamics of Fluids[END_REF].

Moreover, these investigations have shown that the first non-Hadwiger term gives a contribution much smaller than the Hadwiger terms, at least one order of magnitude smaller [START_REF] Laird | Interfacial Free Energy of a Hard-Sphere Fluid in Contact with Curved Hard Surfaces[END_REF]. Hence, morphological thermodynamics can allow for formulating useful approximations in practice even when Hadwiger theorem does not hold rigorously. Fluids adsorbed in porous materials with a complex morphology of pore space provide an interesting ground for further testing the applicability of morphological thermodynamics. This constitutes the main objective of the present work. Some of us (SLZ and WD) have collaborated during a long time with Myroslav Holovko in the study of fluids confined in porous media [START_REF] Holovko | A Highly Accurate and Analytic Equation of State for a Hard Sphere Fluid in Random Porous Media[END_REF][START_REF] Chen | Comment on "A highly accurate and analytic equation of state for a hard sphere fluid in random porous media[END_REF][START_REF] Patsahan | Fluids in Porous Media. III. Scaled Particle Theory[END_REF][START_REF] Holovko | Fluids in Random Porous Media: Scaled Particle Theory[END_REF][START_REF] Chen | Scaled particle theory for multicomponent hard sphere fluids confined in random porous media[END_REF]. It is our great pleasure to dedicate the present article to this special issue for the eightieth anniversary of Professor Holovko.

Our presentation is organized as follows. A brief introduction of morphological thermodynamics and a general equation of state for a fluid confined in various porous media will be presented in the next section. The accuracy of the general equation of state in different confining media will be assessed in Section III. Some conclusions will be presented in the last section.

Theory

A. A brief recall of morphological thermodynamics

The morphological thermodynamics proposed by Mecke, Roth and coworkers is based on the fundamental assumption that the grand potential of a complex interfacial system is a linear combination of four morphological measures, i.e., the system's volume, 𝑉𝑉, the interface area, 𝐴𝐴, the integrated mean curvature, 𝐶𝐶 M , and the integrated Gauss curvatures, 𝐶𝐶 G . From this assumption, one has immediately the following expression for the grand potential,

𝛺𝛺 = -𝑝𝑝 bulk (𝑇𝑇, 𝜇𝜇)𝑉𝑉 + 𝛾𝛾 0 (𝑇𝑇, 𝜇𝜇)𝐴𝐴 + 𝛾𝛾 -1 (𝑇𝑇, 𝜇𝜇)𝐶𝐶 M + 𝛾𝛾 -2 (𝑇𝑇, 𝜇𝜇)𝐶𝐶 G , ( 1
)
where 𝑇𝑇 is temperature, 𝜇𝜇 the chemical potential. The coefficients of the morphological measures are respectively the pressure of the corresponding bulk system, 𝑝𝑝 bulk , the surface tension on a flat surface, 𝛾𝛾 0 , and the bending rigidities, 𝛾𝛾 -1 and 𝛾𝛾 -2 related to the integrated mean and Gauss curvatures. Before applying it to study fluids confined in porous media, it is useful to clarify further some implications of the morphological thermodynamics given in eq.( 1). First, it is to note that the morphological measures are completely decoupled with the properties of the fluid, i.e., the four morphological measures being independent of the fluid state, i.e., 𝑇𝑇 and 𝜇𝜇. Such decoupling allows for separating the study of a complex interfacial system into two simpler tasks, one for determining the morphological measures in the absence of the fluid and the other for determining the coefficients, 𝛾𝛾 0 , 𝛾𝛾 -1 and 𝛾𝛾 -2 , from a simpler system, e.g., the considered fluid near a spherical surface. On the other hand, the independence of the morphological measures on the fluid state does not hold for fluid adsorption in very flexible porous materials since the adsorption can induces large deformation of the materials and modify significantly their morphology. So, the morphological thermodynamics as presented above cannot be applied for the fluid adsorption in flexible porous materials. We can also rewrite eq.( 1) as,

𝛺𝛺 = 𝛺𝛺 bulk + 𝛾𝛾 �𝐴𝐴 , (2)
where

𝛾𝛾 � = 𝛾𝛾 0 (𝑇𝑇, 𝜇𝜇) + 𝛾𝛾 -1 (𝑇𝑇, 𝜇𝜇) 𝐶𝐶 M 𝐴𝐴 + 𝛾𝛾 -2 (𝑇𝑇, 𝜇𝜇) 𝐶𝐶 G 𝐴𝐴 . (3)
Recently, one of us has shown that it is necessary to introduce the concepts of differential and integral surface tensions for strongly confined fluids, i.e., pore size being sufficiently small so that 𝛾𝛾 0 , 𝛾𝛾 -1 and 𝛾𝛾 -2 in eq.( 3) can also depends on the pore size in addition of 𝑇𝑇 and 𝜇𝜇 [START_REF] Dong | Thermodynamics of interface extended to nanoscales by introducing integral and differential surface tensions[END_REF][START_REF] Dong | Nanoscale thermodynamics needs the concept of a disjoining chemical potential[END_REF].

Although 𝛾𝛾 � defined from eq.( 2) is an integral surface tension, the assumption that 𝛾𝛾 0 , 𝛾𝛾 -1 and 𝛾𝛾 -2 are only functions of temperature and chemical potential is too restrictive to account for more complex situations. Therefore, we can anticipate that the morphological thermodynamics can become less and less accurate when confinement becomes stronger and stronger.

B. Adsorption isotherms based on morphological thermodynamics

The morphological thermodynamics extends Gibbs surface thermodynamics by proposing a concrete recipe to account for the curvature contributions in the grand potential. It is straightforward to show that the surface tension given in eq.( 3) satisfies Gibbs adsorption equation, i.e.,

� 𝜕𝜕𝛾𝛾 � 𝜕𝜕𝜕𝜕 � 𝑇𝑇,𝐶𝐶 M ,𝐶𝐶 G = 1 𝐴𝐴 �� 𝜕𝜕𝛺𝛺 𝜕𝜕𝜕𝜕 � 𝑇𝑇,𝐶𝐶 M ,𝐶𝐶 G -� 𝜕𝜕𝛺𝛺 bulk 𝜕𝜕𝜕𝜕 � 𝑇𝑇 � = - 𝑁𝑁 1 -𝑁𝑁 𝑏𝑏 𝐴𝐴 = -𝛤𝛤 , ( 4
)
where 𝛤𝛤 is defined as adsorption, 𝑁𝑁 1 and 𝑁𝑁 𝑏𝑏 are respectively the particle number of the confined fluid and that of the corresponding bulk fluid in the same volume with the same temperature and chemical potential. This equation provides a relation between the density of the confined fluid, 𝜌𝜌 1 = 𝑁𝑁 1 𝑉𝑉 -1 , and that of the corresponding bulk fluid, 𝜌𝜌 𝑏𝑏 = 𝑁𝑁 𝑏𝑏 𝑉𝑉 -1 ,

𝜌𝜌 𝑏𝑏 = 𝜌𝜌 1 + � 𝜕𝜕𝛾𝛾 � 𝜕𝜕𝜕𝜕 � 𝑇𝑇,𝐶𝐶 M ,𝐶𝐶 G 𝐴𝐴 𝑉𝑉 . (5)
Now, from the equation of state of the bulk fluid and the relation given by eq.( 5), we obtain immediately the following equation of state for the confined fluid,

𝜇𝜇(𝜌𝜌 1 ) = 𝜇𝜇 bulk �𝜌𝜌 1 + � 𝜕𝜕𝛾𝛾 � 𝜕𝜕𝜕𝜕 � 𝑇𝑇,𝐶𝐶 M ,𝐶𝐶 G 𝐴𝐴 𝑉𝑉 � . ( 6
)

Results and discussion

The theoretical framework presented in the last section is a quite general one which can be applied for studying a large variety of confined fluids. We illustrate this by considering some benchmark model systems.

A. A hard sphere fluid confined in a slit pore with hard walls

We consider first a hard sphere (HS) fluid confined in a slit pore formed with two parallel hard walls. A schematic presentation of such a system is given in Figure 1. The morphological thermodynamic approach described in Section 2.B is based on two prerequisites: i) an equation of state for the bulk fluid; ii) the surface tension for the considered interface. For a HS fluid confined in a slit pore with two hard walls, scaled particle theory (SPT) [START_REF] Reiss | Statistical mechanics of rigid spheres[END_REF][START_REF] Dong | Scaled particle theory for bulk and confined fluids: A review[END_REF][START_REF] Qiao | Augmented scaled particle theory[END_REF] provides both the equation of state of the bulk fluid and the surface tension of a HS fluid near a flat hard wall. The SPT equation of state for a bulk HS fluid is given by,

𝛽𝛽𝜇𝜇 SPT = ln(𝛬𝛬 3 𝜌𝜌 𝑏𝑏 ) -ln(1 -𝜂𝜂 𝑏𝑏 ) + 7𝜂𝜂 𝑏𝑏 1-𝜂𝜂 𝑏𝑏 + 15 2 � 𝜂𝜂 𝑏𝑏 1-𝜂𝜂 𝑏𝑏 � 2 + 3 � 𝜂𝜂 𝑏𝑏 1-𝜂𝜂 𝑏𝑏 � 3 . ( 7
)
where 𝛬𝛬 is the thermal wavelength, 𝛽𝛽 = (𝑘𝑘 B 𝑇𝑇) -1 (𝑘𝑘 B : Boltzmann constant), 𝜂𝜂 𝑏𝑏 = 𝜋𝜋𝜎𝜎 3 𝜌𝜌 𝑏𝑏 6 -1

is the packing fraction of the bulk fluid. The surface tension for a fluid-wall interface and the adsorption depend on the choice of Gibbs dividing surface. For the present work, we choose the wall surfaces located at 𝑧𝑧 = ± 𝐿𝐿 2 ⁄ as the dividing surfaces (see Fig. 1). For this choice, the SPT gives the following result for the surface tension,

𝜋𝜋𝜎𝜎 2 𝛽𝛽𝛾𝛾 0 SPT = 3𝜂𝜂 𝑏𝑏 1-𝜂𝜂 𝑏𝑏 + 9 2 � 𝜂𝜂 𝑏𝑏 1-𝜂𝜂 𝑏𝑏 � 2 . ( 8
)
From the Gibbs adsorption equation, i.e., eq.( 4), and the results given in eqs.( 7) and ( 8), we obtain straightforwardly the adsorption on one pore wall,

𝜋𝜋𝜎𝜎 2 𝛤𝛤 0 SPT = -𝜋𝜋𝜎𝜎 2 � 𝜕𝜕𝛾𝛾 0 SPT 𝜕𝜕𝜕𝜕 SPT � 𝑇𝑇 = -𝜋𝜋𝜎𝜎 2 � 𝜕𝜕𝛾𝛾 0 SPT 𝜕𝜕𝜂𝜂 𝑏𝑏 � 𝑇𝑇 � 𝜕𝜕𝜕𝜕 SPT 𝜕𝜕𝜂𝜂 𝑏𝑏 � 𝑇𝑇 � = - 3𝜂𝜂 𝑏𝑏 �1-𝜂𝜂 𝑏𝑏 � 1+2𝜂𝜂 𝑏𝑏 . ( 9
)
Since the slit pore is composed of two hard walls, we have the following relation between the packing fraction of the confined fluid, 𝜂𝜂 1 = 𝜋𝜋𝜎𝜎 3 𝜌𝜌 1 6 -1 , and that of the bulk one,

𝜂𝜂 1 = 𝜂𝜂 𝑏𝑏 + 2𝜋𝜋𝜎𝜎 2 𝛤𝛤 0 SPT 6𝐿𝐿 * = 𝜂𝜂 𝑏𝑏 - 𝜂𝜂 𝑏𝑏 �1-𝜂𝜂 𝑏𝑏 � 𝐿𝐿 * �1+2𝜂𝜂 𝑏𝑏 � , (10)
where 𝐿𝐿 * = 𝐿𝐿𝜎𝜎 -1 is the pore width measured with the HS diameter. Since eq.( 10) is a second order polynomial, we obtain the following explicit expression of 𝜂𝜂 𝑏𝑏 in terms of 𝜂𝜂 1 ,

𝜂𝜂 𝑏𝑏 = 1-𝐿𝐿 * +2𝐿𝐿 * 𝜂𝜂 1 +�(1-𝐿𝐿 * +2𝐿𝐿 * 𝜂𝜂 1 ) 2 +4𝐿𝐿 * (1+2𝐿𝐿 * )𝜂𝜂 1 2(1+2𝐿𝐿 * ) . ( 11
)
Now, the adsorption isotherm for a HS fluid confined in a slit pore is given by,

𝛽𝛽𝜇𝜇(𝜂𝜂 1 ) = 𝛽𝛽𝜇𝜇 bulk (𝜂𝜂 𝑏𝑏 ) . ( 12
)
The right-hand-side of eq.( 12) is given by substituting eq.( 11) into eq.( 7). Eq. ( 10) shows clearly that with the increase of 𝐿𝐿 * , the difference between 𝜂𝜂 1 and 𝜂𝜂 𝑏𝑏 becomes smaller and smaller and lim 𝐿𝐿 * →∞ 𝜂𝜂 1 = 𝜂𝜂 𝑏𝑏 , one recovers the bulk result from eq.( 12) for 𝐿𝐿 * → ∞ as one can expect. For point particle, i.e. 𝜎𝜎 = 0, eqs. ( 7) and ( 12) gives the same exact result for the confined and the bulk ideal gas as we expect.

In Figure 2, some results given by eq. ( 12) are presented for a few pore widths and compared to the results of Monte-Carlo (MC) simulations. Since the adsorption for this system is negative, the adsorption isotherms of the confined fluid are always above the isotherm of the bulk fluid. As the surface contribution decreases when the pore width becomes larger, the isotherm of the confined fluid approaches more and more the isotherm of the bulk fluid. The agreement between the results of morphological thermodynamics and those of the MC simulations is very good for 𝐿𝐿 * ≥ 3.5. results of morphological thermodynamics deviate more and more from those of the MC simulations.

For very narrow pores, a disjoining pressure [START_REF] Derjaguin | Direct measurement of molecular forces[END_REF][START_REF] Derjaguin | Surface forces[END_REF] arises in the confined fluid. One of us (WD) has shown that due to the contribution of the disjoining pressure, the surface tension for very narrow pores is no longer equal to that for the fluid near one single wall [START_REF] Dong | Thermodynamics of interface extended to nanoscales by introducing integral and differential surface tensions[END_REF]. However, the contribution of the disjoining pressure is not accounted for by the morphological thermodynamics. We believe this is the main reason for the failure of morphological thermodynamics in describing accurately the strongly confined fluids at high densities in very narrow slit pores.

Slit pore is a model extensively used for studying confined fluids. Many simulations have been carried out for fluids confined in a slit pore. Labik and Smith reported NVT-ensemble Monte Carlo simulation results for hard spheres in a hard slit pore [START_REF] Labik | Computer simualtion of the chemical potential and adsorption isotherm of hard spheres in a hard slit-like pore[END_REF] and our simulation results are in good agreement with theirs. Smith and coworkers [START_REF] Pospisil | Integral euqation and computer simulation studies of hard spheres in a slit pore[END_REF] and Alejandre et al [START_REF] Alejandre | Effect of pore geometry on a confined hard sphere fluid[END_REF] have used integral equations to study a hard sphere fluid in a slit pore. Many investigations on a variety of fluids confined in a slit pore have been made with the help of density functional theory (DFT) and many of them aim at determining the fluid distribution inside the pore (see e.g., [START_REF] Yu | Amodified fundamental measure theory for spherical particles in microchannels[END_REF], an exhaustive review is beyond the scope of the present article). The theoretical approaches based on DFT or integral equations requires first determining the one-and two-body distribution functions. Our approach in this work focuses only on the thermodynamic properties, which requires only the equation of state of the considered fluid in the bulk phase and the surface tension for the considered interface. When SPT is used for these properties, we obtain a totally analytical result for the adsorption isotherms of a hard sphere fluid confined in a slit pore.

B. A hard sphere fluid confined in an ordered or a disordered hard sphere matrix

The porous matrix model with quenched matrix particles proposed by Madden and Glandt [START_REF] Madden | Distribution Functions for Fluids in Random Media[END_REF] accounts for more characteristics of porous media than the simple slit pore model, e.g., pore connectivity, curved pore surface, variation of pore size inside a porous material. Figure 4 illustrates three types of porous matrices: i) ordered hard sphere matrix with matrix particles fixed on a lattice;

ii) disordered hard sphere matrix with matrix particles quenched from an equilibrium configuration of a fluid [START_REF] Madden | Distribution Functions for Fluids in Random Media[END_REF]; iii) overlapping hard sphere matrix with matrix particles distributed totally randomly in space. By construction, the morphological thermodynamics cannot distinguish the matrices of types i and ii. When the dividing surface is chosen as the surface of HS matrix particles, the fluidmatrix interface area of type i or ii matrices can be calculated easily and we have,

𝐴𝐴 = 𝜋𝜋𝜎𝜎 0 2 𝑁𝑁 0 , ( 13
)
where 𝑁𝑁 0 and 𝜎𝜎 0 are respectively the number and diameter of the matrix particle. The fundamental assumption of the morphological thermodynamics for treating the type i and type ii matrices is to describe the adsorption as the sum of adsorption around each matrix particle.

SPT gives the following result for the adsorption around the surface of one matrix particle,

𝜋𝜋𝜎𝜎 0 2 𝛤𝛤 𝑠𝑠 SPT = -𝜋𝜋𝜎𝜎 2 � 𝜕𝜕𝛾𝛾 SPT 𝜕𝜕𝜇𝜇 SPT � 𝑇𝑇,𝜎𝜎 0 = -𝜋𝜋𝜎𝜎 2 � 𝜕𝜕𝛾𝛾 SPT 𝜕𝜕𝜂𝜂 𝑏𝑏 � 𝑇𝑇,𝜎𝜎 0 � 𝜕𝜕𝜇𝜇 SPT 𝜕𝜕𝜂𝜂 𝑏𝑏 � 𝑇𝑇,𝜎𝜎 0 � = - 3𝜂𝜂 𝑏𝑏 �1-𝜂𝜂 𝑏𝑏 � 1+2𝜂𝜂 𝑏𝑏 �𝜏𝜏 2 + 𝜏𝜏�1-𝜂𝜂 𝑏𝑏 � 1+2𝜂𝜂 𝑏𝑏 + �1-𝜂𝜂 𝑏𝑏 � 2 3�1+2𝜂𝜂 𝑏𝑏 � � , ( 14
)
where 𝜏𝜏 = 𝜎𝜎 0 𝜎𝜎 ⁄ is the size ratio between the matrix and fluid particles (𝜎𝜎 0 : matrix particle diameter). Besides the contribution given by eq.( 14), there is also an additional contribution to the adsorption given by,

𝜋𝜋𝜎𝜎 0 2 𝛤𝛤 0 = �𝜙𝜙 0 HS -1�𝜌𝜌 𝑏𝑏 𝑉𝑉 𝑁𝑁 0 , ( 15
)
where 𝑉𝑉 is the volume of the matrix sample and 𝜙𝜙 0 HS the geometric porosity for a hard sphere matrix,

𝜙𝜙 0 HS = 1 - 𝜋𝜋𝜎𝜎 0 3 𝑁𝑁 0 6𝑉𝑉 = 1 -𝜂𝜂 0 𝜏𝜏 3 , ( 16
)
with 𝜂𝜂 0 = 𝜋𝜋𝜎𝜎 3 𝑁𝑁 0 𝑉𝑉 -1 6 -1 . Now, the relation between the packing fraction of the confined fluid and that of the bulk fluid with the same chemical potential is given by,

𝜂𝜂 1 = 𝜙𝜙 0 HS 𝜂𝜂 𝑏𝑏 - 3𝜂𝜂 0 𝜂𝜂 𝑏𝑏 �1-𝜂𝜂 𝑏𝑏 � �1+2𝜂𝜂 𝑏𝑏 � �𝜏𝜏 2 + 𝜏𝜏�1-𝜂𝜂 𝑏𝑏 � 1+2𝜂𝜂 𝑏𝑏 + �1-𝜂𝜂 𝑏𝑏 � 2 3�1+2𝜂𝜂 𝑏𝑏 � � . ( 17
)
Unlike the case of a slit pore, it is not possible to solve eq.( 17) for obtaining an analytical expression of 𝜂𝜂 𝑏𝑏 in terms of 𝜂𝜂 1 like eq.( 11). Nevertheless, it is still quite easy to obtain the adsorption isotherm of the confined fluid with the following procedure. With a given 𝜌𝜌 𝑏𝑏 , one obtains a value of chemical potential, 𝜇𝜇, and the fluid density of the confined fluid having the same chemical potential is given by eq.( 17), thus the relation between 𝜇𝜇 and 𝜂𝜂 1 is found. When the size of the fluid HS shrinks to zero, i.e., 𝜎𝜎 = 0, eq.( 17) becomes 𝜌𝜌 1 = 𝜙𝜙 0 HS 𝜌𝜌 𝑏𝑏 and this allows for obtaining the exact result for an ideal gas confined in a HS matrix, i.e., 𝛽𝛽𝜇𝜇 = ln�𝛬𝛬 3 𝜌𝜌 1 𝜙𝜙 0 HS ⁄ �. So, without the additional contribution to the adsorption given in eq.( 15), it is impossible to recover the exact ideal gas results in the limit of point particles. A more detailed discussion about this point is presented in Appendix A.

The results for a HS fluid confined in two ordered HS matrices with different matrix particle sizes and matrix densities are presented in Fig. We consider next the effect of matrix disorder on the chemical potential of the confined fluid.

The results for a HS fluid confined in a disordered HS matrix are presented in Fig. 6. By construction, the morphological thermodynamics does not take the effect of matrix disorder into account. However, our Monte-Carlo simulation results given in Fig. 6 are performed for disordered HS matrices (see Appendix B for simulation methods and computation conditions). The very good agreement between the results of the morphological thermodynamics and the simulation ones show that the effect of matrix disorder on the chemical potential is negligibly small. This conforms also to the same finding of our recent work [START_REF] Qiao | Connect the Thermodynamics of Bulk and Confined Fluids: Confinement-Adsorption Scaling[END_REF]. In Fig. 6, the results of the morphological thermodynamics are also compared to those of SPT2b1 theory we proposed previously [START_REF] Holovko | Fluids in Random Porous Media: Scaled Particle Theory[END_REF], which is an approach developed specifically for disordered porous matrices while the morphological thermodynamics is a general approach applicable for a larger variety of inhomogeneous fluids. III of [START_REF] Patsahan | Fluids in Porous Media. III. Scaled Particle Theory[END_REF]).

C. A hard sphere fluid confined in an overlapping hard sphere matrix

The application of morphological thermodynamics to this case follows the general procedure described above. So, we need first determine the porosity of the matrix for a point particle, 𝜙𝜙 0 , which is required for determining ideal gas contribution to the chemical potential.

𝜙𝜙 0 OHS = 1 𝑉𝑉𝑉𝑉 𝑁𝑁 0 ∫ 𝑑𝑑𝒓𝒓 ∏ ∫ 𝒒𝒒 𝑖𝑖 𝑁𝑁 0 𝑖𝑖=1 e -𝛽𝛽 ∑ 𝑢𝑢 fm ��𝒓𝒓-𝒒𝒒 𝑗𝑗 �� 𝑁𝑁 0 𝑗𝑗=1 = 1 𝑉𝑉 𝑁𝑁 0 �𝑉𝑉 - 4𝜋𝜋𝑅𝑅 0 3 3 � 𝑁𝑁 0 = exp �- 4𝜋𝜋𝑅𝑅 0 3 𝜌𝜌 0 3 � = exp(-𝜂𝜂 0 𝜏𝜏 3 ) , ( 18
)
where 𝒓𝒓 is the position vector of the point particle and 𝒒𝒒 𝑖𝑖 the position vector of the i-th matrix particle. The interaction between the point particle and a matrix particle is given by,

𝑢𝑢 fm ��𝒓𝒓 -𝒒𝒒 𝑗𝑗 �� = � ∞ , �𝒓𝒓 -𝒒𝒒 𝑗𝑗 � < 𝑅𝑅 0 0 , �𝒓𝒓 -𝒒𝒒 𝑗𝑗 � ≥ 𝑅𝑅 0 . ( 19
)
with 𝑅𝑅 0 = 𝜎𝜎 0 2 ⁄ . The thermodynamic limit is taken, i.e., lim 𝑁𝑁 0 →∞,𝑉𝑉→∞ 𝑁𝑁 0 𝑉𝑉 ⁄ = 𝜌𝜌 0 , when going to the third equality of eq.( 18). The interface area is given by,

𝐴𝐴 OHS = 𝑑𝑑�𝑉𝑉-𝜙𝜙 0 OHS 𝑉𝑉� 𝑑𝑑𝑅𝑅 0 = 𝑉𝑉4𝜋𝜋𝑅𝑅 0 2 𝜌𝜌 0 exp �- 4𝜋𝜋𝑅𝑅 0 3 𝜌𝜌 0 3 � = 𝑉𝑉4𝜋𝜋𝑅𝑅 0 2 𝜌𝜌 0 𝜙𝜙 0 OHS . ( 20
)
The integrated mean and Gauss curvatures are given respectively by,

𝐶𝐶 -1 OHS = 1 2 𝑑𝑑𝐴𝐴 𝑑𝑑𝑅𝑅 0 = 𝑉𝑉4𝜋𝜋𝑅𝑅 0 𝜌𝜌 0 �1 -2𝜋𝜋𝑅𝑅 0 3 𝜌𝜌 0 �exp �- 4𝜋𝜋𝑅𝑅 0 3 𝜌𝜌 0 3 � = 𝑉𝑉4𝜋𝜋𝑅𝑅 0 𝜌𝜌 0 �1 - 3 2 𝜂𝜂 0 𝜏𝜏 3 � 𝜙𝜙 0 OHS , (21)
𝐶𝐶 -2 OHS = 1 2 𝑑𝑑 2 𝐴𝐴 𝑑𝑑𝑅𝑅 0 2 = 𝑉𝑉4𝜋𝜋𝜌𝜌 0 �1 -12𝜋𝜋𝑅𝑅 0 3 𝜌𝜌 0 + 8𝜋𝜋𝑅𝑅 0 6 𝜌𝜌 0 2 �exp �- 4𝜋𝜋𝑅𝑅 0 3 𝜌𝜌 0 3 � = 𝑉𝑉4𝜋𝜋𝜌𝜌 0 �1 -9𝜂𝜂 0 𝜏𝜏 3 + 9 2 𝜂𝜂 0 2 𝜏𝜏 6 � 𝜙𝜙 0 OHS . (22)
Finally, the morphological thermodynamics combined with SPT gives the following relation between the density of the confined fluid to that of the bulk fluid which has the same chemical potential,

𝜂𝜂 1 = 𝜙𝜙 0 OHS 𝜂𝜂 𝑏𝑏 - 3𝜂𝜂 0 𝜂𝜂 𝑏𝑏 �1-𝜂𝜂 𝑏𝑏 � �1+2𝜂𝜂 𝑏𝑏 � 𝜙𝜙 0 OHS �𝜏𝜏 2 + 𝜏𝜏�1-𝜂𝜂 𝑏𝑏 � 1+2𝜂𝜂 𝑏𝑏 �1 - 3 2 𝜂𝜂 0 𝜏𝜏 3 � �+ �1-𝜂𝜂 𝑏𝑏 � 2 3�1+2𝜂𝜂 𝑏𝑏 � � + �1-𝜂𝜂 𝑏𝑏 � 2 3�1+2𝜂𝜂 𝑏𝑏 � �1 -9𝜂𝜂 0 𝜏𝜏 3 + 9 2 𝜂𝜂 0 2 𝜏𝜏 6 �� . ( 23
)
The results for a HS fluid confined in an overlapping hard sphere matrix are presented in Fig. 7.

The results given by morphological thermodynamics are again in excellent agreement with the simulation ones under all the conditions we have studied.

D. A hard sphere fluid confined in a hard sponge matrix

The last system we consider is a hard sphere fluid confined in a hard sponge model proposed by two of us (SLZ and WD) with Q. H. Liu [START_REF] Zhao | Fluids in Porous Media. I. A Hard Sponge Model[END_REF]. Figure 8 gives a schematic presentation of this model.

In this case, the interaction potential between a fluid particle and the matrix is a non-additive nbody one given by,

𝑣𝑣 fm �𝒓𝒓, 𝒒𝒒 1 , 𝒒𝒒 2 , . . . , 𝒒𝒒 𝑁𝑁 0 � = -k B 𝑇𝑇ln �1 -e -𝛽𝛽 ∑ 𝑢𝑢 fm ��𝒓𝒓-𝒒𝒒 𝑗𝑗 �� 𝑁𝑁 0 𝑗𝑗=1 � , (24)
where 𝑢𝑢 fm ��𝒓𝒓 -𝒒𝒒 𝑗𝑗 �� is given by eq.( 19) and now 𝜎𝜎 0 is the diameter of a spherical cavity in the sponge matrix. Although this fluid-matrix interaction potential is non-additive, it is still possible to obtain an analytical result for the porosity, i.e.,

𝜙𝜙 0 HSG = 1 𝑉𝑉𝑉𝑉 𝑁𝑁 0 ∫ 𝑑𝑑𝒓𝒓 ∏ ∫ 𝒒𝒒 𝑖𝑖 𝑁𝑁 0 𝑖𝑖=1 e -𝛽𝛽𝑣𝑣 fm �𝒓𝒓,𝒒𝒒 1 ,𝒒𝒒 2 ,...,𝒒𝒒 𝑁𝑁 0 � = 1 𝑉𝑉𝑉𝑉 𝑁𝑁 0 ∫ 𝑑𝑑𝒓𝒓 ∏ ∫ 𝒒𝒒 𝑖𝑖 𝑁𝑁 0 𝑖𝑖=1 �1 -e -𝛽𝛽 ∑ 𝑢𝑢 fm ��𝒓𝒓-𝒒𝒒 𝑗𝑗 �� 𝑁𝑁 0 𝑗𝑗=1 � = 1 - 1 𝑉𝑉 𝑁𝑁 0 �𝑉𝑉 - 𝜋𝜋𝜎𝜎 0 3 6 � 𝑁𝑁 0 = 1 -exp �- 𝜋𝜋𝜎𝜎 0 3 𝜌𝜌 0 6 � = 1 -𝜙𝜙 0 OHS . ( 25
)
The thermodynamic limit is taken, i.e., lim 𝑁𝑁 0 →∞,𝑉𝑉→∞ 𝑁𝑁 0 𝑉𝑉 ⁄ = 𝜌𝜌 0 , when going to the fourth equality of eq.( 25). It is not difficult to see that if the cavity diameter is the same as that of the matrix particle in an OHS matrix and their number densities are also equal, the following relations hold between the interface area, integrated mean and Gauss curvatures of the hard sponge matrix and those of the OHS matrix,

𝐴𝐴 HSG = 𝐴𝐴 OHS , ( 26
)
𝐶𝐶 -1 HSG = -𝐶𝐶 -1 OHS , ( 27
)
𝐶𝐶 -2 HSG = 𝐶𝐶 -2 OHS . ( 28
)
The relation given in eq.( 27) reflects the simple fact that the confined fluid is adsorbed on a convex surface in a OHS matrix but on a concave surface in a hard sponge matrix. In figure 9, the results of the morphological thermodynamics for a HS fluid in a hard sponge matrix are presented along with simulation results under different conditions. The accuracy of our approach based on the morphological thermodynamics is again remarkable.

D. Contribution of surface curvatures

The comparison of the respective results for OHS and hard sponge matrices indicates that the contribution of the surface curvature terms to the adsorption isotherm can be quite small. More detailed analyses show that when the matrix particle is three times larger than the fluid particle, the contribution of the surface curvatures is negligible. Figure 10 shows that even for the size ratio (𝜏𝜏 = 𝜎𝜎 0 𝜎𝜎 ⁄ ) equal to 1, 2 or 3, the contribution of the surface curvatures is quite moderate. Bryk et al have studied the adsorption of hard spheres confined between two uniaxial cylinders by using DFT. They compared the results for such a confined fluid with those for a HS fluid confined in a slit pore of two flat walls and found that the surface curvature effect is quite small [START_REF] Bryk | Adsorption of fluids in pores formed between two hard cylinders[END_REF]. Our results presented in Fig. 10 are consistent with their finding. To the best of our knowledge, there exist no experimental technique for measuring the interface curvatures of any porous materials. So, our finding here provides a simplified procedure to interpret the experimental results for adsorption isotherms by neglecting the contribution of surface curvature terms.

Conclusion

In the present work, we propose a general approach based on the morphological thermodynamics for determining the chemical potential of a fluid confined in a large variety of porous media, from a simple slit pore to a random hard sponge matrix. Our approach requires an equation of state of the considered fluid in a bulk phase and the surface tension of the fluid on a wall with a much simpler morphology than the complex porous medium under consideration. For the hard sphere fluid confined in the various hard porous media considered in this work, scaled particle theory gives both the equation of state of the bulk fluid and the surface tension for a HS fluid on a hardsphere wall. These are enough for constructing a totally analytical approach for all porous media considered in the present work. The comparison with simulation results show that the overall accuracy of our approach is excellent. Moderate discrepancies are found only for very narrow slit pores (𝐿𝐿 ≤ 3𝜎𝜎) and at high fluid densities.

Although previous investigations have shown that non-Hadwiger terms (high-order curvature terms not included in morphological thermodynamics) do not vanish rigorously [START_REF] Laird | Interfacial Free Energy of a Hard-Sphere Fluid in Contact with Curved Hard Surfaces[END_REF][START_REF] Blokhuis | Existence of a Bending Rigidity for a Hard-Sphere Liquid near a Curved Hard Wall: Validity of the Hadwiger Theorem[END_REF][START_REF] Urrutia | Bending Rigidity and Higher-Order Curvature Terms for the Hard-Sphere Fluid near a Curved Wall[END_REF][START_REF] Hansen-Goos | Communication: Non-Hadwiger Terms in Morphological Thermodynamics of Fluids[END_REF], their contributions are usually smaller by one order of magnitude. In the present work, we find that even the integrated mean and Gauss curvature terms have a negligible contribution to the chemical potential of the confined fluid when the surface curvature is not too large. So, even the contribution of the integrated mean and Gauss curvatures can be neglected in many cases. This simplify significantly the treatment of the experimental results for the adsorption isotherms in practice since the experimental technique is currently lacking for measuring the interface curvatures inside a porous material. Thus, it is in principle possible to elaborate an experimental method for determining the surface tension for the interface between a fluid and the pore wall inside a porous material, with the help of its adsorption isotherms. The well-known Gibbs adsorption equation implies that a nonzero adsorption should lead to a nonzero surface tension. So, this result seems to be in contradiction with that of eq.(A5) which

shows the surface tension is zero. Now, we will show that there is in fact no contradiction. It is recently revealed that two surface tensions can arises, one is the differential surface tension and the other is the integral surface tension [START_REF] Dong | Thermodynamics of interface extended to nanoscales by introducing integral and differential surface tensions[END_REF][START_REF] Dong | Nanoscale thermodynamics needs the concept of a disjoining chemical potential[END_REF]. Moreover, it is the integral surface tension which satisfies a generalized Gibbs adsorption equation [START_REF] Dong | Thermodynamics of interface extended to nanoscales by introducing integral and differential surface tensions[END_REF]. The integral surface tension includes a contribution from the differential surface tension and another contribution from the disjoining pressure. We calculate now the disjoining pressure and show it is indeed nonzero for the system considered here. The pressures of where ℓ � = 𝑉𝑉 𝐴𝐴 ⁄ .

Finally, we show that accounting adequately for the integral surface tension and the related different matrix samples.

As for the slit pore, each simulation for a given matrix sample includes about 2 × 10 5 Monte-Carlo cycles for preparing the system to equilibrium and 10 6 production cycles for determining the chemical potential. The simulation code can be found on GitHub (https://github.com/qiaochongzhi/MC-for-ConfinedFluid).

Figure 1 .

 1 Figure 1. Schematic presentation of a fluid of hard spheres (blue spheres) confined in a slit pore with two walls (grey) separated by a distance of 𝐿𝐿 (surface normal along the z-direction). The dashed lines mark the closest accessible plans for the hard-sphere centers, ℓ = 𝐿𝐿 -𝜎𝜎 being the accessible width (𝜎𝜎: fluid hard sphere diameter).

Figure 2 .

 2 Figure 2. Chemical potential of a hard sphere fluid confined in a slit pore as a function of fluid density from morphological thermodynamics, i.e., eq. (12) (continuous line) and NVT-ensemble Monte-Carlo simulations (symbols, see Appendix B for the details about simulation method and computation conditions). 1) 𝐿𝐿 * = 3.5 (red); 2) 𝐿𝐿 * = 5.0 (blue); 3) 𝐿𝐿 * = 7.5 (yellow); 4) 𝐿𝐿 * → ∞ (black).

Figure 3 Figure 3 .

 33 Figure3shows the results for some narrower pores, 𝐿𝐿 * ≤ 3.0. Up to the moderate fluid density, i.e., 𝜌𝜌 1 𝜎𝜎 3 ≤ 0.6, morphological thermodynamics gives accurate results. For higher densities, the

Figure 4 .

 4 Figure 4. Schematic presentation of a fluid of hard spheres (blue spheres) confined in different porous matrices composed quenched matrix particles (grey). a) ordered porous matrix with hard sphere matrix particles fixed on a simple cubic lattice; b) hard sphere matrix quenched from an equilibrium liquid configuration; c) overlapping hard sphere matrix with matrix particles distributed totally randomly in space.

  . The analytical approach by combining morphological thermodynamics and SPT gives excellent results for the chemical potential, as evidenced by the comparison with the results of the Monte-Carlo simulations we performed. It is to note that we define the density of a fluid confined in a porous matrix with respect to the sample volume which contains fluid and matrix particles. Due to this definition, the range of fluid densities considered for fluids confined in porous matrices appears smaller than that considered for a fluid confined in a slit pore. For example, in Figs.2 and 3the considered fluid density goes up to 0.8 while the fluid density considered in Fig.5goes only to 0.6. In fact, to make a more plausible comparison between the fluid density in a porous matrix and that in a slit pore, we should take into account the porosity of the porous matrix. For the matrix with 𝜎𝜎 0 𝜎𝜎 ⁄ = 10 and 𝜌𝜌 0 𝜎𝜎 3 = 5.7804 × 10 -4 in Fig.5, the porosity is given by 𝜙𝜙 0 = 1 -𝜋𝜋𝜌𝜌 0 𝜎𝜎 0 3 6 ⁄ = 0.69734. Dividing the fluid density defined with respect to the sample volume by the porosity give the fluid density defined with respect to the void volume. So, 𝜌𝜌 1 𝜎𝜎 3 = 0.6 corresponds in fact to a fluid density defined with respect to void volume equal to 𝜌𝜌 1 𝜎𝜎 3 𝜙𝜙 0 ⁄ = 0.86.

Figure 5 .

 5 Figure 5. Chemical potential of a hard sphere fluid confined in an ordered porous matrix, with hard sphere matrix particles fixed on a simple cubic lattice, as a function of fluid density. Matrix to fluid particle size ratio: 𝜎𝜎 0 𝜎𝜎 ⁄ = 5, 10 ; Matrix density: 𝜌𝜌 0 𝜎𝜎 3 = 0.001, 5.7804 × 10 -4 ; Red lines: Morphological thermodynamics combined with SPT; Black squares: NVT-ensemble Monet-Carlo simulation (see Appendix B for the details about simulation method and computation conditions).

Figure 6 .

 6 Figure 6. Chemical potential of a hard sphere fluid confined in a disordered hard sphere matrix as a function of fluid density. a) Matrix to fluid particle size ratio: 𝜎𝜎 0 𝜎𝜎 ⁄ = 2 ; Matrix density: 𝜌𝜌 0 𝜎𝜎 3 = 0.05 ; b) Matrix to fluid particle size ratio: 𝜎𝜎 0 𝜎𝜎 ⁄ = 3 ; Matrix density: 𝜌𝜌 0 𝜎𝜎 3 = 0.0149208, 0.019397; c) Matrix to fluid particle size ratio: 𝜎𝜎 0 𝜎𝜎 ⁄ = 5; Matrix density: 𝜌𝜌 0 𝜎𝜎 3 = 0.0044421, 0.0057473; Red lines: Morphological thermodynamics combined with SPT; Blue dashed lines: SPT2b1 theory [16]; Black squares: Monte-Carlo simulation (Results in (a) from our NVT-ensemble Monte-Carlo simulations with the details of simulation method and computation conditions given in Appendix B, those in (b) and (c) from the µVT-ensemble Monte-Carlo simulations given in TableIIIof[START_REF] Patsahan | Fluids in Porous Media. III. Scaled Particle Theory[END_REF]).

Figure 7 . 1 ; 3 ; 5 ;

 7135 Figure 7. Chemical potential of a hard sphere fluid confined in an overlapping hard sphere matrix as a function of fluid density. Red lines: Morphological thermodynamics combined with SPT; Blue dashed lines: SPT2b1 theory [16]; Black squares: Monte-Carlo simulation (Results in (a) from our NVT-ensemble Monte-Carlo simulations with the details of simulation method and computation conditions given in Appendix B, those in (b) and (c) from µVT-ensemble Monte-Carlo simulations of [15]). a) Matrix to fluid particle size ratio: 𝜎𝜎 0 𝜎𝜎 ⁄ = 1; Matrix density: 𝜌𝜌 0 𝜎𝜎 3 = 0.1590 (the red line and the blue dashed lines nearly overlap each other); b) Matrix to fluid particle size ratio: 𝜎𝜎 0 𝜎𝜎 ⁄ = 3; Matrix density: 𝜌𝜌 0 𝜎𝜎 3 = 0.020683, 0.031328; c) Matrix to fluid particle size ratio: 𝜎𝜎 0 𝜎𝜎 ⁄ = 5; Matrix density: 𝜌𝜌 0 𝜎𝜎 3 = 0.006128, 0.0092824.

Figure 8 .

 8 Figure 8. Schematic presentation of a fluid of hard spheres (blue spheres) confined in a hard sponge

Figure 9 .

 9 Figure 9. Chemical potential of a hard sphere fluid confined in a hard sponge matrix as a function of fluid density. a) Cavity to fluid particle size ratio: 𝜎𝜎 0 𝜎𝜎 ⁄ = 5 ; Cavity center density: 𝜌𝜌 0 𝜎𝜎 3 = 0.020685, 0.012855 ; b) Cavity to fluid particle size ratio: 𝜎𝜎 0 𝜎𝜎 ⁄ = 7 ; Cavity center density: 𝜌𝜌 0 𝜎𝜎 3 = 0.006128, 0.003809; c) Cavity to fluid particle size ratio: 𝜎𝜎 0 𝜎𝜎 ⁄ = 10; Cavity center density: 𝜌𝜌 0 𝜎𝜎 3 = 0.001816, 0.001129 ; Red lines: Morphological thermodynamics combined with SPT; Blue dashed lines: SPT2b1 theory [16]; Black squares: µVT-ensemble Monet-Carlo simulations (results from Fig. 4 of [16]).

Figure 10 .

 10 Figure 10. Comparison of the results for the chemical potential of a confined fluid given by the full morphological thermodynamics (full orange lines) and by its simplified version without the contribution of curvature terms (dashed green lines). a) a HS fluid in a HS matrix with size ratio 𝜎𝜎 0 𝜎𝜎 ⁄ = 2 and matrix density 𝜌𝜌 0 𝜎𝜎 3 = 0.05 ; b) a HS fluid in a HS matrix with size ratio 𝜎𝜎 0 𝜎𝜎 ⁄ = 3 and matrix density 𝜌𝜌 0 𝜎𝜎 3 = 0.0149208, 0.019397; c) a HS fluid in an overlapping HS matrix with size ratio 𝜎𝜎 0 𝜎𝜎 ⁄ = 1 and matrix density 𝜌𝜌 0 𝜎𝜎 3 = 0.1590.

  A6) and (A7) lead immediately to the following nonzero adsorption,

  the confined and bulk fluids are given respectively by, obtains the following result for the disjoining pressure,𝛱𝛱 = 𝑝𝑝 -𝑝𝑝 bulk = k B 𝑇𝑇 𝑒𝑒 𝛽𝛽𝛽𝛽 𝛬𝛬 3 �𝜙𝜙 0 HS -1� . (A11)Disjoining pressure was discovered by Derjaguin in 1930's for a fluid confined between two closely approached flat solid surfaces[START_REF] Derjaguin | Direct measurement of molecular forces[END_REF][START_REF] Derjaguin | Surface forces[END_REF]. To the best of our knowledge, nonzero disjoining pressure has never been found for a fluid confined in a fluid confined in a porous matrix. The integral surface tension is given by, that this integral surface tension and the adsorption given in eq.(8) satisfies the following generalized Gibbs adsorption equation[START_REF] Dong | Thermodynamics of interface extended to nanoscales by introducing integral and differential surface tensions[END_REF][START_REF] Dong | Nanoscale thermodynamics needs the concept of a disjoining chemical potential[END_REF], i.e.,

Acknowledgement

This work is supported by the National Natural Science Foundation of China (Nos. 22008063 and 22178072), the Guangxi Science and Technology Base and Talent Special Project (No. AD21220017). CZQ thanks the National Natural Science Foundation of China for the financial support (project No. 22008063). HRJ is grateful for the Eiffel scholarship of the French government and the CSC scholarship of the Chinese government.

Appendix A. An ideal gas confined in a hard sphere matrix

In this appendix, we show how surface thermodynamics can be applied for an ideal gas adsorbed in a hard sphere matrix and how the surface tension can be defined in this case. As in the main text, we consider a grand canonic ensemble. The fluid-matrix interaction is given by,

where 𝒓𝒓 𝑖𝑖 is the position vector of ith fluid particle and 𝒒𝒒 𝑗𝑗 the position vector of jth matrix particle.

The interaction potential between a fluid particle and a matrix particle is given by,

The partition is given by, Eq.(A4) shows that the grand partition function does not depend on the configuration of the porous matrices. So, the disorder of the matrix configurations does not have any influence on the thermodynamics of this system. Moreover, the grand potential does not depend on the surface area of the matrix particle. Thus, the immediate consequence of this is that the differential surface tension is zero, i.e.,

At the first glance, this indicates that surface thermodynamics does not apply for such a system.

However, it is also straightforward to show that the adsorption in this hard sphere matrix is not zero.

The number of the confined fluid inside the matrix is given by,

However, the number of a bulk ideal gas occupying the same volume and having the same temperature and the same chemical potential is given by, adsorption allows for introducing properly the porosity, 𝜙𝜙 0 , into the adsorption isotherms, i.e., the chemical potential as a function of the density of the confined fluid. From eq.(A7), one has the following equation of state for the bulk ideal gas,

From eqs.(A6) -(A8), we obtain,

where 𝜌𝜌 bulk = 𝑁𝑁 bulk 𝑉𝑉 -1 , 𝜌𝜌 cf = 𝑁𝑁 cf 𝑉𝑉 -1 . Substituting eq.(A15) into eq.(A14), we obtain,

When going to the last equality of eq.(A16), eq.( A6) is used. This achieves the proof that the matrix porosity enters into the isotherm of adsorption, i.e., the relation between the chemical potential and the density of the confined fluid, through the adsorption due to the disjoining pressure.

Appendix B. Simulation method and computational conditions

In the main text, accompanying the presentation of the results given by morphological thermodynamics, Monte-Carlo simulation results are also presented to assess the accuracy of theoretical approach in each case. Except those results with their references being indicated, all the other ones are obtained from our own simulations. In this appendix, we summarize the simulation methods and the computation conditions for each porous medium.

Fluid in a slit pore

The Monte-Carlo simulation results presented in Figs. 2 and3 are obtained by our own simulation in a NVT-ensemble. The chemical potential is calculated by using the test particle method based on Widom's potential distribution theorem [START_REF] Widom | Some topics in the theory of fluids[END_REF]. The two pore walls are respectively placed at

The pore wall has a square shape of the size 10𝜎𝜎 × 10𝜎𝜎 . The periodic boundary condition is used in the two directions parallel to the pore walls. For each simulation run, 2 × 10 5

Monte-Carlo cycles are first performed to prepare the system to equilibrium, then 10 6 MC cycles are performed to calculate the chemical potential.

Fluid in various porous matrices

For the ordered porous matrix (results shown in Fig. 5), we studied only the case with matrix particles placed on a simple cubic lattice. A single matrix particle is placed at the center of the cubic simulation box. The periodic boundary condition is applied in three directions to generate a simple cubic lattice. The simulation box has respectively the size of 20𝜎𝜎 × 20𝜎𝜎 × 20𝜎𝜎 and 36𝜎𝜎 × 36𝜎𝜎 × 36𝜎𝜎 for the cases of 𝜎𝜎 0 = 5𝜎𝜎 and 𝜎𝜎 0 = 10 𝜎𝜎.

For the disordered HS matrix (results shown in Fig. 6), a matrix sample is first generated with For the overlapping hard sphere matrix (results shown in Fig. 7a), 159 matrix particles with the same diameter as the fluid particle, i.e. 𝜎𝜎 0 = 𝜎𝜎, are placed totally randomly in the simulation box of the size, 10𝜎𝜎 × 10𝜎𝜎 × 10𝜎𝜎. The average over matrix configurations is realized with 10


Effects of Shock-Induced Separation on Boundary Layer Transitional Mechanisms

Nikhil Mahalingesh * , Sébastien Piponniau † and Pierre Dupont ‡ Aix Marseille Univ, CNRS, IUSTI, Marseille, France

The transitional mechanisms of a natural boundary layer as well as one interacting with shock-waves are investigated at supersonic speeds. Experiments are performed over a flat plate and a 6°compression ramp at a Mach number of 1.65, using Pitot probes and hot-wire anemometry. Comparisons are drawn between the two configurations focusing on the evolution of unstable modes of the boundary layer. Measurements are made inside the laminar boundary layer to study the primary and secondary instabilities of the natural transition process. While the flat plate exhibits the classical mechanisms of transition for a zero pressure gradient boundary layer, the compression ramp shows an altered transition process where the secondary instabilities play a less important role.

I. Nomenclature

II. Introduction

The transitional mechanisms of laminar boundary layers have been an important field of study for many decades in the discipline of aerodynamics [START_REF] Schubauer | Laminar boundary-layer oscillations and stability of laminar flow[END_REF]. Understanding this flow phenomenon is of great importance due to the significant advantage of lower skin-friction drag associated with laminar boundary layers, especially at supersonic and hypersonic speeds. Further, wall-heat fluxes from a transitioning boundary layer at very high speeds may be larger than even those encountered with turbulent boundary layers [START_REF] Benay | Shock wave/transitional boundary-layer interactions in hypersonic flow[END_REF]. Consequently, laminar-to-turbulent transition has been extensively studied from theoretical [START_REF] Mack | Boundary-layer linear stability theory[END_REF][START_REF] Lees | Investigation of the stability of the laminar boundary layer in a compressible fluid[END_REF], experimental [START_REF] Laufer | Stability and transition of a supersonic laminar boundary layer on an insulated flat plate[END_REF][START_REF] Kendall | Wind tunnel experiments relating to supersonic and hypersonic boundary-layer transition[END_REF] and numerical [START_REF] Sandham | Direct simulation of breakdown to turbulence following oblique instability waves in a supersonic boundary layer[END_REF] points of view, with the first studies as far back as the early 20th century [START_REF] Tollmien | Über die entstehung der turbulenz. nachr. ges. wiss. göttingen 21-24[END_REF]. The classical self-excited Tollmien-Schlichting (TS) waves of incompressible boundary layers become three-dimensional at supersonic speeds (i.e. oriented with an angle with respect to the stream-wise direction) and hence are referred to as oblique modes [START_REF] Fasel | Direct numerical simulation of transition in supersonic boundary layers: oblique breakdown[END_REF]. At low supersonic speeds, the preferred scenario for transition involves the breakdown of stream-wise streak-like structures created by the non-linear interaction of oblique modes [START_REF] Chang | Oblique-mode breakdown and secondary instability in supersonic boundary layers[END_REF]. It was found that pairs of oppositely oriented oblique waves generate rapidly growing stream-wise structures [START_REF] Shadloo | Laminar-turbulent transition in supersonic boundary layers with surface heat transfer: a numerical study[END_REF]. These highly unstable streaks breakdown locally to create turbulent spots, leading to turbulence [START_REF] Mayer | Direct numerical simulation of complete transition to turbulence via oblique breakdown at Mach 3[END_REF]. At hypersonic speeds, "Mack" modes are the most dominant characterized by very high frequencies [START_REF] Mack | Boundary-layer linear stability theory[END_REF]. This transition scenario is quite different from the breakdown process followed by the growth of oblique modes and streaks. Although extensive research has been carried out on both the oblique modes and the streak-like structures, including characterizing their spatial and temporal dynamics, their effects and influence on the length and time-scales of various other flow phenomena involving laminar boundary layers are still not well understood.

One such flow phenomenon is the interaction of boundary layers with shock-waves. Shock-wave boundary layer interaction (SBLI) is a classical flow phenomenon of high speed aerodynamics. SBLIs are encountered in various applications such as engine inlets, transonic wings, rocket nozzles, high-speed compressors, etc. and hence, this topic has received a lot of attention over the past 70 years [START_REF] Délery | Shock-wave boundary layer interactions[END_REF][START_REF] Dolling | Fifty years of shock-wave/boundary-layer interaction research: what next?[END_REF]. When the shock strengths are high, the interaction tends to separate the boundary layer and a characteristic length of separation is obtained [START_REF] Dupont | Space and time organization in a shock-induced separated boundary layer[END_REF]. This length of interaction was found to scale non-linearly with shock strength, Reynolds and Mach number, with further dependence on the type of flow geometry such as oblique shock reflections and compression ramps [START_REF] Souverein | A scaling analysis for turbulent shock-wave/boundary-layer interactions[END_REF]. Further, such interactions are typically associated with intense low-frequency oscillations of shock-waves and large-scale "breathing" of the separation bubble [START_REF] Erengil | Unsteady wave structure near separation in a Mach 5 compression rampinteraction[END_REF][START_REF] Piponniau | A simple model for low-frequency unsteadiness in shock-induced separation[END_REF]. The mechanism of this unsteadiness has been the subject of much debate, where the scientific community was divided on the origins of this low-frequency dynamics, between the upstream boundary layer and downstream separated region [START_REF] Clemens | Low-frequency unsteadiness of shock wave/turbulent boundary layer interactions[END_REF][START_REF] Touber | Large-eddy simulation of low-frequency unsteadiness in a turbulent shock-induced separation bubble[END_REF]. It has to be noted that most of these studies were extensively focused on turbulent SBLIs and relatively less attention has been given to the interaction of laminar boundary layers with shock-waves [START_REF] Doerffer | Unsteady effects of shock wave induced separation[END_REF][START_REF] Babinsky | Shock wave-boundary-layer interactions[END_REF].

Recent interest in transitional interactions sought to characterize and establish an experimental and numerical database for various flow parameters and geometries [START_REF] Doerffer | Transition Location Effect on Shock Wave Boundary Layer Interaction: Experimental and Numerical Findings from the TFAST Project[END_REF]. These studies found a strong correlation between the location of boundary layer transition and the size of the separated region [START_REF] Giepman | A parametric study of laminar and transitional oblique shock wave reflections[END_REF]. Additionally, a discrepancy was found in the frequency of large-scale unsteadiness between experimental and numerical studies [START_REF] Diop | On the length and time scales of a laminar Shock Wave Boundary Layer Interaction[END_REF][START_REF] Larchevêque | Low-and Medium-Frequency Unsteadinesses in a Transitional Shock-Boundary Reflection with Separation[END_REF]. Detailed investigations into the dynamics of transitional SBLIs showed that low-frequency unsteadiness may be directly influenced by the breakdown of coherent structures at reattachment [START_REF] Sansica | Instability and low-frequency unsteadiness in a shock-induced laminar separation bubble[END_REF]. Further, Large Eddy Simulations (LES) of transitional SBLIs found that the amplitude of low-frequency unsteadiness varied significantly with the size of the separated region, which in turn was dependent on the amplitude of initial disturbances in the incoming boundary layer [START_REF] Larchevêque | Low-and Medium-Frequency Unsteadinesses in a Transitional Shock-Boundary Reflection with Separation[END_REF]. On the experimental side, high-resolution Laser Doppler Anemometry (LDA) and hot-wire measurements showed that specific frequency scales are exponentially amplified over the interaction, reaching a maximum at the location of transition, suggesting that they might be related to boundary layer transitional mechanisms [START_REF] Diop | High resolution LDA measurements in transitional oblique shock wave boundary layer interaction[END_REF]. However, these experiments could not characterize the exact mechanism of transition over the SBLI and required further studies to complete the picture.

The current paper hopes to advance the collective knowledge on the transitional mechanisms of supersonic laminar boundary layers. In particular, the paper focuses on the differences between the natural transition process over a flat plate and the accelerated transitional mechanisms over a compression ramp. The following section details the experiments that used high-resolution hot-wire anemometry as the primary diagnostic tool. The results of the natural boundary layer are discussed first and followed by the compression ramp experiments. The current paper focuses on the mechanism of transition and the associated differences between the two geometries.

III. Experimental Methodology

The experimental facilities of the supersonic wind tunnel at the IUSTI laboratory of Aix-Marseille University are used for performing the measurements. It is a closed-circuit continuous-operation wind tunnel with hypo-turbulent free-stream conditions. Experiments can be conducted up to 4 hours with nearly constant operating conditions. A classical converging-diverging nozzle accelerates the flow to a Mach number of 1.65 and a free-stream velocity of 464 m/s. The total pressure of the free-stream was set to 0.8 atm, with a total temperature of 305 K, resulting in a unit Reynolds number of 11.01 million 𝑚 -1 . Turbulence intensity of 0.05% in terms of mass-flux fluctuations, 0.03% in terms of velocity fluctuations and 0.11% in terms of pressure fluctuations, was measured using hot-wire anemometry at the exit of the nozzle. Such low values of turbulence intensities are necessary to avoid bypass transition of the laminar boundary layers [START_REF] Van Ingen | The eN method for transition prediction[END_REF]. Downstream of the nozzle, the test section has a height of 105 mm and span-wise width of 170 mm.

The two geometrical configurations of the flat plate and the compression ramp are similar in design and construction. Both have a sharp leading edge and a thickness of 5 mm, placed at a height of 25 mm from the floor of the test section. The picture and schematic of the flat plate geometry is shown in Fig. 1. The compression ramp has a 6°deflection at 115 mm from the leading edge as shown in Fig. 2. A supersonic Pitot probe, with a flat tip of 0.3 mm in height with a pressure tap of 0.15 mm, is used for measurements of mean flow quantities. The probe is traversed every 0.05 mm in the wall-normal direction to obtain boundary layer profiles, and longitudinal exploration over the compression ramp are made at a height of 5 mm from the wall. Isentropic compressible flow equations and normal shock relations are applied to the Pitot measurements and steady-state flow quantities are obtained.

A hot-wire anemometer setup with a tungsten and platinum wire of 2.5 𝜇m diameter is used with Dantec Dynamics' Streamline amplifier. The anemometer is operated in constant temperature mode with a symmetrical bridge. The overheat ratio is set to 0.8 and the effective frequency bandwidth is approximately 150 kHz (depending on the wire characteristics). Measurements are sampled at 400 kHz for 10 seconds using the National Instruments NI6133 acquisition system with 14 bits of resolution.

Boundary layer hot-wire probes are used for measurements at the wall (Fig. 3b). A number of studies have shown that such probes can be placed inside laminar boundary layers to obtain reliable measurements , and even at supersonic flow conditions [START_REF] Schubauer | Laminar boundary-layer oscillations and stability of laminar flow[END_REF][START_REF] Laufer | Stability and transition of a supersonic laminar boundary layer on an insulated flat plate[END_REF][START_REF] Kendall | Wind tunnel experiments relating to supersonic and hypersonic boundary-layer transition[END_REF]. A small resin at the corner of the prongs is applied so that they can be placed flush on the wall while resting on the resin. This is to ensure the vibrations and oscillations of the probe and its support system are damped. Figure 4 illustrates the measurement points made inside the laminar boundary layer, where the wire was at a constant wall-normal distance of 0.4 mm from the wall. From the illustration, it can be seen that the boundary layer thickness is lower than 0.5 mm over the first half of the flat plate. As a result, boundary layer measurements could not be performed here as the probe was outside or close to the edge of the boundary layer. However, hot-wire probes cannot be placed inside separated regions with reverse flow as the wire cannot distinguish between positive and negative velocities, and furthermore, the wire is easily broken. Hence, classical flat hot-wire probes are traversed over the SBLI (Fig. 3a). The probe in the external flow is sensitive to fluctuations radiated by the boundary layer and associated flow structures of the interaction [START_REF] Smits | Turbulent shear layers in supersonic flow[END_REF][START_REF] Diop | Transition mechanism in a shock wave boundary layer interaction[END_REF]. Additionally, fluctuations from the flow on the side walls and the roof of the test section also influence the hot-wire measurements, but the amplitude of these "artifacts" decay rapidly moving away from their respective sources. Given the large size of the test section (relative to the boundary layer thickness on the side walls and roof), a probe placed in the mid-span and close to the laminar boundary layer has minimal influence from these artifacts. The instantaneous voltage of the hot-wire is influenced by fluctuations of total temperature, mass-flux and Mach number of the flow. For Mach numbers higher than 1.4 and high overheat ratios, non-dimensional fluctuations of mass-flux can be directly obtained from the hot-wire measurements using the King's law coefficient [START_REF] Morkovin | Fluctuations and hot-wire anemometry in compressible flows[END_REF][START_REF] Dupont | Etude expérimentale des champs turbulents dans une couche limite supersonique fortement chauffée[END_REF]. Information regarding the time-scales of the boundary layer and the interaction can thus be obtained from the measurements of the hot-wire anemometer. Fourier analyses are applied to the measured time-signals and estimates of Power Spectral Density (PSD) are determined using Welch's method. Frequencies up to 200 kHz are resolved and fully converged spectral statistics are obtained with 2 16 discrete Fourier points.

IV. Results

The first part of the results section focuses on the measurements over the flat plate and the second part details the results of the compression ramp.

A. Natural Boundary Layer

Measurements of the boundary layer profile at different stream-wise locations on the flat plate are shown in Fig. 5. The profiles are represented in similarity coordinates, where non-dimensional wall-normal distance is shown on the horizontal axis, scaled with respect to the boundary layer thickness (𝛿 99 ) and non-dimensional velocity is shown on the vertical axis, scaled with respect to the external velocity outside the boundary layer (𝑢 𝑒 ). The measurements are compared with the compressible Blasius and turbulent boundary layer profiles, which are obtained numerically for the Mach and Reynolds number of the current experiments. The stream-wise coordinate is expressed in terms of Reynolds number based on the compressible boundary layer thickness (𝑅𝑒 𝛿 * ), which considers compressibility and Mach effects. Such a scaling parameter has been previously used to study transitional laminar boundary layers [START_REF] Mack | Boundary-layer linear stability theory[END_REF][START_REF] Laufer | Stability and transition of a supersonic laminar boundary layer on an insulated flat plate[END_REF]. At low Reynolds numbers, the measured profiles exhibit a large linear behavior near the wall, characteristic of laminar boundary layers. Comparisons with theoretical Blasius profiles yield good agreements. However, at higher Reynolds numbers, the measurements deviate significantly from the Blasius profile and are closer to the turbulent profile. To understand the state of the boundary layer in more detail, table 1 shows the incompressible shape factor (𝐻 𝑖 ) of the profiles, and comparisons are drawn with theoretical shape factors for laminar and turbulent boundary layers. The shape factor at lower Reynolds numbers are close to the Blasius profile, but it reduces significantly at higher Reynolds numbers, but does not reach the shape factor of a fully turbulent profile. This suggests that the boundary layer is still in its transitional stages. It is to be noted that the lack of measurement points close to the wall has an influence on the calculated shape factor. Nevertheless, these results suggest that the beginning of transition significantly changes the boundary layer profile, making it "fuller".

Moving to the time-scales of the boundary layer, Fig. 6 shows the contours of PSD at different stream-wise locations on the flat plate. The PSD is represented in pre-multiplied and normalized format to highlight the dominant frequencies at each stream-wise location. Non-dimensional angular frequency (𝜔 * ) is shown on the vertical axis in logarithmic scale. Frequency is normalized according to Eq. 1, where 𝑈 ∞ is the free-stream velocity and 𝑅𝑒 𝑢 is the unit Reynolds number of the free-stream. This normalization is typically used to identify boundary layer modes [START_REF] Mack | Boundary-layer linear stability theory[END_REF].

𝜔 * = 2𝜋 𝑓 𝑈 ∞ 𝑅𝑒 𝑢 (1)
The contour map on the bottom corresponds to the measurements made inside the boundary layer (Fig. 4) and the one on the top corresponds to measurements made outside the boundary layer at a constant height of 3.85 mm from the wall. Inside the boundary layer, two time-scales are most dominant. At lower Reynolds numbers, the PSD is concentrated at high frequencies (HF) around 2 × 10 -5 ≤ 𝜔 * ≤ 4 × 10 -5 . These scales are typically associated with the compressible version of Tollmien-Schlichting (TS) waves, which are 3D and oblique in nature for supersonic boundary layers [START_REF] Laufer | Stability and transition of a supersonic laminar boundary layer on an insulated flat plate[END_REF][START_REF] Kendall | Wind tunnel experiments relating to supersonic and hypersonic boundary-layer transition[END_REF]. According to stability theory, these scales are the natural eigen modes of the boundary layer, arising due to the linear growth of small disturbances [START_REF] Mack | Boundary-layer linear stability theory[END_REF]. These oblique modes have been documented extensively by numerical and experimental studies [START_REF] Fasel | Direct numerical simulation of transition in supersonic boundary layers: oblique breakdown[END_REF][START_REF] Chang | Oblique-mode breakdown and secondary instability in supersonic boundary layers[END_REF][START_REF] Kosinov | Experiments on the stability of supersonic laminar boundary layers[END_REF][START_REF] Demetriades | Growth of disturbances in a laminar boundary layer at Mach 3[END_REF].

Further downstream on the plate, a medium frequency (MF) scale around 4 × 10 -6 ≤ 𝜔 * ≤ 2 × 10 -5 , appears alongside the HF scale. These scales are characteristic of stream-wise streak-like structures associated with the non-linear stage of the transition process [START_REF] Sandham | Direct simulation of breakdown to turbulence following oblique instability waves in a supersonic boundary layer[END_REF]. These scales are created when the amplitude of the oblique modes reach a certain threshold, triggering the formation of secondary scales [START_REF] Teramoto | Large-eddy simulation of transitional boundary layer with impinging shock wave[END_REF][START_REF] Giepman | High-resolution PIV measurements of a transitional shock wave-boundary layer interaction[END_REF]. The presence of these scales indicate the beginning of the non-linear regime of the transition process, which is soon followed by breakdown and turbulence [START_REF] Mayer | Direct numerical simulation of complete transition to turbulence via oblique breakdown at Mach 3[END_REF].

Outside the boundary layer, the hot-wire probe picks up radiated fluctuations from the boundary layer. It can be seen that for the beginning of the flat plate, the hot-wire records mostly the noise of the wind tunnel. No coherent frequency scale can be observed and relatively low energy is spread over a broad range of high frequencies. However, near the end of the flat plate, a clear time-scale can be observed corresponding to the MF scale seen inside the boundary layer. It seems that the HF scale observed inside the boundary layer is not seen outside, only the MF scales are observed at roughly the same location as seen inside. These results suggest that the HF oblique modes do not radiate fluctuations outside the boundary layer unlike the MF streak-like structures.

To shed light on the growth rates of these scales, Fig. 7 shows the non-dimensional amplitude of mass-flux fluctuations along the flat plate, from the hot-wire measurements. The PSD is integrated over a short band of frequencies associated with the two time-scales of interest to obtain band-passed standard deviation (𝜎 𝑥 ), according to Eq. 2, where 𝑒 is the mean voltage of the hot-wire probe. The band-passed standard deviation at each location is compared with the noise of the wind tunnel (𝜎 0 ) over the same band of frequencies to obtain non-dimensional amplitude (𝑛), as shown in Eq. 3). The noise of the wind tunnel is measured at the exit of the converging-diverging nozzle and serves as the input to the laminar boundary layer. This normalization helps distinguish the presence of coherent time-scales associated with the boundary layer from the wind tunnel noise and the electronic noise of the hot-wire system. Measurements outside the boundary layer show that the HF scales have a nearly constant amplitude over the flat plate. As there is no significant amplification and no distinct signature of a coherent scale from the frequency spectrum (Fig. 6), this corresponds to the free-stream noise of the wind tunnel. The amplitude suggests that the noise has a higher level compared to the free-stream, but constant nonetheless. The MF scales also exhibit a similar constant amplitude up to 𝑅𝑒 𝛿 * ≈ 3300, after which they show a steady exponential growth. Measurements made inside the boundary layer show that both time-scales are continuously amplified. The HF scale starts with a higher relative amplitude compared to the MF scale. However, the growth rate of the MF scale is higher than the HF scale. Once the two curves cross, the HF scale grows with the same rate of the MF scale, reaching a relative amplitude of more than 2 decades over the initial value. For context, the broadband turbulence intensity (in terms of mass-flux fluctuations) at the end of the flat plate is 0.21% outside the boundary layer and 5.65% inside the boundary layer. It is also interesting to observe the "dissipation" of the amplitude measured outside the boundary layer relative to what was measured inside. The amplitude of the MF scales outside the boundary layer is approximately half of the value measured inside. It is to be noted that the stream-wise coordinate of the measurements made outside are projected to the wall based on characteristic angles and hence, there should not be an offset based on the location of the sources.

𝜎 𝑥 = 𝑓 2 ∫ 𝑓 1 𝑃𝑆𝐷 (𝑒) 2 𝑑𝑓 (2) 𝑛 = 𝜎 𝑥 𝜎 0 (3)
These results suggest that the boundary layer undergoes transition through the growth of two time-scales; the first step involves the growth of self-excited HF modes corresponding to the 3D oblique waves. While these modes dominate the frequency space inside the boundary layer, their presence cannot be detected outside the boundary layer. Once these modes reach a critical amplitude, they trigger the growth of secondary instabilities corresponding to MF stream-wise streak-like structures. The presence of these MF scales indicates the beginning of the non-linear regime of the transition process, resulting in modifications to the mean boundary layer profile, gradually making it more turbulent. Additionally, these MF scales have a higher growth rate compared to the HF scales and their presence can also be detected outside the boundary layer, albeit with lower amplitude.

Thus, a critical Reynolds number can be determined based on the location where the secondary instabilities are first observed outside the boundary layer. This value was found to be 𝑅𝑒 𝛿 * ≈ 3314 or 𝑅𝑒 𝑥 ≈ 1.4 × 10 6 (based on the stream-wise coordinate). This critical value denotes the beginning of the transition process. The location where a fully turbulent profile can be obtained will be further downstream of this critical Reynolds number.

For reference, the location of the corner of the compression ramp expressed in terms of this Reynolds number is 𝑅𝑒 𝛿 * = 3117. Thus, the location of the corner is in the part of the boundary layer where the profile is still fully laminar (Fig. 5) and the secondary instabilities are not "yet" dominant inside the boundary layer (Fig. 6). This ensures that the compression ramp interacts with a canonical laminar boundary layer and not one that has begun the transition process.

B. Shock-wave Boundary Layer Interaction (SBLI)

Compression Waves

S C R

Reattachment Shock

Laminar BL Turbulent BL

L int

Separated Shear Layer The interaction between the laminar boundary layer and the compression ramp is illustrated in Fig. 8, where the points 𝑆, 𝐶 and 𝑅 represent separation, corner of ramp and reattachment of the boundary layer respectively. Separation is characterized by weak compression waves and a stronger flow deflection at reattachment leads to a coalescing of compression waves into a shock-wave further away from the wall. The recirculating region is long and thin, shown by the gray region between the points 𝑆, 𝐶 and 𝑅. Such large aspect ratios of the separated region was also found by previous studies [START_REF] Giepman | A parametric study of laminar and transitional oblique shock wave reflections[END_REF][START_REF] Larchevêque | Low-and Medium-Frequency Unsteadinesses in a Transitional Shock-Boundary Reflection with Separation[END_REF][START_REF] Diop | High resolution LDA measurements in transitional oblique shock wave boundary layer interaction[END_REF]. Far downstream of reattachment, the boundary layer is expected to be fully turbulent. Revisiting the topic of critical Reynolds number from the previous section, given the large aspect ratio of the separated region, Reynolds number at separation would be a better parameter to estimate the state of the boundary layer upstream of the interaction. Reynolds number at separation of the boundary layer was found to be 𝑅𝑒 𝛿 * = 2888, which is clearly upstream of the critical Reynolds number estimated in the previous section.

Figure 9 shows the external longitudinal measurements from the Pitot probe, which are made a constant height of 5 mm above the wall. The stream-wise coordinate is expressed in non-dimensional form with 𝑋 * = 0 corresponding to the mean location of separation, 𝑆 and 𝑋 * = 1 corresponding to the location of the corner, 𝐶 (Fig. 8). The average location of separation is determined by identifying the peak in the gradient of longitudinal pressure. As the measurements are made at a height of 5 mm from the wall, 𝑋 * is determined using projected coordinates, along characteristics. The left vertical axis shows non-dimensional pressure coefficient (𝑐 𝑝 ) over the interaction, and a two-step pressure rise can be observed, typical of separated SBLIs. The overall pressure rise agrees quite well with the inviscid step (shown by the dashed-dotted line), with a slight underestimation from the experiments, as the total pressure loss across the reattachment shock cannot be accurately determined.

Fig. 9 Stream-wise evolution of non-dimensional pressure over the compression ramp.

The right vertical axis of Fig. 9 shows the pressure measurements expressed in coefficient of free-interaction (𝐹), which takes into account the effect of Mach and Reynolds number (through the skin-friction coefficient (𝑐 𝑓 ), according to Eq. 4) to provide universal values for separation and plateau pressure of SBLIs [START_REF] Chapman | Investigation of separated flows in supersonic and subsonic streams with emphasis on the effect of transition[END_REF]. It can be seen that the experimental measurements agree quite well with typical values of separation (𝐹 = 0.8, indicated by the dashed line) and plateau pressure (𝐹 = 1.5, indicated by the dotted line) for laminar separated SBLIs [START_REF] Babinsky | Shock wave-boundary-layer interactions[END_REF]. These results validate the canonical nature of the current SBLI experiments. The length of interaction, 𝐿 𝑖𝑛𝑡 , defined as the stream-wise distance between the mean location of separation, 𝑆 and the location of the corner, 𝐶 (Fig. 8), was measured to be 62𝛿 * 𝑠 , where the displacement thickness at separation was taken as reference. A more detailed analyses of the interaction lengths for different Reynolds numbers as well as comparisons with oblique shock reflections have been published by the authors in the proceedings of the Turbulent Shear Flow Phenomenon (TSFP) conference [START_REF] Mahalingesh | Transitional shock-wave boundary layer interaction over a compression ramp[END_REF].

𝑝 -𝑝 0 𝑞 0 = 𝐹 (𝑋 * ) √ 2𝑐 𝑓 𝑀 2 -1 (1/2) (4)
The frequency space of the hot-wire measurements provide a more detailed outlook on the time-scales of the interaction. External longitudinal explorations by the hot-wire probe are made at a height of 23 mm above the wall, which was chosen to avoid interactions between the probe support system and the ramp geometry. Figure 10 shows the evolution of normalized pre-multiplied PSD over the compression ramp, similar to Fig. 6. Frequency is expressed in non-dimensional Strouhal number, 𝑆𝑡 𝐿 = 𝑓 𝐿 𝑖𝑛𝑡 /𝑈 ∞ (where 𝐿 𝑖𝑛𝑡 is the length of interaction) on the left vertical axis in logarithmic scale. The right vertical axis shows non-dimensional angular frequency expressed according to boundary layer scaling (Eq. 1). The conversion factor between the two non-dimensional frequency scales is 𝑆𝑡 𝐿 /𝜔 * = 3.5 × 10 -5 . The horizontal axis represents the projected non-dimensional stream-wise coordinate (𝑋 * ) similar to Fig. 9.

Fig. 10 Longitudinal evolution of normalized pre-multiplied power spectral density over the compression ramp.

The separation bubble (0 ≤ 𝑋 * ≤ 1) seems to be dominated by a broad range of high frequency scales centered around 𝑆𝑡 𝐿 ≈ 1 (2 × 10 -5 ≤ 𝜔 * ≤ 1 × 10 -4 ). These scales do not seem to be related to the interaction and might be the measurement of high frequency noise of the wind tunnel. Similar noise was measured over the flat plate outside the boundary layer (Fig. 6). This will be further addressed in the following paragraphs.

Reattachment shock can be clearly identified downstream of the corner at 𝑋 * = 1.4. The shock signal seems to be characterized by a wide range of low frequencies centered around 𝑆𝑡 𝐿 ≈ 0.09 (indicated by the dashed white line). A signature of this low frequency unsteadiness can also be observed just downstream of reattachment. Surprisingly, this LF scale does not seem to be dominant at the separation of the boundary layer, suggesting that the amplitude of unsteadiness might be quite low. Similar low frequency unsteadiness was found by experimental investigations of transitional oblique shock reflections at similar Mach and Reynolds number conditions [START_REF] Diop | High resolution LDA measurements in transitional oblique shock wave boundary layer interaction[END_REF][START_REF] Bonne | Analysis of the two-dimensional dynamics of a Mach 1.6 shock wave/transitional boundary layer interaction using a RANS based resolvent approach[END_REF]. However, numerical simulations found lower Strouhal numbers of 𝑆𝑡 𝐿 ≈ 0.04 [START_REF] Larchevêque | Low-and Medium-Frequency Unsteadinesses in a Transitional Shock-Boundary Reflection with Separation[END_REF][START_REF] Sansica | Instability and low-frequency unsteadiness in a shock-induced laminar separation bubble[END_REF]. At higher Mach numbers, experiments on transitional SBLIs over compression ramps found low-frequency unsteadiness at 𝑆𝑡 𝐿 ≈ 0.025 [START_REF] Threadgill | Transitional shock boundary layer interactions on a compression ramp at mach 4[END_REF]. This discrepancy remains an open question in the field of transitional SBLIs.

Additionally, an intermediate scale can be seen around 𝑆𝑡 𝐿 ≈ 0.15 (𝜔 * ≈ 6 × 10 -6 ) downstream of reattachment, annotated by the dotted white line. The frequency of these scales are approximately in the same range of frequencies associated with streak-like structures observed over the natural boundary layer (Fig. 6). Recent experimental and numerical studies on transitional compression ramps at low hypersonic speeds, found similar streak-like structures in the boundary layer at reattachment [START_REF] Lugrin | Transitional shockwave/boundary layer interaction experiments in the R2Ch blowdown wind tunnel[END_REF]. These time-scales were found to be essential for the breakdown process leading to turbulence [START_REF] Lugrin | Transition scenario in hypersonic axisymmetrical compression ramp flow[END_REF]. Although the Mach number of these investigations are higher than the current experiments, the transitional mechanisms involved similar amplification of oblique modes and stream-wise streak-like structures. Although these scales are observed after reattachment, no concentration of PSD can be seen at very high frequencies, suggesting the absence of the typical cascade to turbulence. The results from Fig. 11 and 10 suggest that the boundary layer is not fully turbulent at reattachment, but still undergoing transition. Fig. 11 shows the longitudinal evolution of non-dimensional mass-flux fluctuations over the compression ramp. The amplitude is again expressed according to Eq. 3 and calculated for three time-scales: low frequency (LF) corresponding to the unsteadiness found at reattachment, medium frequency (MF) found downstream of reattachment and high frequency (HF) found over the separation bubble. The three time-scales exhibit a general trend of an overshoot due to the reattachment shock, with subtle differences upstream and downstream. LF seems to be the most dominant at separation (𝑋 * = 0) and reaches a lower plateau value till reattachment. MF remains fairly constant up to reattachment with a small increase over the second half of the bubble. Downstream of the reattachment shock, both the LF and MF have a higher amplitude relative to upstream of the interaction, suggesting a localized amplification due to the shock-wave. Their amplitude gradually decreases moving downstream and reaching the initial value found upstream of the interaction. The HF scale shows a decreasing amplitude over the separation bubble and continues to decrease after reattachment, with a relatively lower overshoot (compared to the LF and MF scales), thus confirming that the HF scales observed over the interaction are mostly related to the noise of the wind tunnel.

The overarching theme of these results is the absence of steady exponential growth either before or after reattachment, indicating the absence of canonical transitional mechanisms as seen in Fig. 7. But perhaps the boundary layer undergoes transition in another scenario. As the boundary layer has experienced a strong adverse pressure gradient due to the interaction, the mean velocity profile is expected to be severely decelerated downstream of reattachment and not to follow the shape of canonical boundary layer profiles. Such a velocity profile is not suited for the steady growth of conventional modes or structures associated with natural transition.

In order to shed more light on the state of the boundary layer after reattachment, hot-wire measurements are made inside the boundary layer at various locations downstream of the interaction using the same boundary layer probe shown in Fig. 3b. Figure 12 shows PSD of the hot-wire measurements made inside the boundary layer after reattachment. The PSD is represented in pre-multiplied form without normalization to highlight the evolution of relative amplitudes between different stream-wise locations. Frequency is shown in logarithmic scale on the horizontal axis, with non-dimensional angular frequency, 𝜔 * (boundary layer normalization) on the bottom and Strouhal number based on the length of interaction, 𝑆𝑡 𝐿 (SBLI normalization) on the top. It is important to note that the measurements made outside the boundary layer have four orders of magnitude lower amplitude compared to the measurements made inside the boundary layer. Further, these results highlight the strong contrast between the frequency scales observed inside and outside the boundary layer. Outside the boundary layer, most of the energy is concentrated at relatively lower frequencies. Evidence of MF scales are found, represented by the dashed vertical line, similar to Fig. 6 and10. Additionally, evidence of remnant LF unsteadiness can also be found around 𝑆𝑡 𝐿 ≈ 0.09, annotated by the dashed-dotted vertical line. Both of these scales decay quite rapidly moving further away from the interaction. No coherent high frequency scales can be identified. The last measurement point, which is nearly four length scales downstream of the interaction, measures very low amplitudes, close to the free-stream noise of the wind tunnel.

Inside the boundary layer, most of the energy is concentrated at higher frequencies. Evidence of boundary layer oblique modes (2 × 10 -5 ≤ 𝜔 * ≤ 4 × 10 -5 ) can be seen, in the same range as observed from the natural boundary layer measurements (Fig. 6). Also, a very high frequency scale (VHF) can be observed between 5 × 10 -5 ≤ 𝜔 * ≤ 1 × 10 -4 , which was not found inside the natural boundary layer. The black vertical line at 𝜔 * = 1.8 × 10 -4 represents the theoretical cut-off of the hot-wire anemometer for reference. Moreover, both of these scales reduce in amplitude moving downstream, with the oblique modes decaying faster than the VHF scales. This suggests that the boundary layer is still recovering, up to three interaction lengths downstream of reattachment. On first look, it appears that the MF scales are not found inside the boundary layer. By limiting the vertical axis, MF scales can be identified in the same range of frequencies (𝜔 * ≈ 5 × 10 -6 ) found in the natural boundary layer (Fig. 13, annotated by the vertical dashed line). The HF and the VHF scale are also shown for reference using the dotted and dashed-dotted line respectively. Comparing Fig. 12 and 13, it can be seen that the amplitude of HF scales are more than one order of magnitude higher than the MF scales inside the reattached boundary layer. Also, the amplitude of these MF scales decay quite rapidly and cannot be observed at the next measurement point. These results suggest an alternate breakdown scenario of laminar boundary layers, bypassing the growth of secondary instabilities.

On the other hand, no evidence of LF unsteadiness can be observed inside the boundary layer (Fig. 13), possibly due to their very low amplitude with respect to the MF, which cannot be resolved by the hot-wire anemometer.

Figure 14 shows the evolution of broadband mass-flux fluctuations (RMS) after reattachment inside and outside the boundary layer. Both measurements tell the story of a relaxing and recovering boundary layer, with the RMS gradually decreasing over large distances. However, the level of RMS measured inside the boundary layer is greater than 1%. Further, the RMS decreases to approximately half of its original value over a distance of four interaction lengths. The measurements made outside the boundary layer are analogous to what was seen in Fig. 11, with very low amplitude fluctuations (< 0.1%), close to the free-stream noise of the wind tunnel.

It seems that unlike the transition process of the natural boundary layer, HF oblique modes are the most dominant time-scales inside the boundary layer during the end stages of transition. This is evident by their high amplitude inside the boundary layer downstream of reattachment. Even with such a high amplitude, their presence cannot be detected outside the boundary layer. MF scales are observed inside the boundary layer as well, but having an amplitude one order of magnitude lower than the HF scales. Their existence is apparent by their measurement outside the boundary layer. However, these MF scales are not dominant and decay quite rapidly moving downstream. The presence of new VHF scales imply that the boundary layer is undergoing transition through breakdown of oblique modes, bypassing the rapid growth of MF scales. This transition scenario is clearly different from the natural process over the flat plate. The current experiments detail the transitional mechanisms over the compression ramp for only one set of Mach and Reynolds numbers. Moreover, the results shown here correspond to a relatively small adverse pressure gradient (i.e. 6°c ompression ramp). The behavior of these modes when encountered with higher shock strengths cannot be extrapolated from the current results. Further experiments are required to understand the effects of these flow parameters on these boundary layer time-scales. Moreover, the influence of these transitional modes on the length and time-scales of the SBLI are not fully understood. Parametric studies are necessary to uncover these subtle interdependencies between flow conditions, boundary layer modes and the spatial-temporal organization of transitional SBLIs.

V. Conclusion

Experiments are performed to investigate the transitional mechanisms of supersonic boundary layers. Specifically, comparisons are drawn between the natural transition scenario over a simple flat plate and the altered transition mechanism through shock-wave boundary layer interaction over a compression ramp. The canonical nature of the laminar boundary layer and SBLI has been verified using Pitot probe measurements. Hot-wire anemometer is used to obtain temporal information, including measurements inside boundary layers and over separated shear layers of the SBLI.

It was found that the natural boundary layer undergoes transition in a procedural manner, involving the development of consecutively amplified modes. First the growth of self-excited oblique modes was observed, which in turn triggered the growth of secondary instabilities or stream-wise streak-like structures. Once a critical Reynolds number was reached, these secondary instabilities are amplified exponentially due to non-linear growth, subsequently modifying the mean velocity profile and having large amplitude fluctuations.

Conversely, the boundary layer adopts a different mechanism to undergo transition when interacting with shock-waves. While the same instabilities of the natural transition process are found in the SBLI, their procedural and ordered growth of disturbances are not found. Instead, the boundary layer seems to undergo breakdown of oblique modes, without the rapid growth of secondary instabilities.

These results highlight the subtle contrast in the transitional mechanisms of the laminar boundary layer when interacting with shock-waves. Further experiments and analyses are required to understand the interplay between the influence of flow conditions such as Mach number, Reynolds number and shock strength on the boundary layer transitional modes. Moreover, the consequent effect of these instabilities on the length and time scales of the SBLI require detailed parametric investigations.

  𝑐 𝑓 = Skin-friction coefficient 𝑐 𝑝 = Pressure coefficient 𝛿 99 = Boundary layer thickness [𝑚] 𝛿 * = Compressible displacement thickness [𝑚] 𝑒 = Mean voltage of the hot-wire probe [𝑉] 𝑓 = Frequency [𝐻𝑧] 𝐻 𝑖 = Incompressible shape factor 𝐿 𝑖𝑛𝑡 = Length of interaction [𝑚] 𝑀 = Mach number 𝑛 = Non-dimensional amplitude 𝜔 * = Non-dimensional angular frequency 𝑝 = Static pressure [𝑃𝑎] 𝑞 = Dynamic pressure [𝑃𝑎] 𝑅𝑒 𝛿 * = Reynolds number based on the compressible displacement thickness 𝑅𝑒 𝑢 = Unit Reynolds number [𝑚 -1 ] 𝑅𝑒 𝑥 = Reynolds number based on the stream-wise coordinate (𝜌𝑈) ′ /(𝜌𝑈) = Non-dimensional fluctuations of mass-flux 𝑆𝑡 𝐿 = Strouhal number based on the length of interaction 𝜎 𝑥 = Band-passed standard deviation at a location 𝑥 𝜎 0 = Band-passed standard deviation at the free-stream of the wind tunnel 𝑈 ∞ = Free-stream velocity [𝑚/𝑠] 𝑢 𝑒 = External velocity outside the boundary layer [𝑚/𝑠] 𝑋 * = Non-dimensional stream-wise coordinate

  (a) Picture of the flat plate. (b) Schematic of the flat plate.

Fig. 1

 1 Fig. 1 Flat plate geometrical configuration.

Fig. 2

 2 Fig. 2 Compression ramp geometrical configuration.

Fig. 3

 3 Fig. 3 Pictures of the two types of hot-wire prongs.

Fig. 4

 4 Fig. 4 Illustration of the hot-wire measurements made inside the boundary layer.

Fig. 5

 5 Fig. 5 Boundary Layer velocity profiles at different Reynolds numbers.

Fig. 6

 6 Fig. 6 Stream-wise evolution of normalized pre-multiplied power spectral density of the laminar boundary layer.

Fig. 7

 7 Fig. 7 Stream-wise evolution of fluctuations of the laminar boundary layer.

Fig. 8

 8 Fig. 8 Illustration of the SBLI over the compression ramp.

Fig. 11

 11 Fig. 11 Longitudinal evolution of non-dimensional fluctuations of mass-flux over the compression ramp.

Fig. 12

 12 Fig. 12 Pre-multiplied power spectral density downstream of the SBLI.

Fig. 13

 13 Fig. 13 Pre-multiplied power spectral density inside the boundary layer, downstream of the SBLI.

Fig. 14

 14 Fig. 14 Non-dimensional fluctuations of mass-flux, downstream of the SBLI.

Table 1 Incompressible shape factor of boundary layer profiles.

 1

	𝑅𝑒 𝛿 *	𝐻 𝑖
	Blasius	2.69
	2221	2.64
	2714	2.83
	3161	2.70
	3448	2.13
	3684	2.10
	Turbulent 1.42

Acknowledgments

This work is part of the European TEAMAero project (Towards Effective Flow Control and Mitigation of Shock Effects in Aeronautical Applications). The authors would like to acknowledge the support from the European Union's Horizon 2020 research and innovation programme under grant agreement No. EC grant 860909.