Manel Wannassi

Isabelle Raspo

Numerical study of non-isothermal adsorption of Naphthalene in supercritical CO 2 : behavior near critical point

Keywords: Supercritical fluids, Adsorption, Piston effect, Numerical analysis

In this study, adsorption in a model binary mixture is investigated near the critical point in a side-heated cavity. The diverging behavior of the equilibrium constant and the Piston effect are taken into account and their influence on the adsorption process is pointed to. The modeling is based on numerical integration of the differential equations, considering the Navier-Stokes equations coupled with the energy and mass diffusion balances. By means of this model, the temperature, density and adsorbed concentration profiles are drawn at different times. Some fundamental concepts about the system's response to the heating are illustrated.

The results reveal that the adsorption process is influenced by the combined effect of several parameters, such as the gravity and the proximity to the critical point. In particular, the adsorbed amount exhibits a reversed dependency on the wall heating very close to the critical point, which confirms the complexity of such a process in binary systems near critical conditions.

Highlights

• Non-isothermal adsorption in near-critical binary mixtures was investigated by numerical simulations.

• The adsorption behavior near solvent's critical point has been analyzed.

• The effect of divergent properties and the piston effect were highlighted.

• A strong dependence to temperature and pressure variations in the vicinity of the critical point was depicted.

A c c e p t e d M a n u s c r i p t

Introduction

The supercritical state was first reported in 1822 by Baron Gagniard de la Tour [START_REF] De La Tour | Exposé de quelques résultats obtenu par l'action combinée de la chaleur et de la compression sur certains liquides, tels que l'eau, l'alcool, l'éther sulfurique et l'essence de pétrole rectifiée[END_REF], but only one hundred years later, supercritical techniques have received increased attention and have been used in analytical and on an industrial scale. This state is achieved when the temperature and the pressure of a substance is set over their critical values. So the properties A c c e p t e d M a n u s c r i p t 3 of a supercritical fluid range between those of a liquid and a gas and the distinction between the liquid and the gas phases is not possible. Some of the properties of a supercritical fluid are more liquid-like, whereas others are more gas-like.

Moreover, very close to the critical point, some properties diverge and others tend to zero. In fact, a small raise in pressure remarkably increases the fluid density and this effect diminishes with increasing distance from the critical point. On the other hand, a supercritical fluid has a higher diffusion coefficient and lower viscosity and surface tension than a liquid solvent, which leads to a more favorable mass transfer. Supercritical fluids exhibit very interesting qualities with regard to their physicochemical properties as well as ecology and economy. They are used as an alternative to organic liquid solvents in several applications such as cleanings [START_REF] Aslanidou | A novel approach for textile cleaning based on supercritical CO2 and Pickering emulsions[END_REF][START_REF] Ventosa | Copper cleaning in supercritical CO(2) for the microprocessor interconnects[END_REF][START_REF] Ito | Performance of air filters cleaned by supercritical carbon dioxide[END_REF]. Adsorption technologies using supercritical fluids have been also focused due to their potential applications including analytical extractions, activated carbon regeneration and soil remediation. Several studies have investigated the supercritical adsorption characteristics of many systems [START_REF] Srinivasan | Supercritical Fluid Desorption from Activated Carbon[END_REF][START_REF] Madras | Supercritical Fluid Regeneration of Activated Carbon Loaded with Heavy Molecular Weight Organics[END_REF][START_REF] Macnaughton | Supercritical Adsorption and Desorption Behavior of DDT on Activated Carbon Using Carbon Dioxide[END_REF][START_REF] Tan | Supercritical Regeneration of Activated Carbon Loaded with Benzene and Toluene[END_REF][START_REF] Tan | Desorption of Ethyl Acetate from Activated Carbon by Supercritical Carbon Dioxide[END_REF][START_REF] Andrews | Supercritical fluid extraction of aromatic contaminants from a sandy loam soil[END_REF][START_REF] Erkey | Supercritical Fluid Extraction of Heavy Molecular Weight Organic Contaminants from Soil[END_REF][START_REF] Kothandaraman | Supercritical extraction of polynuclear aromatic hydrocarbons from soil[END_REF]. When adsorption is concerned, thermodynamic and kinetic aspects should be involved to know more details about its performance and mechanisms.

In the framework of isothermal supercritical adsorption, there have been numerous publications in literature dealing with the modeling of adsorption equilibrium using the most common adsorption isotherm models, i.e. the Langmuir, the Freundlich and the Redlich-Peterson models [START_REF] Madras | Supercritical Fluid Regeneration of Activated Carbon Loaded with Heavy Molecular Weight Organics[END_REF][START_REF] Macnaughton | Supercritical Adsorption and Desorption Behavior of DDT on Activated Carbon Using Carbon Dioxide[END_REF][START_REF] Macnanghton | Proc. 3rd Int. Syrup. On Supercritical Fluids[END_REF][START_REF] Humayun | High-resolution adsorption isotherms of supercritical carbon dioxide on activated carbon[END_REF][START_REF] Aranovich | Adsorption isotherms for microporous adsorbents[END_REF][START_REF] Lochmüller | Adsorption isotherms on silica for methanol and 1-hexanol modifiers from supercritical carbon dioxide[END_REF]. All the experimental conditions used correspond to thermodynamic states relatively beyond the critical point because the adsorption equilibrium is influenced by the system temperature, pressure and by the supercritical fluid properties in the vicinity of the critical point. In contrast, supercritical adsorption systems close to the solvent's critical point have received much less explicit attention in the open literature. The experimental studies in this area are scarce. A thermodynamic analysis of near critical binary A c c e p t e d M a n u s c r i p t 4 mixtures was established by Afrane and Chimowitz [START_REF] Afrane | Adsorption in near-critical binary solvent mixtures: thermodynamic analysis and data[END_REF]. The authors studied the adsorption thermodynamics of dilute solutes adsorbing from high pressure supercritical fluid using the Henry's law. However, set under high pressures, the results showed an extremely weak dependence to pressure and to the composition of the supercritical solvent phase. In chromatography, the proximity to the critical point was early reviewed by Van Wasen et al. [START_REF] Van Wasen | Physicochemical Principles and Applications of Supercritical Fluid Chromatography (SFC)[END_REF]. The authors pointed out the unusual behavior of equilibrium partition coefficients in the near-critical region. Many other works also showed interesting features of data in this region [START_REF] Yonker | Mechanism of Solute Retention in Supercritical Fluid Chromatography, Supercritical Fluids: Chemical and Engineering Principles and Applications[END_REF][START_REF] Shim | Adjustable Solute Distribution Between Polymers and Supercritical Fluids[END_REF]; in particular, papers by Schmitz et al. [START_REF] Schmitz | Chromatography with Mobile Phases in the Liquid and the Supercritical State[END_REF] and Klesper and Schmitz [START_REF] Klesper | Gradient Methods in Supercritical Fluid Chromatography[END_REF] provided striking evidence of the highly nonlinear behavior of equilibrium coefficients with respect to pressure and temperature variations, as the critical point of the fluid phase is approached. We believe that an adequate explanation of the thermodynamic basis of these phenomena in adsorption process taking into account both temperature and pressure effects is necessary.

And it is also important to show the influence of the divergent character of thermodynamic properties and transport coefficients in near-critical systems on adsorption system behavior. This is precisely the aim of this paper. For this purpose, adsorption of a model solute from supercritical CO 2 was investigated in a small side-heated cavity by means of 2D numerical simulations. Naphthalene was chosen as a model solute because its phase equilibria with CO 2 has been thoroughly studied [START_REF] Raspo | Heterogeneous reaction induced by the piston effect in supercritical binary mixtures[END_REF][START_REF] Raspo | Fast Mass Transfer at a solid-supercritical fluid interface by Piston Effect[END_REF]. There are extensive data available for this system that have been confirmed. The first section of the paper is devoted to the mathematical modeling of the problem and the numerical method used for the simulations. The modeling of the adsorption reaction at the solid boundaries is exposed in details. Then, the effect of the mass fraction and the proximity to the critical point are discussed for wide temperature and pressure conditions. The results show a strong dependence to temperature and pressure variations when the critical point is approached. We ended up with the effect of Damköhler number on the adsorbed mass fraction.

