Amine Ammar

Emmanuelle Abisset-Chavanne

Francisco Chinesta
email: francisco.chinesta@ec-nantes.fr

Roland Keunings

Flow modelling of quasi-Newtonian fluids in two-scale fibrous fabrics Advanced simulations

Keywords: Permeability, Computational homogenization, Quasi-Newtonian fluids, Generalized Newtonian fluids, Extensional viscosity, Composite materials, Model order reduction, Proper generalized decomposition, Parametric solutions

Permeability is the fundamental macroscopic material property needed to quantify the flow in a fibrous medium viewed as a porous medium. Composite processing models require the permeability as input data to predict flow patterns and pressure fields. In a previous work, the expressions of macroscopic permeability were derived in a double-scale porosity medium for both Newtonian and generalized Newtonian (shear-thinning) resins. In the linear case, only a microscopic calculation on a representative volume is required, implying as many microscopic calculations as there are representative microscopic volumes in the whole fibrous structure. In the non-linear case, and even when the porous microstructure can be described by a unique representative volume, a large number of microscopic calculations must be carried out as the microscale resin viscosity depends on the macroscopic velocity, which in turn depends on the permeability that results from a microscopic calculation. An original and efficient offlineonline procedure was proposed for the solution of non-linear flow problems related to generalized Newtonian fluids in Francisco Chinesta

Introduction

Isothermal flows of complex fluids in complex microstructures can be simulated by solving the momentum and mass balance equations and a suitable rheological constitutive model. For inertialess incompressible flows, these balance equations read,

∇ • σ = 0, (1)
and

∇ • v = 0, (2)
respectively. Here, σ is the Cauchy stress tensor and v the velocity field, both defined at time t at each point within the fluid domain f . When considering porous media, the domain is assumed fully saturated, with the fluid phase occupying the region f whereas the remaining part s = f is occupied by a solid phase assumed at rest. An appropriate constitutive equation must be postulated to describe the fluid's rheology. There are many possible choices, the most usual ones being related to Newtonian, generalized Newtonian and viscoelastic fluids, as summarized below.

-For a Newtonian fluid, the constitutive equation reads σ = -pI + τ = -pI + 2ηD, [START_REF] Binetruy | Flows in Polymers, Reinforced Polymers and Composites. A multiscale approach[END_REF] where p is the pressure field that can be interpreted as the Lagrange multiplier associated with the incompressibility constraint, I is the identity tensor, τ the extra-stress tensor, η the constant fluid viscosity and D the rate of strain tensor, i.e. the symmetric part of the velocity gradient, 2D = ∇v + (∇v) T . -For a generalized Newtonian fluid, the constitutive Eq.

3 remains formally unchanged but the viscosity η now depends on the strain rate γ usually expressed from the second invariant of the rate of strain tensor, i.e. γ = √ 2D : D. The simplest of such models is the power-law (shear-thinning) viscosity given by

η = K γ n-1 , (4)
where K and n are the consistency and power-law index, respectively. The value n = 1 corresponds to a Newtonian fluid. -Viscoelastic models combine viscous and elastic behaviours and results in a history-dependent extrastress tensor τ . Both integral and differential models exist, but the latter are more popular in flow simulations. They have the deceptively simple generic form

Dτ Dt = F(τ , ∇v), (5)
meaning that the Lagrangian or material time derivative of the extra-stress tensor is given as a model-dependent function F of the local extra-stress and velocity gradient evaluated along fluid trajectories.

In all these cases, suitable boundary conditions must be specified at the fluid domain boundary = ∂ f . For inelastic fluids, they are of two types: prescribed velocity on D and prescribed tractions on N , with D ∪ N = and D ∩ N = ∅. Inlet extra-stresses must also be specified for viscoelastic fluids.There is a variety of discretization techniques for solving numerically such a flow problem. The interested reader can refer to [9] and the numerous references therein. Advanced computational techniques for differential and integral viscoelastic models are detailed in [START_REF] Binetruy | Flows in Polymers, Reinforced Polymers and Composites. A multiscale approach[END_REF][START_REF] Keunings | Simulation of viscoelastic fluid flow, in Fundamentals of Computer modeling for Polymer Processing[END_REF][START_REF] Keunings | Finite Element Methods For Integral Viscoelastic Fluids[END_REF][START_REF] Keunings | Micro-macro Methods for the Multi-Scale Simulation of Viscoelastic Flow using Molecular Models of Kinetic Theory[END_REF][START_REF] Owens | Computational rheology[END_REF].

Describing the internal structure of heterogeneous materials at a variety of length scales requires a high-enough resolution for capturing the size, orientation and distribution of the material constituents, i.e. the heterogeneity and its effects on the velocity field. Fibrous media with a large number of fibres and inter-fibre volumes cannot be described by considering each of these individually, which would lead to intractable boundary value problems to solve for the flow. Thus, another continuum approach at a coarser level is called for. One possibility consists in separating the scales (when this can be performed) and associating to each point X ∈ a representative volume at the finest scale ω(X) (with boundary ∂ω(X)) from which upscaling can be carried out.