A c c e p t e d M a n u s c r i p t

Mathematical modeling

Problem under investigation

The problem we consider is that of a dilute solute (Naphthalene in this case, named species 2) in supercritical CO 2 (named species 1). The physical properties of each pure compound are given in Table 1. The Naphthalene-CO 2 mixture is enclosed in a square cavity of height H=1mm and subjected to the earth gravitational field g. The cavity vertical walls are made of activated carbon (see Fig. 1). The activated carbon was chosen as a model adsorbent for this problem allowing as considering an adsorption reaction at the solid-fluid interface.

Here, we emphasize that the chosen mixture as the adsorbent material is only generic since the aim of this study is to qualitatively investigate the influence of the proximity to the critical point on an adsorption reaction. Initially, the fluid is considered in thermodynamic equilibrium at a constant temperature i T slightly above the mixture critical temperature

Governing equations

The mathematical model is based on the 2D time-dependent and compressible Navier-Stokes equations, coupled with energy and mass diffusion equations including the supplemental Peng-Robinson equation of state. In order to reduce computational costs, a low Mach number approximation is used [START_REF] Paolucci | On the filtering of sound from the Navier-Stokes equations[END_REF]. This approximation is valid since Mach numbers about 10 -4 are obtained. Thus, the total pressure is split into two parts: a homogeneous thermodynamic part ( ) th P t , which appears in the equation of state and in the energy equation and only depends on time t, and a non-homogeneous dynamic part ( , , ) dyn P x y t , appearing in the momentum equation and which varies with time and space. In this study, the dynamic pressure is strongly smaller than the thermodynamic part. Consequently, the total pressure is little different from the thermodynamic pressure and the evolution of th P governs that of the total pressure. In [START_REF] Accary | An adaptation of the low Mach number approximation for supercritical fluid buoyant flows[END_REF], a modification of the Low Mach number Approximation was proposed to account for the strong stratification of fluids near the critical point. We tested this modification and we noted that, for the present problem, the results obtained with and without the modification were the same. Therefore, the original approximation [START_REF] Paolucci | On the filtering of sound from the Navier-Stokes equations[END_REF] was used for the simulations reported in this paper. The dimensionless formulation was obtained using cm T as characteristic temperature, i ρ as characteristic density, ( )

1 / i cm R M T
ρ (with R is the perfect gas constant (R=8.3145 J mol -1 K - 1 )) as characteristic pressure, H as characteristic length, the time scale of the piston effect as characteristic time, ( ) ,

i i vi i C D µ λ
) . Thus, the governing equations in a dimensionless form are:

( ) 0 . = ∇ + ∂ ∂ V ρ ρ t (1) ( ) 1 1 1 . . Re 3 Re dyn P t Fr ρ ρ ρ ∂ + ∇ = -∇ + ∆ + ∇ ∇ + ∂ V V V V V (2)
A c c e p t e d M a n u s c r i p t

( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) 0 0 , * * * * * 0 2 1 2 1 0 , * 21 2 . 1 .
. 1 RePr

1

.

1 v th th w vi v th th w vi T C P T P T t C T C P T U U P T V V C T D w Le w ρ ρ ρ ρ γ γ λ γ ρ γ θ   ∂ ∂   + ∇ = - - - ∇     ∂ ∂         ∂   + ∇ ∇ - - + - - - ×         ∂         ∇ ∇ - V V (3) ( ) ( ) * 21 2 1 .
. 1

w w D w t Le ρ ρ ρ γ ∂ + ∇ = ∇ ∇ ∂ - V (4)
V is the velocity of components u and v in the x-and y-directions respectively, w is the mass fraction, γ is the ratio of the isobaric and isochoric specific heats calculated from the equation of state (see Appendix A) with 0 γ and 0 v C corresponding to the values for a perfect

gas 0 0 1 ( 1.3 , 3 / ) v C R M γ = =
. The value of vi C for the initial state was taken from the NIST (see Table 2). The dimensionless numbers are respectively, the Mach number Ma , the Reynolds number Re, the Froude number Fr, the Prandtl number Pr and the Lewis number Le and are defined as:

( )

2 0 21
, Re , , ,

PE i PE PE i vi i i i i vi i V V H V C Ma Fr Pr Le c g H C D ρ µ γ λ µ λ ρ γ = = = = = where 0 0 1 ( / ) cm c R M T γ = is the sound speed and PE PE H V t =
is the characteristic velocity of the piston effect. ( )

* 2 / k k vi cm U U M C T = and ( ) * 2 / / k k i V V M ρ = for k=1,2.
The expressions of k U and k V calculated using the Peng Robinson equation of state are given in Appendix B.

The following relationship is used for ( )

w θ :
( )

1 2 1 (1 ) M w w M θ = -- with 1
M and 2 M (kg mol -1 ) are respectively, the molecular weight of CO 2 and Naphthalene (see Table 1).

The superscript (*) refers to dimensionless parameters.

For thermal conductivity λ (W.m -1 .K -1 ), the following correlation is used [START_REF] Arai | Supercritical fluids: Molecular Interactions, Physical Properties and new Applications[END_REF]:

( ) ( ) ( ) ( ) 0 , , , e c T T T λ ρ λ λ ρ λ ρ = + + ∆
The first term ( ) 0 T λ corresponds to the limit of small densities and is expressed as follow: The second term ( ) e λ ρ is the excess property and is expressed as follow:

( )

c c T T exp C T λ ρ ρ ρ -   ∆ = - + × × - -   -  
C if ρ ρ ρ ρ - -  = × <   = × >  
The binary mass diffusion coefficient, D 21 (m 2 s -1 ), is calculated with the Wilke-Chang equation [START_REF] Higashi | Diffusion coefficients of naphthalene and dimethylnaphthalene in supercritical carbon dioxide[END_REF] The thermodynamic state of the mixture is described by the Ping-Robinson equation of state in the framework of the one-fluid theory. We can then compute the thermodynamic pressure as follow:

( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) * 2 2 2 * * * 2 , 1 / 1 2 / / th T w a T w P b w w b w w b w w ρθ ρ ρ θ ρ θ ρ θ = - - + - (5)
Where:

( ) ( )( ) ( ) ( ) ( ) 2 * * * * 2 1 12 2 , 1 2 1 , a T w a T w a T w w a T w = - + - + ( ) ( ) ( ) 2 * * * * 2 1 12 2 1 2 1 , b w b w b w w b w = - + - + ( ) ( ) ( ) 2 * 1 1 1 1 1 1.487422 1 1 / c i cm c cm c T a T T T T T ρ β ρ   = + -     ( ) ( ) ( ) 2 * 1 2 2 2 2 2 2 1.487422 1 1 / c i cm c cm c M T a T T T T M T ρ β ρ   = + -     * 0.253076 i j cj b ρ ρ        =  with j=(1,2)
A c c e p t e d M a n u s c r i p t

( ) ( ) ( ) ( ) 1 , a T a T a T k = - ( ) * * * 1 12 1 2 12 2 1 1 , 2 M b b b l M   = + -     2 0.37464 1.54226 0.26992 j j j β ω ω = + - (j=1,2)
with ω the acentric factor (Table 1). The binary interaction parameters 12 k and 12 l are determined so as to minimize the error between the calculated and experimental solubility data. These two parameters are temperature dependent and they are obtained through these formulae [START_REF] Nicolas | The Sanchez-Lacombe Lattice fluid model for the modeling of solids in supercritical fluids[END_REF]: ' '' 12 12 12 308. [START_REF] Aranovich | Adsorption isotherms for microporous adsorbents[END_REF]

= = = - = -
As it can be noted, the equations (1) -( 5) are coupled for a given time step. This coupling can be reduced by using an explicit scheme to evaluate the convective terms in Eq. (3). But the energy source term involving . V ∇ must be implicitly evaluated because it accounts for the piston effect, namely the thermoacoustic effect responsible for fast heat transfer near the liquid-gas critical point. So in order to decouple the energy equation (Eq. ( 3)) and the Navier- 2)), the velocity divergence, must be calculated using only the thermodynamic variables as explained by [START_REF] Ouazzani | A new numerical algorithm for low Mach number supercritical fluids[END_REF][START_REF] Raspo | Un algorithme faible nombre de Mach pour la simulation des écoulements de fluides supercritiques par des méthodes spectrales[END_REF].

As part of the low Mach number approximation and for the dimensionless equations and if we consider the equation of state written in the general form P th =F (T, ρ, w), the total derivative of F with respect to time t leads to the relation:

, , , th T T w w F d dP F dw dt dt w dt dT F dt T ρ ρ ρ ρ   ∂ ∂   - + -    ∂ ∂     = ∂     ∂   (6)
Moreover, the continuity equation ( 1) can also be written in the following form:

( )

V . ∇ - = ρ ρ dt d (7)
Then, inserting Eq. ( 6) in the energy equation Eq. (3) using Eq. ( 7) for the computation of dρ/dt and Eq. ( 4) for the computation of dw/dt, lead finally to the following expression for the velocity divergence:

( ) ( ) ( ) ( ) ( ) ( ) * * * , , , * 21 2 2 0 0 , , , 1 . , Re Pr . 1 . 1 1 1 th w T w v th th w w vi T w dP F F F T A U V dt T w w T D w Le C F P F P T C T T ρ ρ ρ ρ ρ γ ρ λ θ ρ γ γ ρ ρ     ∂ ∂ ∂       - ∇ ∇ - -           ∂ ∂ ∂           ∇ = ×       × ∇ ∇     -           ∂ ∂ ∂     - - -           ∂ ∂ ∂         V (8) with ( ) ( ) ( ) ( ) * * * * * * 0 2 1 2 1 0 , , 1
v th th w vi C P A U V U U P T V V C T ρ γ   ∂   = - + - - × -     ∂    
A c c e p t e d M a n u s c r i p t

The expressions of ,w

F T ρ ∂     ∂   , , T F w ρ ∂     ∂   and , T w F ρ   ∂   ∂   are reported in Appendix C.
Thanks to Eq. ( 8), we are now able to solve the governing equations in two uncoupled steps, namely the energy equation and the equation of state on the one hand and the Navier-Stokes equations on the other hand. The algorithm is detailed in section 2.4.

In this study, only the gravity effect was considered, the stratification of the fluid was not taken into account, because we have tested several cases with stratification and no effect was observed.

The initial and boundary conditions in dimensionless form are:

CI: 0 t = ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) * 2 2 * * * 1 1 1 1 , 1 / 1 2 / / i i i i i i i i i i i i T w fixed w a w P b w w b w w b w w ε ρ ε θ ε θ θ θ = + = + + = - - + - BC:
No-slip walls were considered so 0 u v = = at the two plates.

At x=0

* 1 h T T ε δ = + + with * / cm T T T δ δ = At x=1 1 co T ε = + At y=0, 1 0 T y ∂ = ∂ (adiabatic walls)
A c c e p t e d M a n u s c r i p t 14

Modeling of the adsorption reaction

At the interfaces x=0 and H, an adsorption reaction of Naphthalene on activated carbon is considered. The choice of such adsorption system can be explained by the extensive use of activated carbon as new-type high-efficiency adsorbent due to its high adsorption capacities and high mass transfer rates. However, the model description can be applied to any other adsorption system. The main objective here is to see how a supercritical mixture behaves in the vicinity of reactive wall. So, we will focus essentially on the fluid side rather than on what happens in the solid itself.

The species diffusion equation, (Eq. ( 4)), can be written in this form:

( ) 21 . w w w u v D w t x y ρ ρ ρ ρ ∂ ∂ ∂ + + = ∇ ∇ ∂ ∂ ∂
where u and v are the velocity of components in the x-and y-directions respectively.

The boundary conditions are:

At the horizontal non-reactive walls

0 w y ∂ = ∂
for y=0 and y=H

At the vertical reactive walls, a first order kinetic adsorption model is used:

21 . a w D K w n ∂ = ∂
where n is the normal to the surface at x=0 and x=H and K a =k a /S ac with k a the adsorption rate constant (m 3 .kg -1 .s -1 ) and S ac the specific surface area of activated carbon (m 2 . kg -1 ).

Assuming small variations of temperature and pressure, K a can be approximated by the firstorder term of its Taylor series in the vicinity of (T i , P i ):

( ) ( ) ( ) , a a a a i i i i p T K K K K T P T T P P T P ∂ ∂     = + - + -     ∂ ∂     (9)
A c c e p t e d

M a n u s c r i p t

The adsorption rate constant a k can be computed by 2 .

a d k K k =
where 2 K is the adsorption equilibrium constant (m 3 .kg -1 ) and d k is the desorption rate constant obtained through the Arrhenius law:

exp d Ea k A RT   = -    
with E a the activation energy and A the pre-exponential factor of Arrhenius.

The first term of Eq. ( 9) can be then considered as:

( ) ( ) ( ) 2 , , d i i i ai a i i ac k T K T P K K T P S = =
The partial derivatives of K a with respect to temperature and pressure are written as follow:

( ) ( )

2 2 2 , 1 d i i i a a a p p p ac ac k T K T P K k E LnK T S T S RT T   ∂ ∂ ∂       = = +         ∂ ∂ ∂         ( ) ( ) 2 2
, 1

d i i i a a T T T ac ac k T K T P K k LnK P S P S P ∂ ∂ ∂       = =       ∂ ∂ ∂      
Thus Eq. ( 9) becomes:

( ) ( ) ( ) ( )( ) 2 2 2 , 1 , , a a ai i i i i i i P T E LnK LnK K T P K T P T T T P P P RT T P     ∂ ∂       = + + - + -         ∂ ∂           (10)
The partial derivatives of K 2 with respect to temperature and pressure are obtained by equating differentials of the logarithm of the solute fugacity in the fluid and solid phases [START_REF] Jha | Modeling of adsorption equilibria in supercritical fluids[END_REF]:

( ) 2 2 2 2 2 ln m IG ads m p h h H K T RT α - + ∆ ∂   = +   ∂   (11)
A c c e p t e d M a n u s c r i p t