When the same physics is valid and applies at each scale, for example in the case of thermo-mechanical models, homogenization procedures work quite well in both the linear and non-linear cases. Linear and non-linear thermomechanical homogenization, including its computational counterpart, has been the subject of active research over the last three decades. The interested reader can refer to [10, [START_REF] Nguyen | Multiscale continous and discontinous modeling of heterogeneous materials: a review on recent developments[END_REF][START_REF] Castañeda | Nonlinear composites[END_REF] and the references therein. The main issue in computational homogenization is the consideration of nonlinear behaviours that involve heavy computations. Use of model order reduction techniques is an appealing route for improving the efficiency of multi-scale solvers [START_REF] Chinesta | Alleviating mesh constraints: Model reduction, parallel time integration and high resolution homogenization[END_REF][START_REF] Lamari | Routes for efficient computational homogenization of Non-Linear materials using the proper generalized decomposition[END_REF][START_REF] Lopez | Advanced thermal simulation of processes involving materials exhibiting fine-scale microstructures[END_REF].

In the case of flows in porous media, the macroscopic scale corresponds to the scale of the part, much bigger than the heterogeneities scale, whereas the microscopic scale coincides with the scale of heterogeneities. In these circumstances, and for a Newtonian fluid, the macroscopic behaviour is well described by the phenomenological Darcy law that relates the fluid velocity to the pressure gradient at the macroscopic scale, whereas at the microscopic level the flow behaviour is described by the Stokes problem defined in the domain occupied by the fluid (in this work, for the sake of simplicity, we consider the solid phase as rigid and at rest). Even if both descriptions are notably different, appropriate bridges can be defined to associate to each point X in a permeability tensor K(X) relating locally the equivalent macroscopic velocity V(X) and the pressure gradient ∇P (X), by upscaling the microscopic physics that only involves as material parameter the fluid viscosity η.

In [START_REF] Lopez | Flow modeling of linear and nonlinear fluids in two and three scale fibrous fabrics[END_REF], we derived the expressions of the macroscopic permeability in a double-scale porous medium for both Newtonian and shear-thinning fluids. In the linear case, a single microscopic calculation on a representative volume is required, which implies as many microscopic calculations as there are representative microscopic volumes in . In the non-linear case, and even when the microstructure can be described by a unique representative volume, the microscopic calculation must be carried out repeatedly as the microscale fluid viscosity depends on the macroscopic velocity, which in turn depends on the macroscopic permeability obtained from the microscopic calculation. In [START_REF] Lopez | Flow modeling of linear and nonlinear fluids in two scale fibrous fabrics[END_REF], we proposed an original and efficient offline-online procedure that allows one to solve the multi-scale shear-thinning fluid flow problem at a computational cost typical of a macroscopic non-linear solution. This procedure is based on the calculation of the microscopic parametric solution, given the velocity field at each point in the representative volume ω(X), v(x ∈ ω(X)), for any value of the macroscopic velocity field V, v(x ∈ ω(X); V). From the parametric flow kinematics, parametric permeabilities K(V) were successfully and efficiently computed.

The route followed in [START_REF] Lopez | Flow modeling of linear and nonlinear fluids in two scale fibrous fabrics[END_REF] was based on the use of the so-called Proper Generalized Decomposition -PGD -widely described in our recent works [START_REF] Chinesta | An overview of the proper generalized decomposition with applications in computational rheology[END_REF][START_REF] Chinesta | A. PGD-based computa-tional vademecum for efficient design, optimization and control[END_REF][START_REF] Chinesta | The Proper Generalized Decomposition for Advanced Numerical Simulations[END_REF][START_REF] Chinesta | Model order reduction. Encyclopedia of computational mechanics[END_REF]. The main idea of the PGD technique, when focusing on a given microstructure ω(X), consists in introducing the parameters as extra-coordinates and then solving the resulting multidimensional problem by using a separated representation of the problem solution for circumventing the curse of dimensionality that high-dimensional models entail. Within the PGD framework, the microscopic parametric velocity field reads:

v(x ∈ ω(X), V) ≈ Q n=1 R n (x) • S n (V), ( 6
)
where R n and S n are a priori unknown functions and the symbol "•" denotes the Hadamard product.

In the present paper, we extend to some of the richer rheologies encountered in industrial applications the above methodology that was successfully used in our previous works for modelling the flow of Newtonian and generalized Newtonian fluids in complex porous microstructures.

In the field of structural composites, there is a recent trend to substitute thermosets by thermoplastics. The high shear viscosity of thermoplastics, however, makes preform impregnation difficult. Two solutions are currently envisaged: use of thermoplastics with low molecular weight, and injection of a low molecular weight pre-polymer that polymerizes after having impregnated the porous medium. In both cases, the shear viscosity remains reasonably moderate to allow for preform impregnation.

The flow around fibres entails a sequence of converging and diverging flows that are likely to activate significant extensional effects. Neither the Newtonian nor the generalized Newtonian fluid models are capable of describing the behaviour of polymeric liquids in extensional flows. Consideration of richer constitutive laws seems compulsory. In the present paper, we follow the simplest route and adopt the constitutive equation of a so-called quasi-Newtonian fluid. This allows us to model realistic behaviours both in shear and extensional flows while neglecting history-dependent viscoelastic phenomena.

In Section Quasi-Newtonian fluid model, we define the quasi-Newtonian fluid model used in the present work. Then, we describe in Section Upscaling of Quasi-Newtonian behaviour the upscaling procedure that allows one to go from the microscopic scale at which the quasi-Newtonian model applies, to the macroscopic scale characterized by an effective Darcy law. Finally, the predicted permeabilities are discussed in Section Numerical results, and compared to those obtained with a simple generalized Newtonian fluid.

Quasi-Newtonian fluid model

The quasi-Newtonian fluid [START_REF] Ryssel | Comparison of a quasi-Newtonian fluid with a viscoelastic fluid in planar contraction flow[END_REF] differs from the generalized Newtonian fluid in the fact that the only material function involved, namely the viscosity function, is not merely a function of the second invariant of the rate of strain tensor but also depends on the relative rate of rotation of the fluid. Consequently, the quasi-Newtonian fluid is able to show shear-thinning in shear flow and extension-thickening in elongational flow. In what follows, and without loss of generality, we consider 2D planar flows.

The constitutive equation for the quasi-Newtonian fluid reads

σ = -pI + 2η QN D, (7)
where η QN is the effective viscosity of the fluid that accounts for shear-thinning as well as extension-thickening according to the local type of flow. Thus, the first step is to define a descriptor of the local flow type. There is a large variety of such descriptors [START_REF] Thompson | Persistence of straining and flow classification[END_REF]. In the present work, we use the particular descriptor considered in [START_REF] Ryssel | Comparison of a quasi-Newtonian fluid with a viscoelastic fluid in planar contraction flow[END_REF] and successfully validated there with the experimental results reported in [START_REF] Quinzani | Birefringence and laser-Doppler velocimetry (LDV) studies of viscoelastic flow through a planar contraction[END_REF].

First, we introduce the relative rate of rotation W from

W = ω -w, (8)
where ω is related to the flow vorticity (2 = ∇v -(∇v) T ) according to

= • ω, ( 9
)
where is the third-order permutation tensor (also known as the Levi-Civita tensor). The vector ω can also be written in terms of the curl of the velocity,

ω = - 1 2 ∇ × v. ( 10
)
The vector w in Eq. 8 represents the angular velocity of the eigenvectors of the rate of strain tensor D. Thus, if the two eigenvectors of D have their orientations expressed by the unit vectors d 1 and d 2 , then we have

w = d i × Dd i Dt , ( 11
)
for i = 1 or i = 2, and where D• Dt is the material derivative. A simple local descriptor of the type of flow can then be constructed from the second invariant of D, γ = √ 2D : D, and from the norm of W, W , according to

χ = 4 W γ + 2 W . ( 12
)
We have 0 ≤ χ ≤ 2, and more concretely

χ = ⎧ ⎨ ⎩ 0, in planar extension 1, in pure shear 2, in rigid motion . ( 13
)
Now, the shear viscosity η S depends on γ , i.e. η S ( γ ), whereas the extensional viscosity η E depends on ˙ , with 2 ˙ = γ (in 2D), i.e. η E (˙ ). As in [START_REF] Ryssel | Comparison of a quasi-Newtonian fluid with a viscoelastic fluid in planar contraction flow[END_REF], we consider the geometrical weighting of both viscosities, affected by a function f (χ) depending on the local flow type:

η QN ( γ , χ) = (η S ( γ )) f (χ) (η E (˙ )) 1-f (χ) . ( 14
)
Here, the function f (χ) satisfies the following conditions

f (χ) = 1, if χ = 1 0, if χ = 0 , ( 15
)
in order to recover the shear viscosity at locations exhibiting a shear flow and the extensional viscosity where planar extension occurs. We consider f (1 < χ ≤ 2) = 1 for approaching the zero shear rate viscosity in the limit case of rigid rotation.

Upscaling of Quasi-Newtonian behaviour

The upscaling procedure for Newtonian fluids, which allows one to define a bridge between the microscopic Stokes model and the macroscopic Darcy law, is summarized in the Appendix. It is adapted in [START_REF] Lopez | Flow modeling of linear and nonlinear fluids in two and three scale fibrous fabrics[END_REF] and [START_REF] Lopez | Flow modeling of linear and nonlinear fluids in two scale fibrous fabrics[END_REF] to handle a generalized Newtonian fluid. In what follows, we extend the procedure to the quasi-Newtonian fluid. The quasi-Newtonian fluid flow model is solved in the representative volume ω(X), where two phases coexist, i.e. the fluid phase occupying the domain ω f (X) and the solid phase, assumed rigid and at rest, occupying the region ω s (X), with ω f (X) ∪ ω s (X) = ω(X) and ω f (X) ∩ ω s (X) = ∅. The flow model consists of the mass and momentum balances and the quasi-Newtonian constitutive equation discussed in the previous section:

⎧ ⎨ ⎩ ∇ • σ = 0 ∇ • v = 0 σ = -pI + 2η QN D . ( 16
)
It is complemented with the boundary condition v(x ∈ ∂ω(X)) = V, where the macroscopic velocity V is considered as an extra-coordinate in the PGD procedure. Thus, one finally obtains the solution for any possible macroscopic velocity in a certain domain , i.e. V ∈ , where is the ball (or circular disk in 2D) whose radius is set to the maximum expected velocity magnitude. It is assumed that the solid phase in ω(X) has an empty intersection with ∂ω(X), i.e. the entire boundary ∂ω(X) of ω(X) is occupied by the fluid phase. If that intersection is not empty, adequate It is important to note that the flow solver considered in this work (i.e. a standard stable Q2/P 1 finite element Stokes solver for a non-uniform fluid viscosity), produces identical results when applied to the test problems considered in [START_REF] Ryssel | Comparison of a quasi-Newtonian fluid with a viscoelastic fluid in planar contraction flow[END_REF].

Once the non-linear flow problem ( 16) is solved for any macroscopic velocity V ∈ , the parametric fluid velocity v(x, V) is available. From it, γ , ˙ = γ /2 and χ can be evaluated everywhere in the fluid domain x ∈ ω f (X) and for each macroscopic velocity V ∈ : we thus obtain γ (x, V), ˙ (x, V) and χ(x, V).

Then, the parametric quasi-Newtonian viscosity η QN (x, V) can be computed at each point in the fluid domain for each macroscopic velocity. After freezing η QN (x, V) at each point x ∈ ω f (X), the following two parametric linear Stokes problems are solved in ω(X):

∇p i (x, V) = ∇ • 2η QN (x, V)D i (v i (x, V)) ∇ • v i (x, V) = 0 , i = 1, 2 (17)
with the boundary conditions (in the 2D case)

v 1 (x ∈ ∂ω(X), V) = (1, 0) T v 2 (x ∈ ∂ω(X), V) = (0, 1) T . ( 18
)
This yields the localization tensor, from which the parametric effective permeability K QN (V; X) can be identified, as explained in the Appendix for the case of linear fluids.

Numerical results

In order to apply the upscaling procedure described above, we calculate first the microscopic velocity v(x, V) associated to each possible value of the macroscopic velocity V by means of the Proper Generalized Decomposition in order to alleviate the computational cost related to the multidimensional character of the problem. With a frozen parametric quasi-Newtonian viscosity η QN (x, V), the flow problems (17) become linear and their solution yields the parametric permeability K QN (V; X). For additional details on the constructor, the interested reader can refer to [START_REF] Lopez | Flow modeling of linear and nonlinear fluids in two scale fibrous fabrics[END_REF]. In this section, we consider the representative volume depicted in Fig. 1. It is slightly asymmetric in order to induce richer kinematics.

Following the discussion in [START_REF] Ryssel | Comparison of a quasi-Newtonian fluid with a viscoelastic fluid in planar contraction flow[END_REF], the quasi-Newtonian fluid rheology is selected such as to approximate the response of a Giesekus model with relaxation time λ, mobility α, solvent viscosity η s and polymeric contribution η p to the zero shear rate shear viscosity. The Giesekus model has a shear-thinning shear viscosity and an extension-thickening planar elongational viscosity for 0 < α < 0.5. Thus, we specify

η S ( γ ) = η p β S ( γ ) + η s , ( 19
)
Fig. 5 Quasi-Newtonian viscosity η QN associated with the velocity field v(x) related to the macroscopic velocity field V T = (V cos θ, V sin θ), V = 5 mm • s -1 and θ = π/3 Fig. 6 Solution v 1 (x, V) of Eq. 17 associated with the quasi-Newtonian viscosity related to the macroscopic velocity field

V T = (V cos θ, V sin θ), V = 5 mm • s -1 and θ = π/3
and

η E (˙ ) = η p β E (˙ ) + η s ; (20)
with

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ β S ( γ ) = y(y-1+2α) 2α(y 2 +(λ γ ) 2 ) y 2 = 1 2 1 + √ 1 + ν ν = 16α(1 -α)(λ γ ) 2 , ( 21
)
Fig. 7 Solution v 2 (x, V) of Eq. 17 associated with the quasi-Newtonian viscosity related to the macroscopic velocity field

V T = (V cos θ, V sin θ), V = 5 mm • s -1 and θ = π/3 and ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ β E (˙ ) = y-1+2α 2αy y 2 = 1 + 1 2 (z + 1) 2 + ν -(z + 1) ν = 64α(1 -α)(λ ˙ ) 2 z = -4(λ ˙ ) 2 . ( 22
)
As in [START_REF] Ryssel | Comparison of a quasi-Newtonian fluid with a viscoelastic fluid in planar contraction flow[END_REF], the function f (χ) in Eq. 14 is defined as follows:

f (χ) = 3 sin 4 (π/2)χ 1+2 sin 4 (π/2)χ if χ ≤ 1 1 i f1< χ ≤ 2 . ( 23
)
In the simulations, we set η p = 1.422 Pa • s, η s = 0.002 Pa • s, α = 0.15 and λ = 0.03 s. The resulting fluid rheology is depicted in Fig. 2, showing indeed significant shear-thinning as well as extension-thickening.

The quasi-Newtonian flow problem is solved for any possible macroscopic velocity V ∈ . The parametric domain consists of velocities having any direction in [0, 2π) and magnitude [0, V max ], with V max = 15 mm • s -1 . One of these solutions, associated with the macroscopic velocity V of magnitude V = V = 5 mm • s -1 oriented at π/3 with respect to the horizontal axis, is depicted in Fig. 3.

Figure 4 depicts the local descriptor related to the type of flow for this particular choice of the macroscopic velocity, that results from the particularization of the parametric solution χ(x, V). The parametric quasi-Newtonian viscosity η QN (x, V) is shown in Fig. 5, where extensional effects are noticeable.

According to the proposed upscaling procedure, the computed quasi-Newtonian viscosity is then frozen when 8 (bottom) and obtained via the same procedure but by enforcing f (χ) = 1. We see that extensional effects have a significant influence at high flow rates, as they reduce the permeability to about half of that of the purely shear-thinning fluid. There is not, however, any noticeable effect on the principal directions of the permeability tensor.

In order to analyze the influence of the fluid rheology, we now specify η p = 1 Pa•s and η s = 1 Pa•s, all other rheological parameters being left unchanged (Fig. 9). In view of the much larger (constant) solvent viscosity, shear-thinning is almost suppressed. Shear-thinning effects act mainly in the regions where χ ≈ 1, i.e. in the neighbourhood of the fluidsolid interfaces. Thus, one would expect a significant reduction of the permeability in comparison with the previous results. Moreover, a lower difference is expected between the quasi-Newtonian and Generalized-Newtonian fluids. Figure 10 compares both effective permeabilities, in particular K QN 11 versus K GN 11 , confirming the expected behaviour.

Conclusions

In this paper, we addressed the definition of a permeability tensor for quasi-Newtonian fluids flowing in two-scale porous media. The upscaling procedure is based on the existence of a localization tensor and the equality of dissipated powers at the different scales.

For linear Newtonian fluids, an effective geometrical permeability tensor can be defined at the macroscopic scale that only depends on the microscopic geometrical features. This permeability tensor can be computed offline (for a given microstructure) by solving two boundary value problems in the 2D case. For non-linear fluids such as the Generalized Newtonian and quasi-Newtonian fluids, the permeability is no longer purely geometrical.

In order to speed up the macroscopic calculation, we proposed to construct, using the PGD, parametric solutions of microscopic flow problems defined in a representative volume with a given microstructure. Once the permeability is expressed with a parametric dependence on the macroscopic velocity, its particularization becomes computationally very fast, making it possible to solve macroscopic flow problems with a computational complexity that scales with the number of iterations needed by the non-linear macroscopic solution procedure.

Numerical results have been obtained for a quasi-Newtonian model that mimics the rheological behaviour of a Giesekus viscoelastic model, showing significant shearthinning and extension-thickening. The predicted effective permeability tensor was found to be greatly affected by the extensional viscosity.

Compliance with Ethical Standards

Fig. 1 Fig. 2

 12 Fig. 1 Representative volume ω

Fig. 3 3 Fig. 4

 334 Fig. 3 Microscopic velocity field v(x) associated with the macroscopic velocity field V T = (V cos θ, V sin θ), V = 5 mm • s -1 and θ = π/3

Fig. 8 Fig. 9

 89 Fig. 8 Parametric effective permeability for the quasi-Newtonian (top) fluid: K QN 11 (V) (top-left) and K QN 22 (V) (top-right) and Generalized Newtonian fluid (bottom): K GN 11 (V) (bottom-left) and K GN 22 (V) (bottom-right)

Fig. 10

 10 Fig. 10 Parametric effective permeability for the quasi-Newtonian and Generalized Newtonian fluids: K QN 11 (V) (left) and K GN 11 (V) (right)

Conflict of interests

The authors declare that they have no conflict of interest.

Appendix Upscaling of linear Newtonian behaviour: from Stokes to Darcy

Without loss of generality, suppose a statistically homogeneous porous medium where two characteristic scales coexist, the coarsest one in and the finest in ω, as depicted in Fig. 11. For the sake of simplicity, we assume the same microstructure everywhere in , which allows one to consider a unique microscopic domain ω. The fluid is Newtonian and incompressible, with a constant shear viscosity η.

At the macroscopic level, the flow is described by Darcy's law and the continuity equation, where K(X) is the macroscopic permeability tensor at position X. In the light of the aforementioned assumption, K(X) is supposed to be the same everywhere, that is K(X) ≡ K.

At the microscopic scale, inclusions occupy the solid region ω s assumed rigid and at rest, i.e. v(x ∈ ω s ) = 0, separated by the fluid domain ω f where the Newtonian fluid flows. To be consistent with the definition of Darcy's permeability (also called intrinsic or geometric permeability), it is assumed that the fluid occupies the entire volume between solid inclusions.

In [START_REF] Lopez | Flow modeling of linear and nonlinear fluids in two and three scale fibrous fabrics[END_REF], the postulated macroscopic Darcy model was upscaled from the microscopic Stokes model by assuming that:

-There exists a localization tensor M(x, X) such that, when the macroscopic velocity V(X) is prescribed at the boundary ∂ω(X) of ω(X), the microscopic velocity v(x) at each point x ∈ ω(X) solution of the resulting Stokes problem in ω(x), can be obtained from

-The power dissipated at both scales must be the same for any velocity V(X). As shown in [START_REF] Lopez | Flow modeling of linear and nonlinear fluids in two and three scale fibrous fabrics[END_REF], this yields

In view of macroscopic Darcy's law, we thus have

Introducing Eq. 25 in the left-hand side of Eq. 27 and developing the contraction product yields:

with

with i = 1,2,3 and k = 1,2,3, in the general 3D case. Thus, as soon as the localization tensor M is determined, from the macroscopic velocity V we can compute the microscopic velocity v(x) at each position x within the representative volume ω(X), and then evaluate the average of their derivatives according to Eq. 27 in order to identify K -1 and from it K.

For calculating the localization tensor M(x, X) in the general 3D case, we solve three Stokes flow problems in the microscopic domain ω(X):

for three different boundary conditions on ∂ω(X), which is considered fully-contained in the fluid phase, i.e. ∂ω(X) ∩ ω s = ∅:

These are compatible with the incompressibility constraint. Thus, the localization tensor M(x, X) is finally obtain as


Tiago De Oliveira
email: tiago.de-oliveira@cnrs-orleans.fr

Régis Guégan
email: regis.guegan@univ-orleans.fr

Thomas Thiebault
email: thomas.thiebault@cnrs-orleans.fr

Claude Le Milbeau
email: claude.le-milbeau@univ-orleans.fr

Fabrice Muller
email: fabrice.muller@univ-orleans.fr

Vinicius Teixeira
email: vsteixeira@ucs.br

Marcelo Giovanela

Mohammed Boussafir
email: mohammed.boussafir@univ-orleans.fr

Adsorption of diclofenac onto organoclays: Effects of surfactant and environmental (pH and temperature) conditions

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

of organoclays for the adsorption of a recalcitrant PPs in numerous aquatic compartments that can be used as a complement with activated carbon for waste water treatment.

Introduction

Despite the presence of inorganic pollutants such as trace elements and their well-known toxic effects, the priority of the water contaminants has turned recently to emerging micro-pollutants such as pharmaceutical products (PPs) which are more and more frequently found in various aquatic compartments. Among PPs that are non-biodegradable and recalcitrant to conventional wastewater treatments, diclofenac, a non-steroidal anti-inflammatory drug, is one of the few compounds that show harmless actions on both aquatic and terrestrial ecosystems [START_REF] Verlicchi | How efficient are constructed wetlands in removing pharmaceuticals from untreated and treated urban wastewaters? A review[END_REF][START_REF] Coquery | Synthèse des travaux du projet AMPERES[END_REF][START_REF] Sadezky | Proposal of an environmental indicator and classification system of pharmaceutical product residues for environmental management KNAPPE project[END_REF][START_REF] Miege | Fate of pharmaceuticals and personal care products in wastewater treatment plants -Conception of a database and first results[END_REF]. The occurrence of diclofenac in surface, groundwater and even drinking water results to its high consumption in modern countries that stands to 200-2300 µg -1 day -1 habitant -1 [START_REF] Schwaiger | Toxic effects of the non-steroidal anti-inflammatory drug diclofenac Part 1: histopathological alterations and bioaccumulation in rainbow trout[END_REF][START_REF] Kasprzyk-Hordern | Illicit drugs and pharmaceuticals in the environment -Forensic applications of environmental data. Part 1: Estimation of the usage of drugs in local communities[END_REF] which leads thus to discharge with concentration up to 1000 ng L -1 in wastewater plants influents [START_REF] Miege | Fate of pharmaceuticals and personal care products in wastewater treatment plants -Conception of a database and first results[END_REF].

Adsorption represents an effective and simple technique to remove both inorganic and organic micropollutants. Clay minerals were for a long time known for their outstanding adsorption properties.

However, although being used in drinking water treatment, these layered materials turn out to be ineffective for the sequestration of persistent nonpolar hydrophobic contaminants and/or negatively charged organic species such as diclofenac. The chemical modification of clay mineral through the use of cationic surfactant generates organoclay adsorbents that combine both properties of the inorganic layered material and hydrophobic environment with the intercalation of the organic cation. Thus, these materials have shown particular properties for the adsorption of many emerging pollutants and were found to be more effective in some cases than activated carbon. Recent studies have highlighted an improvement of the adsorption capacity using organoclays as adsorbents that can potentially be an alternative for the remediation of micro-pollutants [START_REF] Younker | Impact of salinity and dispersed oil on adsorption of dissolved aromatic hydrocarbons by activated carbon and organoclay[END_REF][START_REF] Lee | Removal of anionic metals by amino-organoclay for water treatment[END_REF][START_REF] Balme | Structure, orientation and stability of lysozyme confined in layered materials[END_REF][START_REF] Lepoitevin | BSA and Lysozyme adsorption on Montmorillonite: influence of the interlayer cation[END_REF][START_REF] Guegan | Confinement effects on water structure in membrane lyotropic phases[END_REF][START_REF] Guegan | Self-assembly of a non-ionic surfactant onto a clay mineral for the preparation of hybrid layered materials[END_REF][START_REF] Guégan | Sandwich organization of non-ionic surfactant liquid crystalline phases as induced by large inorganic K4Nb6O17 nanosheets[END_REF][START_REF] Guegan | Nonionic organoclay: A 'Swiss Army knife' for the adsorption of organic micro-pollutants?[END_REF][START_REF] Zhu | Sorption characteristics of CTMA-bentonite complexes as controlled by surfactant packing density[END_REF][START_REF] Lee | Sorption of hydrophobic organic compounds onto organoclays[END_REF][START_REF] Costa Filho | Preparação e caracterização de organobentonita modificada com ionenos alifáticos[END_REF]. Indeed, these materials adsorb a large range of pollutants, such as pesticides [START_REF] Cruz-Guzman | Adsorption of the herbicide simazine by montmorillonite modified with natural organic cations[END_REF][START_REF] Rodriguez-Cruz | Modification of clay barriers with a cationic surfactant to improve the retention of pesticides in soils[END_REF], herbicides [START_REF] Carrizosa | Organosmectites as sorbent and carrier of the herbicide bentazone[END_REF], phenolic compounds [START_REF] Yapar | Effect of adsorbent concentration to the adsorption of phenol on hexadecyl trimethyl ammonium-bentonite[END_REF][START_REF] Yilmaz | Removal of phenol using TDTAB-and HDTABbentonites[END_REF] and various pharmaceutical products [START_REF] Zhang | Sorption and desorption of carbamazepine from water by smectite clays[END_REF][START_REF] Akcay | The equilibrium and kinetics studies of flurbiprofen adsorption onto tetrabutylammonium montmorillonite (TBAM)[END_REF]. The adsorption efficiency of organoclays mainly depends on both the chemical nature and the structural organization of the intercalated surfactants [25][START_REF] Park | Application of organoclays for the adsorption of recalcitrant organic molecules from aqueous media[END_REF]. Thus, surfactants showing a long alkyl chains such as hexadecyltrimethylammonium (HDTMA) create an appropriate organic environment within the inorganic frame for the adsorption of alkanes whereas modifiers such as benzyl decyltrimethylammonium (BDTA) or benzyltrimethylammonium (BTMA) owing aromatic rings show an excellent affinity with aromatic compounds [25][START_REF] Park | Application of organoclays for the adsorption of recalcitrant organic molecules from aqueous media[END_REF]. Numerous works focused on the improvement of the adsorption of organic contaminants by using cationic organoclays where the the density of surfactant reached an amount of maximum of 1 time that of the cation exchange capacity (CEC) of clay mineral allowing one to create a hydrophobic environment without any strong steric effect that may restrict the sequestration of pollutants. Thus, such organoclays with 1 CEC surfactant density are revealed to be an appropriate sorbent and were proposed as realistic material for the treatment of water. While the incorporation of surfactant at high concentration (i.e. > 1CEC) generates a large hydrophobic network with an arrangement of the organic cations in bilayer within the interlayer space that may enhance or reduce the adsorption properties, rare studies focused on the sorption of organic compounds onto organoclays prepared at such high surfactant density.

Thus, this work stresses out the adsorption of diclofenac onto two organoclays prepared with two different long alkyl tail cationic surfactants: the BDTA and the HDTMA. These two organic cations showing different chemical groups, the BDTA displaying an aromatic ring that may contribute to the adsorption driving forces, were adsorbed at a concentration of 4 times that of the cation exchange capacity (CEC) of the clay mineral where adsorbents display a large interlayer space and a consequent organic environment, which may favor the further adsorption of diclofenac. The extent of the adsorption and the identification of the main interactional mechanisms were determined through the combined information obtained by a set of complementary techniques, including X-ray diffraction (XRD), elemental analyses (CHNS), gas chromatography coupled with mass spectrometry (GC/MS) and Fourier transform infrared spectroscopy (FTIR). The latter spectroscopic technique gives critical information of possible re-arrangement of surfactants and allows us to corroborate the adsorption data, but also to follow, through an original way, the density of both diclofenac and the organic cations. In addition to the discussion of the nature of the cationic surfactant for the adsorption of diclofenac and its characterization with the use of FTIR, this study also aimed to evaluate the removal of diclofenac for different pH and temperature experimental conditions for a proper determination of the adsorption isotherms and thermodynamics of this micro-pollutant onto the organoclay adsorbents.

Materials and Methods

Adsorbents

Wyoming sodium montmorillonite (Mt), obtained from the Source Clay Minerals Repository, University of Missouri (Columbia, MO) was used in this study as a starting material. The structural formula can be expressed as: (Ca 0.12 Na 0.32 K 0.05 )[Al 3.01 Fe(III) 0.41 Mn 0.01 Mg 0.54 Ti 0.02 ][Si 7,98 Al 0.02 ]O 20 (OH) 4 .

Mt was fractioned to < 2 μm by gravity sedimentation, purified and Na + exchanged [START_REF] Guegan | Nonionic organoclay: A 'Swiss Army knife' for the adsorption of organic micro-pollutants?[END_REF]. This Mt clay mineral shows a cation exchange capacity (CEC) of 76.4 meq per 100 g clay.

Both benzyldimethyltetradecyl ammonium (BDTA) and hexadecyltrimethylammonium (HDTMA) cationic surfactants were supplied from Sigma Aldrich Chemical. Aqueous solutions of BDTA and HDTMA with Millipore deionized water were prepared at ambient temperature for concentration at 4 times the cation exchange capacity (CEC) of Mt clay. These surfactant solutions were poured in aqueous Mt dispersions before being stirred at 300 rpm for 24 h. At the end of this procedure, the solutions were centrifuged at 3000 rpm for 20 min and the supernatants were removed. The resulting organoclays are dried at 100°C for 24 h before analysis.

Experimental techniques

The concentration diclofenac before and after being in contact with the adsorbents was obtained by GC/MS analysis using Trace GC Ultra chromatograph and TSQ Quantum XLS mass spectrometer (Thermo Scientific). Moreover, elemental analyses (Carbon, Nitrogen) on powder form samples were performed by using a Thermo Scientific flash 2000 analyzer. The concentration diclofenac before and after being in contact with the adsorbents was obtained by GC/MS analysis after a derivatisation of diclofenac with a mixture of MTBSTFA/pyridine (2:1, v/v). Analyses were performed on a Trace GC Ultra gas chromatograph (GC) fitted with a Thermo Trace Gold TG-5 MS capillary column (60 m, 0.25 mm i.d., 0.25 μm film thickness) and coupled to a TSQ Quantum XLS mass spectrometer equipped with an AS 3000 autosampler (both from Thermo Scientific) with the same methods as described by elsewhere [START_REF] Thiebault | Clayey-sand filter for the removal of pharmaceuticals: Percolation experiments[END_REF]: the temperature of the column was held at 50°C for 3 min, increased from 50 to 120°C at 30°C.min -1 , and from 120 to 310°C at 3°C.min -1 with a final isothermal hold at 310°C for 21 min. 2 μL of sample was injected in splitless mode at 280°C. Helium was the carrier gas (1 mL.min -1 ). The mass spectrometer was operated in EI mode at 70eV, from m/z 50 to 500.

Fourier transform infrared (FTIR) measurements in the range 650-4000 cm -1 , were recorded using a Thermo Nicolet 6700 FT spectrometer equipped with a Deuterated Triglycine Sulfate (DTGS) detector and a Nicolet Continuum microscope. The powder samples were spread over a NaCl window of the microscope. The analyzed sample area was a square of side 100 m chosen under the microscope 15X Infinity Reflechromat objective. The analyses were performed in transmission mode and each spectrum corresponded to the average of 256 scans collected at 2 cm -1 resolution.

The d 001 spacing's of the layered materials (organoclays and diclofenac/organoclay complexes) were determined by the first 00l reflection from the X-rays patterns which were recorded in a conventional θθ Bragg-Brentano configuration by using a Thermo Electron ARL'XTRA diffractometer equipped with a Cu anode (CuK α = 1.5418 Å) coupled with a Si(Li) solid detector. The diffractograms on dry samples (100°C for 24 h) were performed between 1 and 24° (2θ) with an angular and time steps of 0.04° and 10 s, respectively.

Adsorption of diclofenac

Diclofenac (2-[(2,6-dichlorophenyl)amino]benzeneacetic acid of which pKa is 4.1) purchased from Sigma Aldrich Chemical, and assumed to have a purity > 98 % was used under sodium salt form which shows a solubility better than 10 g L -1 . Batch adsorption experiments of diclofenac onto organoclays were conducted in duplicate using at least 10 initial aqueous solutions ranging from 10 mg L -1 to 2 g L -1 .

The solid to liquid ratio was kept constant, where 100 mg of organoclays were used for 50 mL of diclofenac aqueous solutions in centrifuge tubes and the pH was adjusted by a mix of HCl acidic and NaOH basic solutions. Samples were shaken on a rotary shaker at 50 rpm during 24 h in order to reach the equilibrium final concentrations and the centrifuged at 5000 rpm for 25 min. Both supernatants and sorbents after contact with the diclofenac were removed and analyzed through GC/MS. The amount of adsorbed diclofenac was calculated by the difference between the initial and equilibrium final concentrations that allowed us to determine the adsorption isotherms which were corroborated with CHNS analyses. Before their XRD and FTIR characterizations, the entire resulting organoclays and Mt after adsorption of diclofenac were dried at 100°C for 48 h for CHNS, XRD and FTIR sorbents characterizations.

Adsorption isotherm

The resulting adsorption isotherms of the diclofenac onto organoclays and Mt were fitted using Langmuir model expressed as (1) [START_REF] Levan | Binary Langmuir and Freundlich isotherms for ideal adsorbed solutions[END_REF]:

q 𝑒 = 𝑞 𝑚𝑎𝑥 𝐾 𝐿 𝐶 𝑒 1 + 𝐾 𝐿 𝐶 𝑒
where q e is the equilibrium diclofenac amount adsorbed on organoclays (mol g -1 ), C e the equilibrium diclofenac concentration in the resulting solution (mol L -1 ), q max the maximun adsorption capacity of the sorbents (mol g -1 ) and K L the Langmuir adsorption constant which is related to the thermodynamic equation of Gibbs free energy by relation (2):

∆G°= -RT ln K L

where R is the universal gas constant (8.214 J mol -1 K -1 ).However, the adsorption of both surfactants and diclofenac contributes to a change of the chemical nature of the surface as well as a heterogeneity in the distribution of the adsorption sites, and therefore Freundlich and Dubinin-Radushkevish (D-R) model isotherm equations were used. The Freundlich isotherm is an empirical equation to describe heterogeneous systems of which the linear form is (3) [START_REF] Levan | Binary Langmuir and Freundlich isotherms for ideal adsorbed solutions[END_REF][START_REF] Ozcan | Adsorption of Acid Blue 193 from aqueous solutions onto BTMA-bentonite[END_REF]:

ln q 𝑒 = ln 𝐾 𝐹 + 1 𝑛 ln C 𝑒
where K F (L g -1 ) and n are both adsorption Freundlich parameters which respectively indicate the adsorption efficiency and the non-linearity degree between diclofenac and sorbents. The Freundlich parameter n suggest when 1/ n is between 0.1 and 1 a favorable adsorption process [START_REF] Liu | Adsorption of methyl violet from aqueous solution by halloysite nanotubes[END_REF]. The D-R

(3) isotherm allow to acquire complementary thermodynamic parameters. The equation is expressed as (4) [START_REF] Dubinin | Equation of the Characteristic Curve of Activated Charcoal[END_REF][START_REF] Hobson | Physical adsorption isotherms extending from ultrahigh vacuum to vapor pressure[END_REF]:

ln q e = ln q m -βε 2 where ε is the Polayi potential, q m the theoretical saturation capacity of the sorbent and β the constant related to the activity coefficient (mol² J -²) connected to the free energy E of adsorption (kJ mol -1 ) with equation ( 5):

E = 1 √2β
The latter parameter gives information about adsorption mechanism. Indeed, if E < 8 kJ.mol -1 , the adsorption process follows preferentially molecular interactional mechanism, while for E > 8 kJ.mol -1 ion-exchange is envisaged [START_REF] Onyango | Adsorption equilibrium modeling and solution chemistry dependence of fluoride removal from water by trivalent-cation-exchanged zeolite F-9[END_REF].

Results and Discussion

BDTA-Mt and HDTMA-Mt organoclays properties

The BDTA-Mt and HDTMA-Mt organoclays were characterized by both XRD (Fig. 1a) and FTIR (Fig. 1b) analyses. The 00l diffraction patterns of both cationic organoclays shift to lower of 2θ angular values underlining the intercalation of surfactant within the interlayer space of the layered starting clay mineral. Indeed, the dehydrated sodium exchanged clay mineral shows a 00l reflection at 9.2 ° (2θ) leading to an interlayer spacing of 9.7 Å which is in agreement with previous works [START_REF] Lepoitevin | BSA and Lysozyme adsorption on Montmorillonite: influence of the interlayer cation[END_REF][START_REF] Thiebault | Adsorption mechanisms of emerging micropollutants with a clay mineral: Case of tramadol and doxepine pharmaceutical products[END_REF]. As expected, with a concentration of 4 times the CEC for the cationic surfactants, the expansion of the interlayer space reaches 34.5 Å for BDTA-Mt and for HDTMA-Mt organoclays, that matches an organization of the organic cations in normal bilayer or paraffin structure. Indeed, FTIR spectra confirm such surfactant organization and show typical features characteristics of the organic compounds, such as the absorption bands at 2840-2920 cm -1 relative to the symmetric and antisymmetric CH 2 stretching vibrations of the surfactant alkyl chains (Fig. 1b) [START_REF] Lepoitevin | BSA and Lysozyme adsorption on Montmorillonite: influence of the interlayer cation[END_REF][START_REF] Guégan | Sandwich organization of non-ionic surfactant liquid crystalline phases as induced by large inorganic K4Nb6O17 nanosheets[END_REF].

(4)

(5) It is well established that the wavenumbers of the CH 2 stretching bands of hydrocarbon chains are extremely sensitive to the conformational ordering and change in the gauche-trans conformer ratio of the chains which can be used as probe, in correlation with the d 001 spacing variation, for the surfactant arrangement within the silicate layers. Here, the wavenumbers of both symmetric and antisymmetric CH 2 stretching vibrations indicate that the organic cations located in the internal structure show an alltrans conformation [START_REF] Guégan | Sandwich organization of non-ionic surfactant liquid crystalline phases as induced by large inorganic K4Nb6O17 nanosheets[END_REF][START_REF] Guegan | Nonionic organoclay: A 'Swiss Army knife' for the adsorption of organic micro-pollutants?[END_REF].

Effect of the surfactant

The final adsorption properties of organoclay materials mainly depend on the chemical nature of the adsorbed surfactant. The intercalation of cationic surfactants, through ion exchange with the inorganic cations, switches the chemical nature of the starting layered material from hydrophilic to hydrophobic while expands at a wide opening the interlayer space in order to the adsorption of organic compounds.

Nevertheless, the improvement of the sequestration mainly depends on both the chemical nature and the structural organization of the intercalated surfactants. Thus, surfactants showing a long alkyl chains such as HDTMA create an appropriate organic environment within the inorganic frame for the adsorption of alkanes whereas modifiers, such as BDTA show an excellent affinity with aromatic compounds [START_REF] Guegan | Nonionic organoclay: A 'Swiss Army knife' for the adsorption of organic micro-pollutants?[END_REF].

Here, since diclofenac consists of two aromatic rings and display a negative charge in aqueous solution for pH above its pKa, thus it is expected the nature of the intercalated surfactant plays on the magnitude of the adsorption. The equilibrium adsorption isotherms of this PP onto both HDTMA-Mt and BDTA-Mt organoclays performed at room temperature and for a free pH (i.e. pH above the pKa of the diclofenac) display similar lineshapes with two main regimes: (i) a gradual growth of the adsorbed amount of diclofenac as its equilibrium concentration is increased; and (ii) a steady state at high concentration (Fig. 2).

The use of Langmuir equation model allows one to quantify the degree of affinity between the diclofenac and the two adsorbents, which here appears to be identical. Nevertheless, at high concentration of diclofenac, BDTA-Mt organoclay seems to be more efficient with a maximum adsorbed amount 1.9 x 10 -4 mol g -1 whereas HDTMA-Mt shows a plateau with a maximum adsorbed amount of 1.55 x 10 -4 mol g -1 . This slight difference (15%) of adsorbed amount may result to the nature of surfactant where BDTA could interact through molecular interactions (π-π) with diclofenac leading to enhancement of the adsorption. Nevertheless, while the adsorption and intercalation of surfactant give a hydrophobic surface, it decreases the specific surface area of the starting clay mineral. The contribution to the reduction of the accessible surface mainly depends on the nature of the inserted surfactant. This also has to be taken into account and may be another cause for the difference of maximum adsorbed amount of diclofenac onto organoclays.

Equilibrium Diclofenac Concentration (mol.L The equilibrium adsorption isotherms highlight the proper affinity of diclofenac to both organoclays.

XRD patterns of the organoclays (Fig. S1 in supplementary information) after being in contact with diclofenac do not show any difference with d 001 spacing's remaining at the same openings. However, the value of adsorbed diclofenac in the steady state regime reaches an amount that exceeds the maximum possible adsorbed onto the external surface of organoclay platelets, stressing out that diclofenac is intercalated within the interlayer space of both organoclays and should disrupt the arrangement of the organic cations. FTIR spectroscopy gives us important information on the conformation of organic cations by focusing on the CH 2 stretching bands wavenumbers windows in the range 2800-3000 cm -1 (Fig. 3). In addition to determining any conformational change, this quantitative technique (in transmission mode) allows one to follow the evolution of both surfactants and diclofenac amount. Thus, FTIR data allows one to extract similar information as adsorption isotherms obtained through GC/MS, but also to analyze the fate of surfactants. Indeed, the integrated intensity of the CH 2 stretching bands is based on the density of surfactants as well as that of the CO 2 -stretching in the range of 1450-1650 cm -1 follows that of diclofenac. The evolution of the integrated intensity CO 2 -stretching bands related to the diclofenac density shows a lineshape similar to that of the resulting GC/MS adsorption isotherms, confirming the proper adsorption of the PP onto both organoclays (Fig. 4 shows the results for HDTMA-Mt, but similar behavior is obtained for BDTA-Mt, see Fig. S2 in supplementary information).

In contrast, the decrease of the integrated intensity of both CH 2 stretching bands for the organic cations in organoclay as the concentration of diclofenac is increased reveals a reorganization of the surfactant within the interlayer space due to the adsorption of PP which is highlighted with the shift to high frequencies of both symmetric and antisymmetric CH 2 stretching bands (Fig. 5). Nevertheless, the reorganization of the surfactant aggregates depends on the adsorbed amount of diclofenac. In low concentration regime, diclofenac may be adsorbed onto the external surface of organoclay platelets and thus does not or slightly perturb the all-trans conformation and paraffin structure of the organic cations within the interlayer space. Once the starting concentration of diclofenac reaches 0.5 g L -1 , the organic cations undergo a deep rearrangement in lateral layers where the alkyl chains show more gauche conformer (i.e. disorder) and are less packed together within the interlayer space.

Starting Diclofenac Concentration (g.L

Thermodynamics

The temperature may drive to antagonist situations. Indeed, while temperature increases mobility of organic (i.e. reduction of mean free path) compounds making the adsorption or contact with a surface easier, it also reduces the energy of activation of numerous processes. Nevertheless, several experimental studies on the adsorption of different types of pollutants onto clay minerals or soils pointed out an increase of temperature reduces the percentage of adsorption [START_REF] Akcay | The equilibrium and kinetics studies of flurbiprofen adsorption onto tetrabutylammonium montmorillonite (TBAM)[END_REF]. The equilibrium isotherms at 5, 25, 35 and 50 °C resulting of the adsorption of diclofenac onto BDTA-Mt and HDTMA-Mt organoclays were fitted with Langmuir, Freundlich and D-R model equations. The whole adsorption parameters derived from and quality of the fitting procedure (R 2 ) are listed in Table 1.

Moreover, we used an error function (F error ) in order to evaluate which equation models was best suited to describe these processes. A lower result from the error function indicated a smaller difference between adsorption capacity calculated by the model (q cal ) and the experimental (q exp ). F error can be expressed according to the following equation ( 6):

F error = ∑( q i cal -q i exp q i exp p i ) ²
where q i cal is the value of adsorbed diclofenac amount predicted by the fitted model, q i exp is the value of adsorbed diclofenac amount measured experimentally, i indicated the initial diclofenac concentration of the experiment and p is the number of experiments performed.

Equilibrium Diclofenac Concentration (mol.L

1 Langmuir model is based on the assumption of monolayer adsorption on a structurally homogenous adsorbent, where all sorption sites are identical and energy equivalent whereas the empirical Freundlich equation is applicable to adsorption on heterogeneous surfaces, where the interaction between the adsorbed molecules is not limited to the formation of a monolayer. Finally, the D-R equation isotherm represents a more general model for adsorption, since it assumes a heterogeneous surface and variability for the adsorption potential and was applied successfully for fitting the adsorption of organic species onto organoclays. Here, the diclofenac adsorption isotherms are properly fitted by the three equation models used as R 2 values show with a magnitude comprised between 0.859 and 0.999 and F error values are between 6.8 x 10 -5 and 0.271 (Table 1). Based on R 2 values, experimental data seem to be better adjusted to the Langmuir model, however its F error are higher than 0.002 and to those for both Freundlich and D-R equations, ranging from 6.8 x 10 -5 to 7.7 x 10 -4 . Although Langmuir equation properly fitted experimental data, the two latter equations appear to be more suitable for modeling the adsorption of diclofenac onto both organoclays surface. Indeed, organoclays shows a heterogeneous surface with a reorganization of the organic cations once diclofenac is adsorbed leading to a distribution of several adsorption sites that are taken into account in both Freundlich and D-R equation models (Fig. S3-4 in supplementary information).

Nevertheless, the quality of the fits obtained by using a Langmuir equation is far enough to work on their parameters (Fig. 6). The whole ∆G° values derived from the Langmuir fitting procedure are negative for both organoclays and decrease when temperature rises ranging from -15.9 kJ mol -1 for BDTA-Mt at 5°C to -19.8 kJ mol -1 at 50°C, suggesting that the adsorption is spontaneous and more favorable at low temperature [START_REF] Shu | Adsorption removal of Congo red from aqueous solution by polyhedral Cu2O nanoparticles: Kinetics, isotherms, thermodynamics and mechanism analysis[END_REF]. In addition, since the thermodynamics parameter derived from the fitting procedure show similar values for both organoclays, it thus reinforces the idea that the nature of surfactant appears to play a minor role in the adsorption of the anionic PP here. The preparation of the organoclays was achieved at a concentration of 4 times the CEC where organic cations, after 1 CEC, are adsorbed through molecular hydrophobic interaction leading to a positive charge of the surfactant aggregates within the interlayer space. Thus, instead of the nature of surfactant it is much more the positive charge of the surfactant and the proportion of accessible mobilized organic cations that can form ionic organic complexes with diclofenac in a first step, which can be adsorbed in a second step onto organoclays, suggesting a chemisorption process. Indeed, if the adsorbed amount increases as the temperature increases, it corresponds to physisorption whereas the opposite trend generally underlines a chemisorption process [START_REF] Sevim | An investigation of the kinetics and thermodynamics of the adsorption of a cationic cobalt porphyrazine onto sepiolite[END_REF]. Here, the adsorbed amount q m decreases as the temperature increases excluding de facto any physisorption mechanism. The increase in temperature enhances the dynamics of the whole organic compounds making easier adsorption but also the mobility of the alkyl-chains of the surfactant cations. Indeed, it is well known that temperature affects alkyl chains with the introduction of gauche conformer or disorder leading to a perturbation of the adsorption sites in the hydrophobic environment grafted onto the clay mineral. If temperature may favor the probability of contact between diclofenac and the accessible organic cations where the driving force involves a chemisorption mechanism, it also disrupts the hydrophobic adsorption sites generated by the grafted surfactant cations onto the clay surface and here it appears that the balance between these two antagonist effects tips in favor of a decrease of the adsorption efficiency at high temperature. Both ∆G° and the estimated adsorbed amount q max as well as the adsorption efficiency more generally decrease as the temperature rises indicating a chemical adsorption of diclofenac onto organoclays. Moreover, the enthalpy change (∆H°) and entropy change (∆S°) can be calculated with the van't Hoff equation ( 7):

∆G°= ∆H°-T∆S°

the octahedral sheets with a release of Al 3+ or Mg 2+ cations that can substitute to the sodium cations, thus to a complete destructuration of the layered structure of the clay mineral. Here, three values of pH were tested where the diclofenac mainly displays a protonated form at the lowest pH value of 3.5 although negative diclofenac are also present in solution whereas for the two other pH > pKa, diclofenac is negative. Interestingly, at a pH of 3.5

where it was expected to get a low adsorption, the equilibrium adsorption isotherm show a lineshape similar to those at high pH (Fig. 7). Moreover, the parameters derived (with G° < 9 kJ mol -1 ) from the Langmuir fit underline that the adsorption thermodynamically spontaneous and favorable. It is likely the proportion of the diclofenac in anions is associated organic cations and follows a similar adsorption mechanism as it was previously explained but it cannot exclude that further protonated molecules could be adsorbed as well onto the adsorbents surface. Nevertheless, the adsorption efficiency is particularly enhanced at high pH value where the fraction of diclofenac is preponderant. Indeed, the K L values for HDTMA-Mt vary from 1.97 x 10 3 to 5.68 x 10 3 L mol -1 for low and high pH respectively confirming the better affinity of diclofenac with adsorbent at high pH values. The same trend is supported with the K F values which varies from 1.31 x 10 -3 L.g -1 (acidic pH) to 3.26 x 10 -3 L g -1 (basic pH) and the theoretical saturation capacity q m of D-R model. Their coherent evolutions point out the importance of the pH value mainly playing on the charge of diclofenac on the adsorption efficiency of the two organoclays.

Conclusion

The use of both BDTA and HDTMA cationic surfactants with a Na + exchanged Mt for the preparation of organoclays display favorable adsorption for diclofenac defined as an emerging micro-pollutant due to its presence at significant concentrations in numerous natural water environments. Thermodynamics parameters derived from Langmuir, Freundlich and D-R equation fitting procedure suggest that the adsorption onto the two sorbents is mainly driven by strong electrostatic interaction between diclofenac and both surfactants. It appears that the nature of surfactant plays a minor role but is rather the mobilized organic cations, which are main vectors for the adsorption of diclofenac. The association of hydrophobic complexes between diclofenac and surfactant in solution leads to their adsorption/intercalation onto organoclays implying a rearrangement of the surfactant as FTIR spectra stressed out. The magnitude of the adsorption is somehow controlled by the experimental temperature and pH conditions. Indeed, as it was expected, a rise in temperature involves a decrease in the adsorbed amount of diclofenac. The grafting of surfactant prevents any deterioration of the organoclay structure while keeping the arrangement of the surfactant aggregates. In acidic conditions, the proportion of diclofenac showing a negative charge decreases, thus reduces the adsorbed amount onto organoclays. Finally, the prepared organoclays show an adequate surface for the adsorption of diclofenac for various experimental conditions that it is possible to get in analogous natural environments, and thus are revealed as interesting materials in complement to activated carbon for water remediation. However, before any potential application or use as treatment method, further research works need to be undertaken, in particular the desorption of diclofenac with organic solvent and the quantification the adsorption efficiency with reuse organoclays.

Fig. 1 .

 1 Fig. 1. (a) XRD patterns for dehydrated Na-Mt, BDTA-Mt and HDTMA-Mt organoclays. (b) FTIR spectra of Na-Mt, BDTA-Mt and HDTMA-Mt organoclays. The proper intercalation of both surfactants is confirmed by the shift of the first 00l reflection at low angle values (XRD data) and with the CH 2 scissoring and stretching bands of the alkyl tails of the organic cations.

Fig. 2 .

 2 Fig. 2.Effect of the chemical nature of cationic surfactant (BDTA-Mt and HDTMA-Mt) on adsorption of diclofenac at ambient temperature (25°C) and using starting concentration of diclofenac ranging from 10 mg L -1 to 2 g L -1 .

Fig. 3 .

 3 Fig. 3. 3D evolution of the FTIR spectra in the wavenumbers window of CH 2 stretching bands for (a) HDTMA-Mt organoclay in the range 2800-3000 cm -1 and the (b) CO 2 -stretching bands of diclofenac in the range 1450-1650 cm -1 for different starting diclofenac concentrations (0.1, 0.4, 0.6, 1 and 2 g L -1 ).

Fig. 4 .

 4 Fig. 4. Evolution of normalized integrated FTIR intensity of both CO 2 -bands (diclofenac) and CH 2 bands (surfactants) as a function of the starting diclofenac concentration. The lines are guides for the eyes.

Fig. 5 .

 5 Fig. 5. Representation of the change in frequency of the symmetric and asymmetric CH 2 stretching of the HDTMA-Mt and BDTA-Mt organoclays as a function of the starting diclofenac concentration.

Fig. 6 .

 6 Fig. 6. Equilibrium adsorption isotherms of diclofenac onto BDTA-Mt organoclay (prepared at a concentration of 4 times the CEC of the starting Mt) for 4 different temperatures: 5, 25, 35, 50°C. The efficiency of the adsorption is greater at low temperature.

Fig. 7 .

 7 Fig. 7. Equilibrium adsorption isotherms of diclofenac (pKa=4.51) onto HDTMA-Mt organoclay (prepared at a concentration of 4 times the CEC of the starting Mt) for three different pH values of 3.5, 6.0-6.5, and 11.0. The efficiency of the adsorption is enhanced at high pH.

Table 1

 1 Adsorption isotherm constants determined with Langmuir, Freundlich, and Dubinin-Radushkevich model fit for the adsorption of diclofenac onto BDTA-Mt and 2 HDTMA-Mt organoclays for different temperatures and pH experimental conditions.

	3
	4
	Adsorbant			Langmuir				Freundlich				Dubinin-Radushkevich
		T (°C) pH	q max (mol.g -1 )	K L (L.mol -1 )	ΔG° (kJ.mol -1 )	R 2 F error	K F (L.g -1 )	n	R 2	F error	q m (mol.g -1 )	β (mol².J -2 )	E (kJ.mol -1 )	R 2	F error
	BDTA-Mt	5	6.5	3.06 x 10 -4 617	-15.9	0.986 0.065	2.26 x 10 -3 2.34 0.983 1.4 x 10 -5	5.80 x 10 -4 6.18 x 10 -9 8.99	0.978 1.2 x 10 -4
		25	6.5	1.88 x 10 -4 2298	-19.1	0.912 0.271	8.73 x 10 -4 3.35 0.989 3.5 x 10 -5	3.39 x 10 -4 3.69 x 10 -9 11.6	0.996 2.5 x 10 -5
		35	6.5	7.37 x 10 -5 5263	-21.2	0.985 0.002	4.48 x 10 -4 3.49 0.941 2.1 x 10 -5	1.08 x 10 -4 2.07 x 10 -9 15.5	0.980 6.8 x 10 -5
		50	6.5	7.94 x 10 -5 3058	-19.8	0.963 0.016	3.23 x 10 -4 3.87 0.995 4.3 x 10 -4	1.22 x 10 -4 2.34 x 10 -9 14.6	0.999 4.1 x 10 -5
	HDTMA-Mt	5	6.5	3.16 x 10 -4 425	-14.9	0.981 0.263	6.67 x 10 -3 1.60 0.971 1.5 x 10 -5	6.48 x 10 -4 7.46 x 10 -9 8.18	0.982 3.0 x 10 -4
		25	6.5	1.71 x 10 -4 1439	-18.1	0.996 0.011	3.11 x 10 -3 1.87 0.962 3.1 x 10 -5	4.48 x 10 -4 6.45 x 10 -9 8.81	0.991 1.3 x 10 -4
		35	6.5	1.26 x 10 -4 1858	-18.6	0.983 0.056	1.47 x 10 -3 2.25 0.972 2.3 x 10 -5	2.54 x 10 -4 4.88 x 10 -9 10.1	0.962 3.9 x 10 -4
		50	6.5	1.10 x 10 -4 2379	-19.3	0.992 0.008	1.28 x 10 -3 2.23 0.914 3.9 x 10 -4	2.67 x 10 -4 5.60 x 10 -9 9.45	0.966 7.4 x 10 -4
	BDTA-Mt	25	4	1.41 x 10 -4 3062	-19.8	0.917 0.065	5.72 x 10 -4 3.72 0.988 7.7 x 10 -4	1.24 x 10 -4 1.14 x 10 -9 20.9	0.859 2.4 x 10 -4
		25	11	1.98 x 10 -4 2187	-19.1	0.912 0.071	7.33 x 10 -4 3.81 0.989 6.9 x 10 -4	1.77 x 10 -4 1.20 x 10 -9 20.4	0.911 5.1 x 10 -4
	HDTMA-Mt	25	4	1.37 x 10 -4 1976	-18.8	0.949 0.003	1.31 x 10 -3 2.47 0.995 2.9 x 10 -4	2.54 x 10 -4 3.94 x 10 -9 11.2	0.994 6.8 x 10 -4
		25	11	1.76 x 10 -4 5681	-21.4	0.994 0.003	3.26 x 10 -3 2.03 0.980 2.4 x 10 -4	4.31 x 10 -4 4.70 x 10 -9 10.3	0.994 1.1 x 10 -5

Acknowledgements

The authors would like to thank the Région Centre (Projects SMAM & HarPE) as well as the Labex VOLTAIRE for their financial support.

The ∆H° values for BDTA-Mt and HDTMA-Mt are positive and, respectively, equal to 10.4 kJ mol -1 and 11.8 kJ mol -1 , reflecting the endothermic process of the adsorption of diclofenac onto the adsorbents. The positive ∆S° values (97.3 J mol -1 K -1 and 97.9 J mol -1 K -1 for BDTA-Mt and HDTMA-Mt, respectively) suggest an adsorption of diclofenac leads to the randomness in the solid/liquid interface.

The Freundlich constant K F is related to the adsorption capacity of the organoclays and is enhanced with the decrease in temperature, which is coherent with the Langmuir analysis.

Further, the empirical parameter of Freundlich 1/n is between 0.1 and 1 confirming also that the adsorption process is favorable for both organoclays. The efficiency decrease of organoclay to adsorb diclofenac at high temperature may also be understood by the theoretical saturation capacity q m determined with the D-R equation model. Indeed, for BDTA-Mt (but similar feature is observed for HDTMA-Mt), q m exhibits a drastic diminution from 5.80 x 10 -4 mol g -1 to 1.22 x 10 -4 mol g -1 if one changes the temperature from 5°C to 50°C. D-R model stresses out with E > 8 kJ mol -1 (Table 1) in the whole studied pH and temperature range, chemisorption mechanism mainly ensures the adsorption of diclofenac, involving strong electrostatic interactions with the mobilized organic cations whatever their chemical nature giving a certain stability of the resulting composite materials as well as a credit to use organoclays as geosorbent for the adsorption of diclofenac and more generally PPs.

Effect of pH on the adsorption of diclofenac onto BDTA-Mt and HDTMA-Mt organoclays

The efficiency of the adsorption and the magnitude of the parameters derived from the fitting procedure by using Langmuir, Freundlich and D-R equation models depend on numerous parameters, such as the temperature, ionic strength, pH, and both the chemical nature and molecular size of the organic cations. The effect of pH modifies both the density of charge of the chemical form of the diclofenac, but also those of the clay surface and its organoclay derivatives as well as its structure. Indeed, it is well known that low pH leads to an attack of