16 , 2 2 ln m m T K P RT ϑ κ ∞ ∂   = -   ∂   (12)
where As shown by Eq. ( 10), the derivatives of the equilibrium adsorption constant, K 2 , with respect to temperature and pressure are directly involved in the definition of the boundary condition at the solid-fluid interface. So, it is important to assess the effect of temperature and pressure on these derivatives. A strong sensitivity to temperature and pressure can be guessed from Eqs. ( 11) -( 12) since the infinite-dilution residual partial molar enthalpy, ( For pressures close to the critical one, namely P≤ 8.909 MPa, a divergence of the derivative is observed as the critical temperature is approached. Then, from 9 MPa, a behavior change can be depicted with the appearance of a maximum which decreases and is shifted to the high temperature domain when the pressure increases. One can also notice a similar trend for high pressures and high temperatures. Therefore beyond the critical point the effect of the pressure and temperature is no longer noticed.

330

For this reason, a particular attention was paid to the effects of temperature and pressure on the equilibrium adsorption constant in this study because a wide range of temperature and pressure were considered. Moreover, changes near and far from the critical point will help us to explain our results later.

It must be noted that the divergence of the infinite-dilution properties (residual partial molar enthalpy and partial molar volume) of the solute is not specific of Naphthalene but it is a universal behavior for dilute mixtures [START_REF] Chang | The critical dilemma of dilute mixtures[END_REF][START_REF] Chang | Behavior of dilute mixtures near the solvent's critical point[END_REF]. Therefore, the results obtained in this paper for Naphthalene in supercritical CO 2 are likely to be observed in all dilute binary mixtures, at least on a qualitative point of view. Moreover, the nature of the adsorbent material appears ads IG H h h ∆ = -, in Eq. ( 11), and therefore just for the variation of K 2 with respect to temperature. We compared the evolutions of d(LnK2)/dT for soil ( In the framework of the Low Mach number approximation, the boundary condition for the mass fraction w on x=0 and H can then be expressed as:

( ) 21 , ai k th w D K D T P w n ∂ = ∂ with ( ) ( ) ( ) ( )( ) 2 2 2 , 1 , , a k th i i i i i i P T E LnK LnK D T P T P T T T P P P RT T P   ∂ ∂     = + + - + -       ∂ ∂      
In the dimensionless form, it will be written as:

* * * 21 1 1 k D w w Da D n ∂ - = - ∂ (13)
D =1: * * 1 1 1 k w w Da D x ∂ - = ∂ (15)

Numerical method

The numerical integration of the model equations has been carried out using a second order semi-implicit scheme [START_REF] Vanel | A pseudospectral solution of vorticity-stream function equations using the influence matrix technique[END_REF]: the convective terms are evaluated by an Adams-Bashforth scheme, and then the time integration of the resulting differential equations has been done

with an implicit second order backward Euler scheme. The space approximation is performed using the Chebyshev-collocation method with Gauss-Lobatto points. For the computation of the convective terms, the derivatives are calculated in the spectral space and the products are performed in the physical one; the connection between the spectral and the physical spaces is realized through a FFT algorithm. On the other hand, the spectral differentiation matrices are used for the derivatives in the diffusive terms.

The computation of the velocity divergence by Eq. ( 8) allowed a decoupling between the thermodynamic variables T, ρ, P th , * w and that of the dynamic field. Consequently, the discretized equations can be solved in two successive steps: first, the thermodynamic variables are computed through the algorithm proposed by Ouazzani and Garrabos [START_REF] Ouazzani | A new numerical algorithm for low Mach number supercritical fluids[END_REF] and then the Navier-Stokes equations are solved using the modified projection method developed in [START_REF] Hugues | An improved projection scheme applied to pseudospectral methods for the incompressible Navier-Stokes equations[END_REF] and extended to variable density flows. These two steps are detailed bellow.

Computation of the thermodynamic variables (T, ρ, P th,

*

w )

The discretised energy and diffusion equations can be written as Helmholtz equations with time-dependent coefficients. In order to solve them using the diagonalization technique developed in [START_REF] Haldenwang | Chebyshev 3D spectral and 2D pseudospectral solvers for the Helmholtz equation[END_REF] for Helmholtz equations with constant coefficients, the density and the transport coefficients * λ and *

21

D are split into a constant part equal to the initial value and a time-dependent part:

( )

1 1 1 1 n n α α + + = + - for α=ρ , * λ , * 21 D
So, the discretized energy equation for example obtained as Helmholtz equation with constant coefficients is written as follow:

( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) (
)

1 1 1 1 * 1 1 1 , 1 1 1 1 1 1 0 0 , * 1 * 1 1 *
+ + + + + ∆ - = - - ∇ -∇   - + + + ∇       ∂   + - - ∇     ∂       - ∇ ∇ - V ( 16
)
where δt is the time-step and the notation AB(.) means an Adams-Bashforth evaluation of the quantity:

( )

1 1 , 2 - - - = n n n n AB φ φ φ .
The diagonalization process of the Helmholtz operator with constant coefficients is executed only once in a preprocessing stage. After that, at each time step, the solution of Eq. ( 16) is reduced to matrix products, leading to a very efficient solution technique.

A c c e p t e d M a n u s c r i p t

Thus, the solution of the energy Eq. ( 16), the diffusion equation and the equation of state is performed through the following iterative algorithm:

1. The variables T k-1 , P th k-1 , ρ k-1 , * ( 1)
)

1 1 1 1 1 , 1 1 1 1 1 1 0 0 , * 1 * 1 1 * 1 * 1 21 2 1 3 3 1 . 1 Re Pr 2 2 Re Pr 4 . 2 1 . 1 , . ( ) 1 k k k k k k n n n n k k k k k v th th w vi k k k k k i k T T T T t t T T A B T t C P P T V C T w A U V D w w Le ρ γ γ ρ λ δ δ ρ ρ δ γ ρ θ γ - - - - - - - - - - - - - - - - - ∆ - = - - ∇ -∇   - + + + ∇       ∂   + - - ∇     ∂       - ∇ ∇ - V 3. The mass fraction *k
w is obtained by the solution of the Helmholtz equation:

( ) ( ) ( ) ( ) 1 3 3 1 1 .(( 1
) ) 2 2 1 1 4 2 k k k k k k k n n n n k w w w D w t t Le Le w w A B V w t ρ ρ δ δ γ γ ρ ρ δ - - - - - - - - ∆ - = - - ∇ -∇ - -   - + + + ∇    
with the following Robin boundary conditions on adsorbing walls:

( ) 1 * * 1 * * 21 1 1 1 1 1 k k k k k k w D w w Da x Da Da D x - - -   ∂ ∂ + = + -   ∂ ∂   for x=0
( ) 1 * * * 1 1 1 1 1 1 k k k k k w w w Da x Da D x - -   ∂ ∂ - = - -   ∂ ∂   for x=1
4. The couple ( )

k k th P ρ
, is computed from the constraint of global mass conservation and the equation of state. This computation must be performed through an iterative process; ) for φ =T,ρ ,w * ,P th , and the maximum number of iterations is fixed to 250.

Solution of the Navier-Stokes equations

The second step is the solution of the Navier-Stokes equations. At the current time step (n+1), temperature, density and velocity divergence are known. A projection-type algorithm such as those developed for the solution of incompressible Navier-Stokes equations can then be used. In the present work, the original projection method of Hugues and Randriamampianina [START_REF] Hugues | An improved projection scheme applied to pseudospectral methods for the incompressible Navier-Stokes equations[END_REF] was modified to account for variable density flows [START_REF] Raspo | Un algorithme faible nombre de Mach pour la simulation des écoulements de fluides supercritiques par des méthodes spectrales[END_REF][START_REF] Ameur | Numerical simulation of the Poiseuille-Rayleigh-Bénard instability for a supercritical fluid in a mini-channel[END_REF]. The advantage of this method compared to other projection methods is that it allows improving the accuracy on pressure and reducing the slip velocity. It consists in solving the Navier-Stokes equations by three successive steps.

1 st Step: Computation of a preliminary pressure

The preliminary pressure Here again, we have to solve Helmholtz equations with variable coefficients for each velocity component. As for energy and mass diffusion equations, the density ρ n+1 is split into a constant part and a time-dependent part and the following Helmholtz equation with constant coefficient is solved iteratively:

( ) ( ) ( ) *, *, 1 *, 1 1 0 0 1 , 1 1 1 1 1 3 3 Re 2 2 4 . 2 1 1 . 3Re l l n l n dyn n n n n n n n T P t t A B t V g Fr ρ ρ ρ δ δ ρ ρ δ ρ + - + - - + + + ∆ - = - + ∇   - - + ∇     - ∇ ∇ - V V V V V V V ( 19
)
The convergence is achieved when Max(Resu, Resv) < 10 -13 ,with Resφ =Max((φ l -φ l-1 )/φ l-1 ) for φ =u, v. Only 3 or 4 iterations are necessary.

nd Step: Correction step

The converged velocity field * V is then corrected by taking into account the pressure gradient at the current time step (n+1) so that the final velocity field satisfies the continuity equation:

( ) ( ) ( )        Ω = ∇ + + - Ω ∂ = Ω ∂ ∪ Ω - -∇ = - + + - + + + + + + + + n i . on . . in * 0 2 4 3 2 3 1 1 1 1 1 1 1 1 1 1 1 n n n n n n B n n dyn n dyn n n n t n n P P t V V V V V ρ δ ρ ρ ρ ρ ρ δ ( 20
)
This system is solved through the following Poisson problem for the intermediate variable The system's response to the heating of the left vertical wall is explained in this section. The analysis below is based on simulations carried out at various initial temperatures 3. It can be noted that, as the critical temperature is approached, the thermal diffusion and mass diffusion times strongly increase whereas the characteristic time of piston effect decreases. For all the simulations, the left side of the cavity (at x=0) is gradually heated following a cosine low expressed as:

0.5 1 cos ( ) heat heat heat t T if t t t T t T if t t π δ      ∆ - ≤        =       ∆ ≥ 
where T ∆ is the temperature increase and t heat is the heating phase time corresponding to 200 time steps. The influence of T ∆ will be discussed later for different cases ( T ∆ =50, 100, 150 and 200 mK).

General description

In the first part, the case of i T =308.15 K and T ∆ =100mK, is discussed. The Damköhler number is fixed to Da=10 -5 and the initial mass fraction w i corresponds to the solubility of Naphthalene in CO 2 .

The evolution of temperature and density distributions between the walls and at the cavity mid-height are illustrated in Figs. associated with large density gradients (Fig. 5) along the vertical walls and the isothermal cavity bulk. Therefore, the temperature and the density profiles exhibit the same behavior, dominated by the piston effect, as in pure fluid [START_REF] Zappoli | Thermoacoustic and buoyancy-driven transport in a square side-heated cavity filled with near-critical fluid[END_REF][START_REF] Accary | Three-dimensional Rayleigh-Bénard instability in a supercritical fluid[END_REF]. In the rest of the paper, we denote as bulk the fluid region, which does not include the boundary layers. The homogeneous bulk temperature and density fields of a supercritical fluid occur at a time that is much shorter than the thermal diffusion time (t PE =0.23 s and t d =107s). They are the signature of the piston effect which was identified a long time ago as responsible for fast heat transport in near-critical pure fluids [START_REF] Boukari | Critical speeding up in pure fluids[END_REF][START_REF] Onuki | Fast adiabatic equilibration in a single-component fluid near the liquid-vapor critical point[END_REF][START_REF] Zappoli | Anomalous heat transport by the piston effect in supercritical fluids under zero gravity[END_REF]. However, the piston effect plays an important role only for short times because of the disappearance of sharp temperature gradients due to the action of thermal diffusivity.

For an advanced time (30s), the system's response, is then markedly different from that observed for shorter times. A similar trend to that of a perfect gas characterized by equilibrium can be depicted. It must be noted that temperature equilibrium is achieved on a time much shorter than the diffusion time. This is probably a consequence of the evolution of mass fraction in the cavity bulk, since w influences the temperature evolution as a source term of the energy equation.

The mass fraction field exhibits a different aspect as shown in Fig. 6. In order to explain the typical behavior of the supercritical mixture near the two reactive walls, we focus on the 10). Since the partial derivative of K 2 with respect to pressure is negative (see Fig. 2), the pressure term tends to reduce K a and as a result, to diminish the adsorbed amount. On the left heated wall, the effect of the temperature increase is amplified by the diverging derivative of K 2 with respect to temperature (see Fig. 3), leading to an important amount of Naphthalene adsorbed at the warm side.

The phenomenon occurring at the isothermal right wall is totally different. The strong density gradient near the boundary x=H (Fig. 5), generated by the piston effect, goes along with an increase of the amount of Naphthalene neat this reactive wall leading to an increase of the adsorbed quantity. Finally, the strong and homogeneous increase of the pressure, induced by the piston effect in the whole cavity, reduces the adsorbed amount at both reactive walls.

A c c e p t e d M a n u s c r i p t We consider now the case when the mixture is subjected to the gravity. In this case, side heating initiates gravity-driven convection in the fluid phase and the temperature and density fields obtained are completely two-dimensional on the contrary to the mainly 1D solutions previously observed in the absence of gravity. In Fig. 7, instantaneous temperature fields are plotted in the (x,y) plan for several times (1 s, 2.5 s, 4 s and 30 s). During a short time, the cavity bulk is heated rapidly due to the piston effect. As a result, upstream rises near This aspect is different from the perfect gas case where thermal boundary layer is formed only near heated side which leads to a single stream [START_REF] Vadim | Simulation of naturel convection in a sideheated cavity with a near-critical fluid[END_REF]. For the sake of comparison and in order to highlight the effect of Earth's gravity, we analyze in the following the Naphthalene mass fraction evolution with and without gravity (Fig. 8). The difference previously observed between the hot wall and the cold wall has disappeared. Gravity in this case tends to balance both sides. Conversely, the amount adsorbed is remarkably lower at the two sides (hot and cold) compared to the case without gravity. On the other hand, the mass fraction at the cavity centre has increased. These observations can be explained as follow:

A c c e p t e d M a n u s c r i p t 32

For the case without gravity, thermal diffusion process is the only highlighted. The characteristic time of thermal diffusion is much higher than the time scale of the piston effect.

Such a long process allows keeping a high temperature in the vicinity of the heated wall.

Consequently, as seen above, the divergence of the partial derivative of the equilibrium constant with respect to temperature near the critical point is believed to affect the mass fraction rate at the hot side. However, in the case with gravity, density variations generate a plume which expands upwardly (as shown in Fig. 7). This intensive plume moves hot fluid to the top boundary and thus decreases the temperature near the warm side. Therefore, the effect of the derivative with respect to temperature is reduced. Such a phenomenon can explain the reduction of the mass fraction when gravity is taken into account.

(a) (b) Fig. 9. The effect of the solid plate temperature on the mass fraction perturbation without gravity (a) and with gravity (b) for

T i =308.15 K at t=500s
The effect of heating intensity is highlighted in Fig. 9 with and without gravity. For the same proximity to the critical point and the same time, Fig. 9 shows a comparison of the mass fraction profiles at the cavity mid-height obtained for four values of the temperature rise ∆T= 50mK, 100mK, 150mK and 200mK: the mass fraction decreases with increasing heating. In the case without gravity, a stronger heating reduces remarkably the mass fraction throughout the entire volume while keeping the asymmetry of adsorbed amount between the hot and the A c c e p t e d M a n u s c r i p t cold walls (Fig. 9(a)). It can be noted that this asymmetry increases with the heating. A similar trend can be observed with gravity (Fig. 9(b)). This surprising result can be explained by the negative pressure term in the adsorption rate K a (Eq. ( 10)). Indeed, a stronger heating of the wall induces a larger thermodynamic pressure increase by the piston effect and, as a result, a larger value of the pressure term which reduces the parameter K a . The reduction of the adsorbed amount at the walls then leads to a smaller mass fraction in the whole cavity. dLnK2/dp=0 in Eq. ( 10) for T i =308.15 K and t=500s.

This explanation is confirmed by Fig. 10 where a case without the pressure term in Eq.

(10) was tested. We can then observe that with the sole presence of the temperature term, a stronger heating increases the mass fraction in the entire volume (Fig. 10(b)). Consequently, the pressure term plays a major role in the expression of the adsorption rate K a (Eq. ( 10)) near the critical point. The evolution of the mass fraction profiles as a function of heating depicted in Fig. 10(a) is then directly attributable to the piston effect which is responsible for the strong and homogeneous increase of the pressure in the cavity. 2,00 x10 -3 2.11x10 - 3 6,76 x10 -5 6,75 x10 -5 2,17 x10 -3 9.35x10 - 4 2,77 x10 -5 2,76 x10 -5 2,28 x10 -3 7.6751x10 - 4 2,25 x10 -5 2,24 x10 -5 2,30 x10 -3

Effect of initial mass fraction

The effect of the initial mass fraction on the adsorbed amount at the two walls is reported in Tables 4 and5 for cases with and without gravity and for two temperature increases ∆T=50mK and 100mK. The maximum value of w i corresponds to the solubility of Naphthalene in CO 2 . We note that the difference between the left and right sides observed in Fig. 8 without gravity can be quantified for each initial mass fraction and for the two temperature increase. Then, the effect of the gravity is confirmed in Table 5 where equilibrium is established between the two sides and only a slight difference can be observed between the heated and the isothermal plates. A detailed analysis of Tables 4 and5 reveals

Table 4

The effect of the initial mass fraction for T i =308.15 K, ∆T=100 mK and ∆T=50 mK at t=30s and without gravity

Table 5

The effect of the initial mass fraction for T i =308.15 K, ∆T=100 mK and ∆T=50 mK at t=30s and with gravity

A c c e p t e d M a n u s c r i p t 35 that, as could be expected, the mass fraction variation at both sides decreases with the initial mass fraction. To better quantify the difference between the hot and the cold sides, the relative

variation 0 0 / x x H x w w w = = =
is reported in the third column. In the first case, without gravity and for ∆T=100mK, a difference of 1.5% can be depicted for the greater initial mass fraction.

Then this difference evolves to 2.5% for the lowest value of w i . When the heating decreases to 50mK, the difference between the heated and isothermal sides is reduced (0.7% and 1.2% for high and low initial mass fraction respectively). However, it can be noted that, for both values of ∆T, the relative difference of adsorbed amount between the hot and cold sides increases when w i decreases and the increase is larger for the smaller heating, (71% for ∆T=50mK and 66% for ∆T =100mK). Table 5 shows that the gravity not only lowers the difference between the two walls (only 0.27% and 0.36% of difference can be observed for higher and lower initial mass fraction for ∆T=100mK and 0.18% and 0.23% for ∆T=50mK) but also reduces the effect of the initial mass fraction: the increase between the largest and the smallest values of w i is about 31% for ∆T=100mK and 26% for ∆T=50mK. These results then confirm the balancing influence of the gravity previously depicted by Fig. 8. The results presented in this section were obtained with an initial mass fraction w i corresponding to the solubility of Naphthalene in CO 2 at temperature T i and density ρ i . The values of w i and the initial pressure P thi calculated using the Peng-Robinson equation are reported in Table 6. A Similar behavior to that described in section 3.1 was observed for the mass fraction for different initial temperatures T i and temperature increases ∆T. However, when we move away from the critical point, the partial derivatives of the equilibrium constant K 2 with respect to temperature and pressure decrease (Figs. 2 and3) and thus, the mass fraction at the heated and isothermal sides are influenced. Moreover, it must be noted that, since the initial mass fraction w i is fixed to the value corresponding to the solubility at T i and ρ i , the resulting pressure gets higher as the initial temperature moves away from the critical one. As a consequence, the initial pressure belongs to the high pressure range where the derivative of K 2 with respect to pressure is smaller. For all the initial temperatures (T i =307.75 K to T i =318.15 K), the piston effect generated by the boundary heating induces a fast and strong pressure rise in the entire volume before the pressure reaches a steady value. When the system is subjected to the Earth's gravity, convection accelerates the pressure increase (Fig. 11(b)). Yet, with or without gravity, above the critical point the piston effect becomes less effective and the pressure plateau for a given ∆T gradually decreases. For T i =318.15 K, steady value is much smaller The lessening of the piston effect depicted by Fig. 11 is confirmed by Tables 7 and8 which

Proximity of the LCEP

show that the pressure value corresponding to the equilibrium state after t=30s decrease when moving away from the critical point for both temperature increases of 50mK and 100mK and for the two cases with and without gravity effect. It can be also noted that the gravity has a very little influence on the pressure evolution, since very close values are obtained for the pressure with and without gravity. Tables 7 and8 also show how the proximity to the critical point affects the adsorbed amount at the two reactive walls. In the two cases, with or without gravity, a change in the variation of the mass fraction as a function of initial temperature can be observed for the farest value of T i . Indeed, the mass fraction at the heated and isothermal plates regularly decreases when the system moves away from the critical temperature up to T i =311.15 K. Then, a strong increase of the mass fraction is observed at T i =318.15 K. A similar behavior change at the highest temperature can also be noted on the relative gap of A c c e p t e d M a n u s c r i p t mass fraction between the two sides. For example, for ∆T=100mK, the difference between the two reactive walls increases from about 1.4% at T i =307.75 K to 1.6% at T i =311.15 K and then decreases to 1% at T i =318.15 K. This abrupt behavior change can be attributed both to the reduction of the piston effect (leading to a smaller pressure increase) and to the decrease of the derivative of the equilibrium constant K 2 with respect to pressure. These two phenomena lead to a decrease of the negative pressure term in the adsorption rate expression (Eq. ( 10)).

Table 8

The effect of the proximity to the critical temperature ∆T=100 mK and ∆T=50 mK at t=30s and with gravity

Table 7

The effect of the proximity to the critical temperature for ∆T=100 mK and ∆T=50 mK at t=30s and without gravity The consequence of this decrease of the pressure term in Eq. ( 10) at the highest initial temperature is depicted by Fig. 14 which shows the mass fraction variation between the two plates for three temperature increases (∆T=50, 100 and 200mK) and two initial temperatures (T i =311.15 K and 318.15 K). Though the tendency for T i =311.15 K is similar to that shown in Fig. 10 for T i =308.15 K, the behavior of the mass fraction distribution for T i =318.15 K is reversed: a stronger heating of the left side increases the mass fraction. The effect is more remarkable at the heated side than at the isothermal one where the profiles of ∆T=50mK and 100mK merge. This behavior change is due to the competition between the derivatives of the equilibrium constant K 2 with respect to temperature and pressure. Far away from the critical point, the diverging behaviors of the isothermal compressibility and of the volume expansivity disappear leading to smaller variations of the derivatives ( )

2 ln / p K T ∂ ∂ and ( ) 2 ln / T K P ∂ ∂ (see Appendix D).
Moreover, for a given heating intensity, the pressure increase generated by the piston effect is much lower for the highest initial temperature. As a result, the negative pressure term in Eq. ( 10) becomes negligible and only the temperature effect is highlighted and causes the enhancement of mass fraction with heating increase. The results presented up to now were obtained for a Damköhler number fixed to 10 -5 .

With reference to literature studies, it has been found that the Damköhler number for the adsorption of Naphthalene and other solutes such as toluene or benzene can vary from 10 -3 to 10 -14 [START_REF] Aisien | Batch study, equilibrium and kinetics of adsorption of naphthalene using waste tyre rubber granules[END_REF][START_REF] Cabal | Kinetics of naphthalene adsorption on an activated carbon : Comparison between aqueous and organic media[END_REF][START_REF] Long | Adsorption of naphthalene onto the carbon adsorbent from waste ion exchange resin : Equilibrium and kinetic characteristics[END_REF][START_REF] Lucas | Mathematical model of supercritical CO2 adsorption on activated carbon effect of operating conditions and adsorption scale-up[END_REF][START_REF] Shi | Adsorption of naphthalene onto a high-surface-area carbon from waste ion exchange resin[END_REF]. The Damköhler number was estimated using the available data in each research work. For the sake of comparison, three values of the Damköhler number were tested. Similar behaviors to those reported in the previous sections were found for temperature, pressure and mass fraction distribution. Only the mass fraction variations at the heated left side and the isothermal right one are presented in Table 9. The same tendency was found with high mass fraction at the heated plate for all the Damköhler numbers. Whereas the Naphthalene mass fraction is found to be very much smaller for the smallest Damköhler number, increasing Da by a decade results in an increase by a decade of the mass fraction for the larger values of Da.

Conclusion

In this paper we have presented new results and a detailed analysis of adsorption in a model binary dilute mixture, the Naphthalene-CO 2 mixture, very close to the critical point.

The results of this study revealed that sufficiently close to the mixture critical point, the increase of the wall heating remarkably affects the adsorbed amount at the two reactive boundaries and the mass fraction inside the cavity. More precisely, the adsorbed amount, as the bulk mass fraction, is reduced by increasing the wall heating. This peculiar behavior is attributed to the Piston effect, coupled with the divergent character of the derivative of the The effect of the Damköhler number on the mass fraction perturbation for T i =308.15 K and ∆T=100 mK A c c e p t e d M a n u s c r i p t adsorption equilibrium constant with respect to pressure. Far enough from the critical point, the Piston effect weakens and a classical behavior is observed. Our results also showed that this retrograde adsorption is obtained with and without gravity. However, in the presence of gravity, convection induces large thermal plumes along the hot and cold boundaries and tends to reduce the temperature gradients near the two walls leading to more symmetric profiles of the mass fraction. Finally, the effect of the Damköhler number was studied. The same behavior was found for all the values considered.

All the results presented in this paper were obtained for the Naphthalene-CO 2 model mixture. However, we believe that this study can be relevant for many dilute binary mixtures.

Indeed, the phenomena observed are due to the divergence of the solvent transport properties (namely the isothermal compressibility and the thermal expansion coefficient) near the critical point leading to the appearance of the Piston effect and to the divergence of the solute thermodynamic properties (such as infinite dilution partial molar volume). And these divergent behaviors occur in a universal way for large classes of systems. Therefore, similar results should be obtained for all binary dilute mixtures involving a non-volatile solute near the solvent's critical point and this kind of dilute mixtures is relevant for many adsorption processes.

A c c e p t e d M a n u s c r i p t

            ∂ ∂ ∂ ∂     = - - + -             ∂ ∂ ∂ ∂                 ( ) 2
  ∂     - ∂     + -   + - -   + -   -= - + - - + - (B.2)
where ϑ is the molar volume of mixture.

Appendix C

In the expression of the velocity divergence ( . V ∇ ) (Eq. ( 8)), the derivatives are calculated as

    ∂   = - + - -       ∂   -       -   ∂ + - × - × -     ∂     (C.
  + -     ∂   + -- -     ∂   + +         (D.4)

ε defines the dimensionless proximity to the critical point ( 1 εA c c e p t e d M a n u s c r i p t 6 Fig. 1 .

 161 Fig. 1. Physical configuration

  t is the characteristic time of thermal diffusion, m γ the capacity ratio of the mixture (see Appendix A) and PE H t was taken as the characteristic velocity. The transport properties such as the dynamic viscosity µ, the isochoric specific heat capacity C V , the thermal conductivity λ and the diffusion coefficient 21 D were dimensionless, relative to their respective initial values (

V

  are respectively the dimensionless partial molar internal energy and partial molar volume expressed as follow:

  :

  3 mol -1 ) the molar volume of Naphthalene at boiling point, Φ the association factor ( Φ =1 for CO 2 ).

A c c e p t e d M a n u s c r i p t

  Stokes equations (Eqs. (1)-(

2 m h and 2 sh 2 IGh

 222 of adsorption of the solute on the solid plate, are the infinite-dilution partial molar enthalpies of solute in the mobile and stationary phases, respectively, is the enthalpy of the solute in the ideal gas state, volume of the solute in mobile phase and m α , m κ are respectively, the volume expansivity and the isothermal compressibility.The infinite-dilution residual partial molar enthalpy ( 2 2 m IG h h -) of the solute, m α , m κ and , 2 m ϑ ∞ are obtained using the Peng-Robinson equation of state and they are reported in Appendix D.

  infinite dilution partial molar volume of the solute, , 2 m ϑ ∞ , the isothermal compressibility, m κ , and the volume expansivity, m α , diverge near the solvent critical point. This assumption is confirmed by Fig. 2 which shows the K 2 derivative profiles as a function of pressure for different temperatures. The profiles show a sharp minimum which becomes more important when the critical temperature is approached. Then this minimum is shifted to the high pressure domain and becomes less significant away from the critical point especially for T=318.15 K. In this region (high pressure domain), the effect of the temperature on the isothermal derivative is less pronounced. In the same way, the isobaric derivative of the A c c e p t e d M a n u s c r i p t equilibrium constant as a function of temperature is shown in Fig. 3 for different pressures.

Fig 3 .

 3 Fig 3. Equilibrium constant derivative vs. temperature

  J/mol) as adsorbent material. The same variations were observed in both cases with maxima occurring for the same temperatures but with slightly different values: for example, the maximum for P= 9.221 MPa was equal to 0.52383 K -1 for soil and 0.47643 K -1 for activated carbon. These identical variations can be explained by the fact that the evolution of d(LnK2)/dT is governed by the very large values of the infinite-dilution residual partial molar enthalpy, ( 2 2 m IG h h -), and of the volume expansivity, m α , near the critical point. And these two properties are completely independent of the characteristics of the adsorbent material. Consequently, the results presented in this paper for adsorption on activated carbon are also relevant for any other adsorbent material.

  is the Damköhler number defined as the ratio of the characteristic fluidic time scale (diffusion characteristic time) and chemical time scale (adsorption characteristic time): 13) leads to two Robin-type boundary conditions at x=0 and x=1: At x=1, the wall is maintained at the initial temperature, thus, * 21

A c c e p t e d M a n u s c r i p t 20

 20

2 .

 2 The temperature T k is obtained by the solution of the Helmholtz equation:

A c c e p t e d M a n u s c r i p t 22 5. 21 kD are updated. 6 .

 22216 The thermal conductivity *k λ and the diffusion coefficient * The velocity divergence ( ) k T V . ∇ is computed by Eq. (8). The steps 2 to 6 are repeated until convergence is achieved on temperature, thermodynamic pressure, density and mass fraction. The convergence criterion used is Max(ResT, Resρ ,Resw, ResP th ) < 10 -11 ,with Resφ =Max(

  a Poisson equation, derived from the discretized momentum equation, with Neumann boundary conditions obtained by the normal projection of the momentum equation on the boundary: A c c e p t e d M a n u s c r i p t

Fig. 4 .Fig. 6 .

 46 Fig. 4. Profiles of the temperature perturbation at several times for T i =308.15 K and ∆T=100 mK in the case g=0

Fig. 5 . 27 iT

 527 Fig. 5. Profiles of the density perturbation at several times for Ti=308.15 K and ∆T=100 mK

25

 25

.

  The values of the different characteristic times of the problem are given in Table

4 and 5

 5 for several times (t = 0.5 s, 1 s, 2.5 s, 4 s and 30 s) in the absence of gravity. Because of the very small thermal diffusivity near the critical point, the heating of the wall causes ultra-thin boundary layers at the wall-fluid interface. Due to the high isothermal compressibility, the fluid close to the heated side expands upward and converts some of the kinetic energy into thermal energy. This results in compressing adiabatically the rest of the fluid and leading to a quick increase of the thermodynamic pressure which induces a fast and homogeneous heating of the cavity bulk by thermo-acoustic A c c e p t e d M a n u s c r i p t effects (piston effect). During the heating phase, temperature, and therefore density, at the heated wall change at each time step and, as a consequence, temperature and density gradients near this wall increase more and more. As the bulk temperature grows, and since the right wall (at x= H) temperature is maintained at its initial value, a cold boundary layer settles near the right wall where the fluid contracts. The contraction causes an expansion of the bulk and reduces the bulk temperature. The behavior of the fluid is the result of these two competing processes, heating by hot boundary layer and cooling by cold boundary layer. The temperature field is then divided into three distinct zones: two thermal boundary layers

A c c e p t e d M a n u s c r i p t 29 boundary conditions developed in section 2 . 3 .

 2923 The dependence to temperature and pressure of the adsorption rate can be depicted from Eq. (

599Fig. 7 .

 7 Fig. 7. Isosurfaces of the temperature variation for several times (a) 1s, (b) 2.5s, (c) 4s and (d) 30s, for T i =308.15 K and ∆T=100mK in the presence of gravity

  s c r i p t the left warm surface and a hot spot at the left corner of the cavity can be depicted after 1s and then it is convected progressively along the top wall for longer times. As in Fig.4for g=0, a cool boundary layer forms near the right isothermal wall, due to the piston effect. As a result, a jet moving down appears near the right wall from 2.5 s. The hot and cold thermal plumes along the top and the bottom sides develop progressively in time as shown in Fig.7. We can clearly see the homogenous increase of the bulk temperature induced by the piston effect.

Fig. 8 .

 8 Fig. 8. The effect of the gravity on the mass fraction perturbation for T i =308.15 K and ∆T=100 mK at t=500s

Fig. 10 .

 10 Fig. 10. Influence of the piston effect on the mass fraction perturbation without gravity with (a) dLnK2/dp≠0 and (b)

Fig. 11 .

 11 Fig. 11. The effect of the proximity to the critical temperature without (a) and with (b) gravity at ∆T=100 mK.

A c c e p t e d M a n u s c r i p t 36

 36

A c c e p t e d M a n u

  s c r i p tshowing that far enough from the critical point the previously observed effects on temperature and mass fraction will be reduced.The effect on temperature and density can actually be observed in Figs.12 and 13comparing the instantaneous temperature and density fields near and far from the critical point (for T i =307.75 K and T i =318.15 K). For these two initial temperatures, the characteristic time scales of the piston effect are t PE =0.19 s and t PE =1.26 s respectively. On the other hand, the characteristic time of thermal diffusion is t d =115.7 s for T i =307.75 K and t d =49.23 s for T i =318.15 K. Therefore, far from the critical point, the piston effect decreases in favor of the thermal diffusion as clearly shown by the temperature field in Fig.12(b). The boundary layers become thicker and the thermal plumes are larger (Figs 12(b) and 13(b)). The figures also show that, near the critical point (T i =307.75 K), top and bottom plumes reach the opposite plate while the reduction of the piston effect away from the critical point reduces this phenomenon.

A c c e p t e d M a n u s c r i p t 39 753Fig. 12 .Fig. 13 .

 391213 Fig. 12. Isosurfaces of temperature with gravity for (a) T i =307.75 K and (b) 318.15 K for ∆T=100 mK at t=500s.

40 Fig. 14 .

 4014 Fig. 14. Evolution of the mass fraction with the heating intensity for (a) T i =311.15 K at t=500s and (b) T i =318.15 K at t=500s

u s c r i p t 41 3. 4

 414 Influence of the Damköhler number

1 )

 1 The derivatives are calculated using the Peng-Robinson equation for pure CO 2 :In Eqs. (A2) and (A3), ( ) a T and b are the coefficients of the Peng-Robinson equation of state written for mass variable in dimensional form.

  way, the capacity ratio of the mixture, m γ is calculated as follow: are calculated using the Peng-Robinson equation of state for the mixture: in which the derivatives i da dT are calculated as described above for pure component.

2 mϑ

 2 [START_REF] Erkey | Supercritical Fluid Extraction of Heavy Molecular Weight Organic Contaminants from Soil[END_REF]-[START_REF] Kothandaraman | Supercritical extraction of polynuclear aromatic hydrocarbons from soil[END_REF], the volume expansivity, α , the isothermal compressibility, κ , the partial molar volume, and the partial molar residual enthalpy,

Table 1

 1 Pure component properties

Table 2

 2

	Initial parameters
	T i (K) 307.75	ρ ρ ρ ρ i (kg.m -3 ) 470	λ i (W.m -1 .K -1 ) 0.098332532	Cv i (J.kg -1 .K -1 ) 1325.839	µ i (Pa.s) 3.33828 x10 -5	(D 21 ) i (m 2 .s -1 ) 2.19525 x10 -8
	308.15	470	0.096196327	1306.27	3.34016 x10 -5	2.19686 x10 -8
	309.15	470	0.091343811	1269.32	3.34485 x10 -5	2.20090 x10 -8
	311.15	470	0.083368892	1214.86	3.35423 x10 -5	2.20895 x10 -8
	318.15	470	0.066814297	1074.57	3.38700 x10 -5	2.23679 x10 -8

Table 3

 3 Characteristic times (piston effect t PE , thermal diffusion t d , mass diffusion t Md , adsorption t ad )

T i (K) t PE (s) t d (s) t Md (s) t ad (s)

	307.75	0.1999	115.6989	45.5530	45.5530 x10 5
	308.15	0.2367	106.9022	45.5194	45.5194 x10 5
	309.15	0.3333	91.1253	45.4359	45.4359 x10 5
	311.15	0.5404	73.0746	45.2704	45.2704 x10 5
	318.15	1.2621	49.2332	44.7069	44.7069 x10 5

Table 6

 6 Initial mass fraction w i and pressure P thi

	T i (K) 307.75	ρ ρ ρ ρ i (kg.m -3 ) 470	wi 7.6751x10 -3	P thi (MPa) 8.67158154
	308.15	470	4.22x10 -3	8.73970242
	309.15	470	2.11x10 -3	8.90942413
	311.15	470	9.35x10 -4	9.22114740
	318.15	470	7.6751x10 -4	10.3949472

100 mK T i (K) w-w i (x=0) w-w i (x=H) (w x=0 -w x=H )/w x=0 (P th -P th0 )/P th0 307.75

			100 mK
	T i (K) 307.75	w-w i (x=0) 4.09 x10 -4	w-w i (x=H) 4.03 x10 -4	(w x=0 -w x=H )/w x=0 1.3885 x10 -2	(P th -P th0 )/P th0 3.8347 x10 -11
	308.15	3.93 x10 -4	3.88 x10 -4	1.4824 x10 -2	3.8053 x10 -11
	309.15	3.75 x10 -4	3.69 x10 -4	1.6061 x10 -2	3.7326 x10 -11
	311.15	3.71 x10 -4	3.65 x10 -4	1.6067 x10 -2	3.6028 x10 -11
	318.15	4.77 x10 -4	4.72 x10 -4	1.0726 x10 -2	3.1747 x10 -11
			50 mK
	T i (K) 307.75	w-w i (x=0) 4.10 x10 -4	w-w i (x=H) 4.07 x10 -4	(w x=0 -w x=H )/w x=0 6.9607 x10 -3	(P th -P th0 )/P th0 1.9197 x10 -11
	308.15	3.94 x10 -4	3.91 x10 -4	7.4364 x10 -3	1.9046 x10 -11
	309.15	3.75 x10 -4	3.72 x10 -4	8.0633 x10 -3	1.8675 x10 -11
	311.15	3.70 x10 -4	3.67 x10 -4	8.0684 x10 -3	1.8016 x10 -11
	318.15	4.76 x10 -4	4.73 x10 -4	5.3793 x10 -3	1.5852 x10 -11
		3.94 x10 -4	3.93 x10 -4	2.6531 x10 -3	3.8326 x10 -11
	308.15	3.78 x10 -4	3.77 x10 -4	2.7778 x10 -3	3.8036 x10 -11
	309.15	3.58 x10 -4	3.57 x10 -4	2.9648 x10 -3	3.7315 x10 -11
	311.15	3.53 x10 -4	3.52 x10 -4	3.0537 x10 -3	3.6024 x10 -11
	318.15	4.56 x10 -4	4.55 x10 -4	2.3868 x10 -3	3.1752 x10 -11
			50 mK
	T i (K) 307.75	w-w i (x=0) 3.98 x10 -4	w-w i (x=H) 3.97 x10 -4	(w x=0 -w x=H )/w x=0 1.7511 x10 -3	(P th -P th0 )/P th0 1.9199 x10 -11
	308.15	3.82 x10 -4	3.81 x10 -4	1.8280 x10 -3	1.9049 x10 -11
	309.15	3.62 x10 -4	3.61 x10 -4	1.9335 x10 -3	1.8679 x10 -11
	311.15	3.56 x10 -4	3.55 x10 -4	1.9657 x10 -3	1.8020 x10 -11
	318.15	4.59 x10 -4	4.58 x10 -4	1.5165 x10 -3	1.5860 x10 -11

Table 9

 9

Acknowledgments: The authors acknowledge the financial support from the CNES (Centre National d'Etudes Spatiales).

with Ω the computational domain (Ω=]-1,+1[×]-1,+1[), ∂Ω its boundary,

the boundary conditions of the velocity V n+1 , ∂/∂n the normal derivative and . ∇ V is calculated from thermodynamic variables and noted ( )

In Eq. ( 17) , the term ∆V n+1 was decomposed in the boundary condition using the formula:

( ) ( )

and the rotational term was evaluated using an Adams-Bashforth scheme.

nd

Step: Computation of a predicted velocity V *

The predicted velocity field V * is computed implicitly from the momentum equation with the gradient of the preliminary pressure instead of that of the actual pressure 1 + n dyn P . The predicted velocity therefore satisfies the following problem:

Appendix A

The ratio of the isobaric and isochoric specific heats for pure CO 2 , γ , and for mixture, m γ , are calculated from the equation of state as follow:

A c c e p t e d M a n u s c r i p t

The difference of the partial molar internal energies of the two components is expressed by:

with ( )

H T and ( )

H T the perfect gas enthalpy of the two components at T 0 =298.15 K and y the mole fraction of component 2, calculated from the mass fraction by the formula:

( )

Cp and 0 1

Cp are the isobaric heat capacities of components 2 and 1 respectively as perfect gas and their difference is expressed by: (

where 12 12 , , , ,