Giulio Sciarra

Francesco Dell'isola

Kolumban Hutter

A solid-fluid mixture model allowing for solid dilatation under external pressure

 Guarascio 1999).

Introduction

In many engineering applications of binary mixture models of solid-fluid interactions the pore space or the permeabilities are prescribed functions of the spatial coordinates, but they do not evolve together with the temporal changes of the other field variables: the engineering models based on this assumption stem from the studies of Terzaghi [START_REF] Von Terzaghi | Theoretical Soil Mechanics[END_REF]. One example to the contrary is for instance the slow creeping deformation and the percolation of brine through a salt formation from a pressurized cavern filled with a liquid (for a detailed discussion of these phenomena see for instance [5]). Observations indicate that an increase in cavern pressure will not only result in a very slow creeping deformation of the salt dome but equally enhance the diffusion of liquid through the salt. This increase in percolation is not only due to an increase of the partial pressure of the liquid as a reaction to an increase of the cavern pressure but also because of grain boundary openings in the vicinity of the cavern walls.

Many of the presently existing mixture theories treating solid-fluid interactions may in principle be able to cope with the fluid dilatation mechanism: Bowen [START_REF] Bowen | Incompressible porous media models by use the theory of mixtures[END_REF], [START_REF] Bowen | Compressible porous media models by use of the theory of mixtures[END_REF], Müller [START_REF] Müller | Thermodynamics[END_REF], Morland [START_REF] Morland | A simple constitutive theory for a fluid saturated porous solid[END_REF], Rajagopal & Tao [START_REF] Rajagopal | On boundary conditions for a certain class of problems in mixture theory[END_REF], Svensdsen & Hutter [START_REF] Svendsen | On the thermodynamics of a mixture of isotropic materials with constraints[END_REF]; scrutiny has shown, however, that very particular constitutive behaviour must be assumed to achieve it. Terzaghi and Fillunger were aware that their models were unable to describe this dilatancy (see for instance de Boer & Ehlers [START_REF] De Boer | The development of the concept of effective stress[END_REF]), Svendsen & Hutter's [START_REF] Svendsen | On the thermodynamics of a mixture of isotropic materials with constraints[END_REF] mixture theory allows for it, but as shown by [START_REF] Dell'isola | A qualitative analysis of the dynamics of a sheared and pressurized layer of saturated soil[END_REF], [START_REF] Dell'isola | Variations of porosity in a sheared pressurized layer of saturated soil induced by vertical drainage of water[END_REF] rather singular functional forms of the thermodynamic free energies ψ(•, n) as functions of the porosity n are required if an appreciable space opening of the pores is to be achieved. It was then thought that introducing density gradients as independent constitutive variables would regularize the formulation [START_REF] Dell'isola | A variational approach for the deformation of a saturated porous solid. A second gradient theory extending Terzaghi's effective stress principle[END_REF], and indeed it did so.

However, in the problem stated above and in the mentioned papers we were confronted with a further difficulty: the flux boundary conditions between the single constituent body and the mixture. As shown by Hutter et al. [START_REF] Hutter | On interfacial transition conditions in two phase gravity flow[END_REF] the jump condition of momentum requires non-vanishing momentum production terms to exist on the singular surface (see also Svendsen & Gray [21]). An alternative, simpler than this, and not requiring surface balances is to postulate a parametrization how the traction on the single-constituent side of the surface is distributed between the tractions on the mixture-side [18]. This parametrization can be expressed as a scalar parameter (akin to surface fraction) depending on a number of variables, say the constituent densities. Obviously, the parametrization must be such that the work done by these boundary tractions in a simply connected closed circuit in this parameter space is zero for otherwise a perpetuum mobile of the first kind would result (Seppecher, personal communication). The construction of the potential for the boundary tractions restricts the parametrization in our case to the extent that the functional dependence is fixed except for a single constant.

A mixture of two constituents one of which is a fluid can only exist in equilibrium when it is confined, i.e., when a pre-stress is exerted on it; this means that reference states with non-vanishing fluid densities are always pre-stressed in such media. This fact gives rise to the question how the free energy describing the interior behaviour of the mixture must be structured that perturbations about such pre-stressed states are stable. For a quadratic dependence of the free energy upon the constituent densities the analysis shows that the stability condition depends strongly on the coupling term involving the two densities. For isotropic stresses and a homogeneous reference state we will prove that there is always an upper bound of the external pressure beyond which such states become unstable, corroborating Fillunger and Terzaghi's explosion under pressure (de Boer [START_REF] De Boer | Theory of porous media[END_REF]).

One of the significant results determined in this paper is the fact that, depending on the coefficients of the parabolic representation of the free energy, a stable perturbation of a homogeneous pre-stressed reference state can give rise to a decrease in the apparent solid density with an increase of the external pressures. If the true density of the fluid is essentially constant this corresponds to a dilatation of the space occupied by the fluid. This is essentially how a sponge responds to the absorption of water from the outside; however, this property is exactly what is needed to achieve an increase of permeability without phase changes of the salt in the salt cavern problem mentioned at the beginning. The spatially dependent pore opening close to the cavern space is manifest as a boundary layer (in what is called Disturbed Rock Zone) and for this reason we close this paper with a presentation of the higher gradient model corresponding to the analysis of the earlier sections. A one-dimensional linear-elastic problem is finally presented as a first application of the newly introduced model: we prove that close to unstable pre-stressed reference configurations the thickness of the boundary layer (where the apparent mass densities are not constant) at the external interface of the solid-fluid mixture tends to infinity. The regularizing properties of the proposed second gradient model are thus established.

The ultimate problem will have to incorporate visco-plastic components to also account for the creeping deformation of the salt (see Cosenza [5]).

Theoretical concepts

Consider a binary mixture of a solid matrix with connected pores which are filled with a liquid. This arrangement can be thought of as a soil, rock or sponge. Let the two components be referred to as the solid and the fluid and indicate them by the suffices s and f . Let, moreover ρ s , ρ f and v s , v f be the solid and the fluid apparent densities and velocities, respectively, in the mixture. The mixture density and the barycentric velocity are then given by

ρ = ρ s + ρ f ⇒ 1 = ρ s ρ + ρ f ρ =: ξ s + ξ f , (1)
v = ξ s v s + ξ f v f , ( 2
)
in which ξ s and ξ f are the mass fractions of the solid and the fluid, respectively. We conceive this mixture to be non-reactive so that the balances of mass for the constituents reduce to

∂ρ a ∂t + div (ρ a v a ) = 0, (a = s, f ) . ( 3
)
In the ensuing analysis we shall restrict ourselves to purely mechanical processes; temperature will play no role, and so the constituent momentum equations are the only additional dynamical equations to be added to [START_REF] Casal | La Capillarité interne[END_REF]. Instead of a direct application of these laws we shall use the Principle of Virtual Power (PVP) applied to the appropriate energy functional to derive them. Let ψ int be the energy volume density and, as we limit ourselves to spherical states, assume it to be depending on the densities ρ a : ψ int = ψ int (ρ s , ρ f ). Other dependencies could be incorporated, but will not be here for simplicity. The PVP states that the variation of the total energy in the considered body (a mixture) related to its (barycentric) motion equals the power of the external forces. If the exterior of the mixture body is a fluid then the boundary traction exerted on the mixture is a pressure, p ext which must be distributed between the constituents via the phenomenological ansatz

p f = d f p ext , p s = d s p ext , d s + d f = 1, (4)
in which d s and d f are surface fraction parameters which, like the energy itself, depend upon the densities ρ s and ρ f . This parametrization is of constitutive nature and has been proposed already in the past (Morland [START_REF] Morland | A simple constitutive theory for a fluid saturated porous solid[END_REF], Rajagopal & Tao [START_REF] Rajagopal | On boundary conditions for a certain class of problems in mixture theory[END_REF]). Generally the areal fractions are identified with the volume fractions. This identification traces back in mining engineering to the law of Delesse (1848), as quoted by de Boer [START_REF] De Boer | Theory of porous media[END_REF]. In this paper we demonstrate the interrelation of the two via an argument concerning the first law of therodynamics. Thus the PVP yields the statement

d dt B ψ int dV = - ∂B d a p ext n • v a dA, (5)
where summation over doubly repeated indices a is understood and n denotes the outward pointing unit normal to ∂B which is the boundary of the Eulerean region B , the actual placement of the mixture body. We assume the external pressure p ext to be a conservative field. This means that the work performed by p ext on a cyclic quasi-static variation of the state parameters ρ s and ρ f must be path independent for otherwise a perpetuum mobile of the first kind would emerge. This is avoided if a potential ψ ext exists such that1

d dt B ψ ext dV = - ∂B d a p ext n • v a dA. ( 6
)
Employing the Reynolds Transport Theorem on the LHS, the divergence theorem on the RHS of (6) and using (3) yields

∂B -ρ a ∂ψ ext ∂ρ a + ξ a ψ ext + d a p ext n • v a dA + B ρ a ∇ ∂ψ ext ∂ρ a • v a dV = 0. (7)
Since [START_REF] De Boer | The development of the concept of effective stress[END_REF] must hold for all velocity fields v a , using a localization procedure the following restrictions on constitutive equations emerge:

∂ψ ext ∂ρ a = 1 ρ ψ ext + d a ρ a p ext , ( 8
)
ρ a ∇ ∂ψ ext ∂ρ a = 0 =⇒ ∂ψ ext ∂ρ a = p ext c a , (9)
where p ext c α are constants and p ext is inserted for convenience. Each of these statements constitutes two equations and from their exploitation ψ ext and d a can be determined. They can be regarded as a vectorial statement in the 2D state-variable space spanned by ρ s and ρ f . The second one implies that (c s , c f and c being integration constants)

ψ ext (ρ s , ρ f ) = p ext c s ρ s + c f ρ f + c , ( 10
)
so that using this in ( 8) we obtain

d s = ξ s + ρ s [ρ] (1 -ξ s ) = ξ s 1 + ρ f [ρ] , c = -1, ( 11
)
where

c s -c f =: 1 [ρ] . ( 12
)
Let us pause to interpret this result: the requirement that the external pressure does not perform work along a closed trajectory in the state space ρ s , ρ f has led to a restriction on the constitutive equations for d s (and d f ), which divide the external pressure p ext between the partial pressures p s and p f , respectively. The division p a = ξ a p ext is an obvious and allowable possibility ([ρ] → ∞), but it is not exhaustive. Any choice obeying ( 8) and ( 9) must have the form [START_REF] Dell'isola | Variations of porosity in a sheared pressurized layer of saturated soil induced by vertical drainage of water[END_REF] in which the quantity ρ s / [ρ] is a scale parameter controlling how the boundary pressure is divided between the constituents.

Recall that 0 < d s < 1 and 0 < d f < 1 which imply

- 1 ρ f < 1 [ρ] < 1 ρ s . ( 13
)
Evidently the typical scale parameter for the division of the pressure p ext between the constituent pressures involves an inverse density. Obviously the mixture density in the reference configuration is a permissible choice for [ρ], in principle it may be infinitely large; d s then lies within the interval 0 < d s < 1 and is a linear function of ξ s between these two values.

Observe that the condition d a = 0 implies that either ρ a = 0 or ρ-ρ a = [ρ]. The second of these conditions is physically impossible as both ρ a and ρ must be considered as variables. This means that when one part of the external pressure p ext acting on the fluid (or on the solid) is equal to zero then the fluid-solid mixture at the boundary reduces to the solid skeleton (or to a pure fluid) only.

Linearization in the neighbourhood of a pre-stressed reference state

Let B 0 be the mixture body in its reference state with positions labeled by X. Constituent densities and the external pressure in this state are denoted by ρ 0 a (a = s, f ) and p ext 0 , respectively; they will be assumed to be different from zero, and thus will give rise to constituent stress distributions p 0 a . Such a nontrivial pre-stressed situation is a necessary requirement for the mixture material introduced in Sect. 2. Indeed, the model is only meaningful when ρ s > 0 and ρ f > 0; thus d s > 0 and d f > 0, and this necessarily implies, through the parametrization of the constituent pressures that p 0 s = 0 and p 0 f = 0, if p ext 0 = 0. For the fluid to exist at all, a nonvanishing pressure is mandatory. So only nontrivial pre-stressed conditions are possible. Consider a pre-stressed reference configuration of the solid and suppose the displacement of the material particles of the solid constituent to be small so that geometric linearizations are justified throughout. Reference values of field quantities will carry the suffix (•) 0 , perturbed fields will be denoted by a tilde, (•); therefore, in view of [START_REF] Dell'isola | A qualitative analysis of the dynamics of a sheared and pressurized layer of saturated soil[END_REF] we have

ψ ext = p ext 0 c s ρ 0 s + c f ρ 0 f -1 +p ext 0 c s ρs + c f ρf + pext c s ρ 0 s + c f ρ 0 f -1 (14)
+ pext c s ρs + c f ρf .

Alternatively, we express the potential ψ int as a quadratic form about the reference state ρ 0 s , ρ 0 f ; this will be done in the form ψ int = ρϕ int , in which ϕ int is the internal energy per unit mass. What obtains is as follows:

ψ int = ρϕ int ρ s , ρ f (ρ 0 + ρ) γ s ρs + γ f ρf + 1 2 γ ss ρ2 s + γ sf ρs ρf + 1 2 γ ff ρ2 f + ... ( 15
) = ρ 0 γ s ρs + γ f ρf + 1 2 (2γ s + ρ 0 γ ss ) ρ2 s + 1 2 2γ f + ρ 0 γ ff ρ2 f + γ s + γ f + ρ 0 γ sf ρs ρf ,
where an arbitrary constant in ϕ int is irrelevant and where the constants

γ s := ∂ϕ int ∂ρ s ρ 0 s, ρ 0 f , γ f := ∂ϕ int ∂ρ f ρ 0 s, ρ 0 f , γ ss := ∂ 2 ϕ int ∂ρ 2 s ρ 0 s, ρ 0 f , γ ff := ∂ 2 ϕ int ∂ρ 2 f ρ 0 s, ρ 0 f , γ sf := ∂ 2 ϕ int ∂ρ s ∂ρ f ρ 0 s, ρ 0 f ( 16
)
are supposed known when the reference state is known. Depending upon the value of γ sf relative to the values of γ ss , γ ff coupling will be called weak or strong. Generally growing γ sf increases this coupling. Next, we wish to write down the equilibrium equations which follow from the localization of

d dt B ψ tot dV = d dt B ψ int -ψ ext dV = 0. ( 17
)
For static conditions this yields

ρ s ∂ψ int ∂ ρs -ξ s ψ int = ρ s ∂ψ ext ∂ ρs -ξ s ψ ext , ρ f ∂ψ int ∂ ρf -ξ f ψ int = ρ f ∂ψ ext ∂ ρf -ξ f ψ ext . ( 18
)
Using ( 14) and ( 15) and expanding as illustrated above yields the following zeroth and first order problems: zeroth order

γ s = p ext 0 ρ 2 0 + ρ 0 f ρ 2 0 c s -c f p ext 0 , γ f = p ext 0 ρ 2 0 - ρ 0 s ρ 2 0 c s -c f p ext 0 . ( 19
)
first order

R 11 R 12 R 21 R 22 ρs ρf =      ξ 0 s + ρ 0 s ρ 0 f ρ 0 c s -c f ξ 0 f - ρ 0 s ρ 0 f ρ 0 c s -c f      pext , ( 20
)
where

R 11 := ρ 0 s (2γ s + ρ 0 γ ss ) + ρ 0 f γ s -p ext 0 ξ 0 f 2 c s -c f + ξ 0 f ρ0 , R 12 := ρ 0 s γ s + γ f + ρ 0 γ sf -ρ 0 s γ f -p ext 0 ξ 0 s 2 c s -c f - ξ 0 s ρ0 , R 21 := ρ 0 f γ s + γ f + ρ 0 γ sf -ρ 0 f γ s + p ext 0 ξ 0 f 2 c s -c f + ξ 0 f ρ0 , R 22 := ρ 0 f 2γ f + ρ 0 γ ff + ρ 0 s γ f + p ext 0 ξ 0 s 2 c s -c f - ξ 0 s ρ0 . ( 21
)
Equations ( 19) relate γ s and γ f to one another. On the other hand, (20) could be inverted to obtain ρs and ρf in terms of pext

ρs ρf = [R] -1      ξ 0 s + ρ 0 s ρ 0 f ρ 0 c s -c f ξ 0 f - ρ 0 s ρ 0 f ρ 0 c s -c f      pext , ( 22
)
however, this inversion is only possible if the matrix [R] is invertible. We will show later that conditions of invertibility agree with the requirements of the stability of the pre-stressed reference state. This is the problem we shall address now.

Linear stability analysis of pre-stressed reference states

The scope of this section is not to embed the equilibrium properties of our system into a full non linear dynamic stability analysis, we simply wish to know whether a particular given reference state is stable with respect to small perturbations. In this spirit, the reference states described by the density fields ρ 0 s , ρ 0 f and the pressure p ext 0 will now be assumed to be spatially uniform. Our interest is in their stability against linear perturbations of the densities and the pressures. It is expected that such stability properties depend on the functional form of the total energy

ψ tot = ψ int -ψ ext (23)
and its properties around an equilibrium state. In a static situation stability of the equilibrium state requires the function ψ tot to assume its minimum in equilibrium (Dirichlet criterion), so that the Hessian matrix H ψ tot is positive definite in a neighbourhood of the equilibrium,

H = 2γ s + ρ 0 γ ss γ s + γ f + ρ 0 γ sf γ s + γ f + ρ 0 γ sf 2γ f + ρ 0 γ ff , ( 24
)
implying the Rouse-Hurwitz criteria

2γ s + ρ 0 γ ss > 0, det H > 0. ( 25
)
With the help of (19) 1 the first inequality implies

2 p ext 0 ρ 2 0 + ρ 0 f ρ 2 0 c s -c f >-1/ρ 0 f p ext 0 + ρ 0 γ ss > 0,
requiring at worst that γ ss > 0. We will therefore suppose that γ ss > 0 for all cases. The second inequality for stability can be written as

-ρ 2 0 c s -c f 2 p ext 0 2 + 2ρ 3 0 β p ext 0 + ρ 6 0 γ ss γ ff -γ 2 sf > 0, (26)
in which

β := β 0 + c s -c f β 1 , β 0 := γ ss + γ ff -2γ sf , ( 27
)
β 1 := ρ 0 s γ sf -γ ss + ρ 0 f γ ff -γ sf .
The LHS of (26) is a quadratic form P 2 p ext 0 . It represents a set of parabolas (see Fig. 1) 2 with a positive value at the vertex and which are open in the downward direction; the two solutions of the equation P 2 p ext 0 * = 0,

p ext 0 * = - 1 c s -c f 2 -ρ 0 β ± ρ 2 0 β 2 + ρ 4 0 c s -c f 2 γ ss γ ff -γ 2 sf , ( 28
)
are positive and negative irrespective of whether β > 0 or β < 0. The stability region is p ext 0 < p ext 0 C , where p ext 0 C is the root of (28) with the negative square root sign (as the other root is negative for all admissible choices of the involved parameters). For csc f = 0 the two parabolas reduce to the same pair of straight lines. The one with positive slope has an unbounded stability limit, while that with negative slope has a finite stability limit. These limits are denoted by p * c When c sc f is zero or [ρ] → ∞ then inequality (26) reduces to a linear statement in p ext 0 and the parabolas become straight lines. For β 0 > 0 the stability limit for p ext 0 is unbounded, while for β 0 < 0 it is bounded and is given by

p ext 0 C = ρ 3 0 γ ss γ ff -γ 2 sf 2 2γ sf -γ ss -γ ff . ( 29
)
2 Figure 1 Noticing that

β 0 > 0 ⇒ γ ss + γ ff 2 > γ sf , β 0 < 0 ⇒ γ ss + γ ff 2 < γ sf (30)
we deduce that the first case corresponds to a weak coupling of the solid and the fluid phases via γ sf , whilst the second one is related to a strong coupling. So stability exists for all pressures in the presence of weak coupling, whilst for strong coupling stability is restricted to small pre-stresses. These conditions change qualitatively, when c sc f is not equal to zero. The stability limit is now bounded in both cases, β > 0 and β < 0. Under these circumstances instability is always reached if only p ext 0 is sufficiently large. To be more exact equation (28) implies that 3

p ext 0 C = 1 c s -c f 2 ρ 0 β + ρ 2 0 β 2 + ρ 4 0 c s -c f 2 γ ss γ ff -γ 2 sf ( 31
β| B := 1 -ρ 0 s c s -c f B ≥0 γ ss -γ sf + 1 + ρ 0 f c s -c f B ≥0 γ ff -γ sf . ( 32
)
Thus, β| B is a weighted average of γ ss -γ sf and γ ff -γ sf with weights which depend on the division of the external traction onto the solid and the fluid phases. Moreover, β| B is positive (negative) according to whether γ ss -γ sf and γ ff -γ sf are positive (negative). Since strong coupling corresponds to large |γ sf | values, it is evident that it enhances the potential of instability. Defining

I = β| B ρ 0 c s -c f B γ ss γ ff -γ 2 sf ( 33
)
for negative β| B , (31) takes the form

p ext 0 C (cs-cf ) B = ρ 0 β| B c s -c f 2 1 + 1 + 1/I 2 .
Therefore, the larger I 2 , or the more negative I is, the closer to zero will be p ext 0 C (cs-cf ) B . Consequently I can be taken as a measure of instability.

Conditions for pressure induced dilatancy of the solid matrix

Notice that our intention is to find pressurized conditions that yield dilatancy. It is plain, that such states make only sense if they are stable. This was the reason why the stability analysis was presented in the first place. In this spirit, let us return to the system (20) determining the perturbations ρs and ρf , if pext is given.

A formal inversion of this system of equations yields

ρs = det R ρs det [R] = det R ρs ρ 0 s ρ 0 f det [H] , ( 34
) ρf = det R ρf det [R] = det R ρf ρ 0 s ρ 0 f det [H]
, 3 We explicitly remark that when csc f tends to zero expression (31) tends to infinity if β > 0 but to (29) if β < 0.

where

det R ρs = pext ρ 0 s ρ 0 f ρ 3 0 p ext 0 ρ 0 c s -c f 1 -c s -c f ρ 0 s (35) +ρ 3 0 1 + c s -c f ρ 0 f γ ff -ρ 3 0 1 -c s -c f ρ 0 s γ sf , det R ρf = pext ρ 0 s ρ 0 f ρ 3 0 -p ext 0 ρ 0 c s -c f 1 + c s -c f ρ 0 f +ρ 3 0 1 -c s -c f ρ 0 s γ ss -ρ 3 0 1 + c s -c f ρ 0 f γ sf .
These can be computed by using Cramer's rule implemented in maple. For stability det [H] > 0, so the signs of ρs , ρf are dictated by det R ρs and det R ρf , respectively.

To obtain conditions of dilatancy of the solid matrix, induced by an increase of pressure applied on its boundary, one can remark that ρs is negative if and only if det R ρs < 0, and this requires that

1. if c s -c f > 0, p ext 0 < ρ 2 0 1 -c s -c f ρ 0 s γ sf -1 + c s -c f ρ 0 f γ ff c s -c f 1 -c s -c f ρ 0 s , ( 36
) 2. if c s -c f < 0, p ext 0 > ρ 2 0 1 + c s -c f ρ 0 f γ ff -1 -c s -c f ρ 0 s γ sf |c s -c f | 1 -c s -c f ρ 0 s , (37) 3
. if c s -c f = 0, (35) 1 implies that det R ρs is independent of the pressure p ext 0 . Then det R ρs (cs-cf )=0 = pext ρ 0 s ρ 0 f (γ ff -γ sf ), (38)
and, therefore, ρs is negative provided that

γ ff < γ sf . ( 39
)
If the RHS of (36) is negative then a solution with ρs < 0 does not exist for p ext 0 > 0. Should the RHS of (37) be negative, then dilatancy occurs for all p ext 0 > 0. The foregoing analysis shows that dilatancy under external pressure is possible in a material if only the coupling coefficient γ sf is sufficiently large. This, however, does not yet demonstrate that a real material exists such that dilatancy is indeed established.

Second gradient energy describing pore micro-deformations

The model equations derived so far enjoy the following properties: when a uniform external pressure is applied to the mixture body, exhibiting purely spherical stress states, the constituent densities are equally uniform. This is so because the model does not account for the possible formation of a boundary layer (in which apparent densities are spatially variable) at the mixture external interface. From a physical point of view one could state that this absence is a consequence of a lacking description of the microscopic pore deformation. This is a singular behaviour of the first gradient theory.

To cure the first gradient model from such singular features we now develop a second gradient mixture model. Such theories were preveously developped as non-simple mixture models (Müller [START_REF] Müller | Thermodynamics[END_REF], Rajagopal & Tao [START_REF] Rajagopal | On boundary conditions for a certain class of problems in mixture theory[END_REF]). We are not aware of any second gradient PVP-approach for mixtures in the spirit pursued here, but our approach essentially follows Gouin [START_REF] Gouin | Variational methods for fluid mixtures of grade n Continuum Models and Discrete Systems[END_REF], Casal [START_REF] Casal | La Capillarité interne[END_REF] and Seppecher [START_REF] Seppecher | Etude des condition aux limites en théorie du second gradient: cas de la capillarité C[END_REF] (for simple mixtures) and Germain [START_REF] Germain | La méthode des puissances virtuelles en mécanique des milieux continus[END_REF] who use the PVP for single constituent bodies.

General balance equations

Let us now consider the PVP for the general situation in which (i) body forces b a are present, (ii) the internal energy ψ int may also depend on the second deformation gradients and (iii) the action of the exterior to the body is given by tractions t a and double forces τ a (see for more details [START_REF] Germain | La méthode des puissances virtuelles en mécanique des milieux continus[END_REF] and [START_REF] Dell'isola | A variational approach for the deformation of a saturated porous solid. A second gradient theory extending Terzaghi's effective stress principle[END_REF]). For such a case equation ( 5) is generalized in the form

d dt B ψ int dV = B b a • v a dV + ∂B t a • v a + τ a • ∂v a ∂n dA; ( 40
)
τ s and τ f are the double forces acting on the solid and on the fluid, respectively. We select the simplest gradient dependence of ψ int ,

ψ int = ε (ρ a ) + λ s 2 f ss , f ss := ∇ρ s • ∇ρ s , ( 41
)
where λ s is a constant. This corresponds to a gradient dependence for the solid but not for the fluid, which is special. The localization of (40) with λ s = 0 is easily shown to be

∇p a -m a = b a , in B , p a = d a p ext , on ∂B , ( 42
)
where

p a := ρ a ∂ε ∂ρ a -ξ a ε, m a := ∂ε ∂ρ a ∇ρ a -∇ (ξ a ε) , (no sum. over a) (43)
The details of the derivation follow eqs.( 7) (8) ( 9) and are e.g. also given in [START_REF] Dell'isola | A variational approach for the deformation of a saturated porous solid. A second gradient theory extending Terzaghi's effective stress principle[END_REF]. It therefore suffices do deal with the additional gradient dependent term in (40)

I add :=   d dt B ψ int dV   add = λ s 2 d dt B ∇ρ s • ∇ρ s dV = ∂B t aadd • v a + τ a • ∂v a ∂n . ( 44
)
Applying the Reynolds Transport Theorem yields (see Appendix A)

I add = λ s B 1 2 (f ss ) I • ∇ (ξ a v a ) + ∇ρ s • ∂ ∂t ∇ρ s + ∇ ⊗ ∇ρ s (ξ a v a ) dV = -λ s B div 1 2 f ss ξ s I-f ss I-∇ρ s ⊗ ∇ρ s + div (ρ s ∇ρ s ⊗ I) + 1 2 ∇f ss - 1 2 ∇ (ξ s f ss ) • v s dV +λ s ∂B 1 2 f ss ξ s I-f ss I-∇ρ s ⊗ ∇ρ s + div (ρ s ∇ρ s ⊗ I) n • v s + 1 2 f ss ξ f n • v f -(ρ s ∇ρ s ⊗ I) n • ∇v s dA. ( 45
)
Consequently, localization of the complete equation (40) leads to the following boundary value problem:

∇p s -m s -λ s div ρ s ρ s I+ 1 2 f ss I-∇ρ s ⊗ ∇ρ s = b s , ∇p f -m f = b f , in B (46) -p s n + λ s ρ s ρ s I+ 1 2 f ss ξ s I-∇ρ s ⊗ ∇ρ s n = t s , -p f n+ λ s 2 f ss ξ f n = t f , -λ s ρ s (∇ρ s • n) n = τ s , 0 = τ f , on ∂B . ( 47
)
Notice that the gradient effects only enter the field equations of the solid, this obviously because of our restrictive assumption (41). The boundary conditions of the fluid are, however, affected by the gradient terms; they generate an additional pressure. The surface double forces only enter the boundary conditions of the solid constituent, because the free energy does not depend on ∇ρ f .

External action potential

We now address the problem of finding a potential ψ ext for the external action specified on RHS of (40). We are looking for a potential ψ ext that depends on the state parameters ρ s , ρ f , and ∇ρ s , such that

d dt B ψ ext dV = B b a • v a dV + ∂B t a • v a + τ a • ∂v a ∂n dA. ( 48
)
Its existence assures that a cyclic quasi-static variation of these parameters is path-independent. Using the Reynolds Transport Theorem the LHS of (48) becomes

d dt B ψ ext dV = B ∂ ∂t ψ ext + div ψ ext v dV . ( 49
)
Performing the differentiations term by term and using the constituent balances of mass yields

∂ ∂t ψ ext = - ∂ψ ext ∂ρ a ρ a I • ∇v a - ∂ψ ext ∂ρ a ∇ρ a • v a - ∂ψ ext ∂∇ρ s • ∇ρ s I • ∇v s -∇ ⊗ ∇ρ s ∂ψ ext ∂∇ρ s • v s + ∂ψ ext ∂∇ρ s ⊗ ∇ρ s • ∇v s -ρ s ∂ψ ext ∂∇ρ s ⊗ I • ∇ ⊗ ∇v s (50) div ψ ext v = div ψ ext ξ a v a = ξ a ψ ext div v a + v a • ∇ ξ a ψ ext .
Substituting these above allows to write (49) as

d dt B ψ ext dV = B (β a • v a + B a • ∇v a + B s • ∇ ⊗ ∇v s ) dV , ( 51
)
where the quantities β a , B a (a = s, f ) and B s are defined by the following expressions:

β s := - ∂ψ ext ∂ρ s ∇ρ s -(∇ ⊗ ∇ρ s ) ∂ψ ext ∂(∇ρ s ) + ∇(ξ s ψ ext ), β f := - ∂ψ ext ∂ρ f ∇ρ f + ∇(ξ f ψ ext ), B s := - ∂ψ ext ∂ρ s ρ s - ∂ψ ext ∂(∇ρ s ) • ∇ρ s + ξ s ψ ext I- ∂ψ ext ∂(∇ρ s ) ⊗ ∇ρ s , (52)
B f := - ∂ψ ext ∂ρ f ρ f + ξ f ψ ext I, B s := -ρ s ∂ψ ext ∂(∇ρ s ) ⊗ I.
The RHS of (48) and (51) agree with one another if b a , t a , and τ a are given by

b s = β s -div (B s -div B s ), b f = β f -div B f , t s = (B s -div B s )n-div s (B s n), t f = B f n, τ s = (B s n)n, τ f = 0; (53)
in these formulas we consider the following decomposition for the gradient of an n-th order tensor field Ω:

∇Ω = ∇ s Ω + n ⊗ ∂Ω ∂n , (54)
∇ s Ω being the restriction of ∇Ω on ∂B , and we define the surface divergence of an n-th order tensor field Ω, the differential operator div s such that

div s Ω T u = (div s Ω) • u + Ω • ∇ s u, ∀ (n -1) -th order tensor field u (55) and ∂B div s Ω dA = ∂∂B Ων dS , ( 56
)
where ν is the outward normal to ∂∂B the line boundary of ∂B . For smooth surfaces ∂B the integral on the RHS of (56) vanishes. In this case the contact action on the solid and the fluid, the double forces on the solid and the fluid are given by (see Appendix B)

t s = - ∂ψ ext ∂ρ s ρ s + ξ s ψ ext - ∂ψ ext ∂(∇ρ s ) • ∇ρ s + div ρ s ∂ψ ext ∂(∇ρ s ) +ρ s ∂ψ ext ∂(∇ρ s ) • n (tr ∇ s n)- ∂ψ ext ∂(∇ρ s ) • n ∂ρ s ∂n n+ρ s ∇ s ∂ψ ext ∂(∇ρ s ) • n , t f = - ∂ψ ext ∂ρ f ρ f + ξ f ψ ext n, (57)
τ s = -ρ s ∂ψ ext ∂(∇ρ s ) • n n, τ f = 0.
We distinguish the normal and the shear parts of t s and assume that the double force acting on the solid depends linearly on the external pressure: τ s = d D p ext n; this is reasonable, since increasing the pressure increases the pore space and the latter is opened by the action of the double force. In so doing we obtain the following forms of the constitutive relations for the coefficient d D and for the coefficients d a (a = s, f ), appearing in (4) valid in the case of second gradient solid matrices:

d D = - 1 p ext ρ s ∂ψ ext ∂(∇ρ s ) • n , d s = 1 p ext ∂ψ ext ∂ρ s ρ s -ξ s ψ ext -ρ s div ∂ψ ext ∂∇ρ s (58) + ρ s ∂ψ ext ∂∇ρ s • n tr ∇ s n + ∂ψ ext ∂∇ρ s • n ∇ρ s • n , d f = 1 p ext ∂ψ ext ∂ρ f ρ f -ξ f ψ ext .
These formulas simply emerge if one divides τ s , and the components of t a normal to the surface by p ext . Assume that the body forces b s and b f vanish; then equations (53) 1,2 yield the conditions

ρ s ∇ ∂ψ ext ∂ρ s -div ∂ψ ext ∂∇ρ s = 0, ρ f ∇ ∂ψ ext ∂ρ f = 0. ( 59
)
Further investigations will be necessary to generalize the results, found in the §2, about ψ ext implied by the condition d f +d s = 1. We simply remark here that one can find in the subsequent one-dimensional problem, developed as an application of the introduced new model, a form for ψ ext verifying the above constraints.

A one dimensional application

Consider a one-dimensional problem in which the body forces on the solid and on the fluid vanish. Assume constant external pressure and suppose that the derivative of ψ ext with respect to ∇ρ s is constant in B i.e. d dx

∂ψ ext ∂(ρ s,x ) = 0, (60)
where ∂ψ ext /∂(ρ s,x ) := ∂ψ ext /∂∇ρ s • e, where e is the unit vector defining the x direction. With this, eqs.

(57) reduce to

d s p ext = ∂ψ ext ∂ρ s ρ s -ξ s ψ ext + ∂ψ ext ∂(ρ s,x ) d ρ s dx , d f p ext = ∂ψ ext ∂ρ f ρ f -ξ f ψ ext , t shear s = 0, t shear f = 0, τ s = -ρ s ∂ψ ext ∂(ρ s,x ) e, τ f = 0 (61)
and equations (59) become

d dx ∂ψ ext ∂ρ a = 0, a = s, f . ( 62
)
Eqs. ( 60) and (62) imply that ψ ext is a linear function of its arguments ρ a and d ρ s /dx . Thus, with an appropriate normalization, one has

ψ ext = p ext c s ρ s + c f ρ f + k s d ρ s dx -1 , (63)
whilst the constitutive relations defining the coefficients d a (a = s, f ) and d D are

d s = ξ s 1 + (c s -c f )ρ f + ξ f k s d ρ s dx , d f = ξ f 1 -(c s -c f )ρ s -ξ f k s d ρ s dx , ( 64
)
d D = -k s ρ s (n • e),
as easily deducible from (63) and (58) or from (61). Note that constitutive relations for d s and d f differ from relations [START_REF] Dell'isola | Variations of porosity in a sheared pressurized layer of saturated soil induced by vertical drainage of water[END_REF], obtained by a first gradient mixture model, by an additive quantity; this is due to the assumption on the derivative of ψ ext with respect to ∇ρ s : ψ ext simply depends linearly on ρ s and d ρ s /dx . Consider a linearized theory and assume that in the reference configuration the constituent density of the solid (and of the fluid) is not uniform. This hypothesis is necessary to appreciate second gradient effects: if the densities were uniform in the reference configuration, then (47) 3 would imply τ s = 0, i.e., k s = 0, so that there would be no possibility to have a non-vanishing double force on ∂B , acting on the solid skeleton.

We also assume that the coefficient γ s of the linear term of the internal potential energy ε is not constant but a linear function of ρ 0 s ,

γ s = α s ρ 0 s , ( 65
)
where α s is assumed uniform in B ; a justification for this will be given shortly. With these prerequisites we may now perform a perturbation analysis in the vicinity of a pre-stressed reference state with ψ ext given by (63) and ψ int by (41). With an approach entirely analogous to that of §3 we then deduce from the balance laws the following zeroth and first order equations: zeroth order problem

d ρ 0 dx α s ρ 0 s + ρ 0 α s d ρ 0 s dx -λ s d 3 ρ 0 s dx 3 = 0, ρ 0 f d ρ 0 dx γ f = 0, ( 66
) first order problem -λ s d 3 ρs dx 3 + 2α s ρ 0 s + ρ 0 γ ss d ρs dx + α s ρ 0 s + γ f + ρ 0 γ sf d ρf dx + +2α s d ρ 0 s dx ρs + α s d ρ 0 s dx ρf = 0, (67)
2γ f + ρ 0 γ ff d ρf dx + α s ρ 0 s + γ f + ρ 0 γ sf d ρs dx + α s d ρ 0 s dx ρs = 0.
To these ODEs at each perturbation order four boundary conditions must be added: we suppose that the tractions on the solid and on the fluid and the double force are known at x = 0, and we assume that the double force (and therefore d ρ 0 s /dx ) vanishes as x → ∞. The condition (65) corresponds to the idea that the apparent density of the solid (and of the fluid) in the reference state is given by the sum of a constant and an exponentially decreasing term, smaller than zero; if γ s were constant in B then equation (66) 1 would imply that ρ 0 s = const because of the condition as x → ∞. The solution of the zeroth order problem is

ρ 0 s (x ) = C 1 + C 2 exp x x 0 + C 3 exp - x x 0 , ρ 0 (x ) = C 4 ,
where

x 0 := λ s /ρ 0 α s . ( 68
)
This explicitly demonstrates that the double forces are responsable for the exponential decay of ρ 0 s (x ) as one moves away from the mixture surrounding environment, since λ s = 0. The boundary condition at x → ∞ implies that C 2 = 0; the boundary condition on the value of the double force at x = 0 implies that C 3 ≤ 0, so the apparent densities of the two constituents are given by

ρ 0 s (x ) = C 1 -|C 3 | exp - x x 0 , ( 69
)
ρ 0 f (x ) = (C 4 -C 1 ) + |C 3 | exp - x x 0 ;
x 0 is the characteristic decay length of the zeroth order solution. We do not show the explicit expression for C 1 , (C 4 -C 1 ) and |C 3 | which are rather cumbersome; we simply recall that they depend on the external actions and on the interface constitutive parameters (e.g. c sc f ) and can be interpreted respectively as the apparent solid and fluid mass densities far from the mixture-surrounding environment and their maximum variations induced by applied external double forces.

To compute the solution of the first order problem, consider the following non-dimensionalization of the independent and dependent variables:

ξ := x x 0 , r s := ρs C 1 , r f := ρf C 4 -C 1 .
Equations (67) constitute a fourth order differential problem; so it can be expressed as a system of first order differential equations given in the following form:

d Y d ξ = (R + exp(-ξ) A 0 + exp(-2ξ) A 1 ) Y, (70)
where Y, R, A 0 , A 1 are defined by

Y :=     r s r f dr s /d ξ d 2 r s /d ξ 2     , R :=     0 0 1 0 0 0 -a 4 /a 6 0 0 0 0 1 0 0 a 1 -a 3 a 4 0     , A 0 :=     0 0 0 0 -a 5 /a 6 0 a 5 /a 6 0 0 0 0 0 2a 2 -a 3 a 5 a 2 a 6 2(a 3 a 5 -a 2 ) 0     , (71)
A 1 :=     0 0 0 0 0 0 0 0 0 0 0 0 a 2 a 5 0 -a 2 a 5 0     ,
and The matrix A(z ) =z -1 R + A 0 + z A 1 has at most a pole at z = 0 but it is analytic for 0 < |z | < a, a > 0 and the point z = 0 is a singular point of the first kind for the system (73) (see [START_REF] Coddington | Theory of ordinary differential equations[END_REF]), so it fulfills the hypotheses of theorem 3.1, p. 117 and 4.1, p.119 in [START_REF] Coddington | Theory of ordinary differential equations[END_REF]. Therefore, the fundamental matrix of system (73) is represented in terms of a series, convergent in the set 0 < |z | < a,

a 1 := x 2 0 /λ s (2α s C 1 + C 4 γ ss ) , a 2 := α s |C 3 | x 2 0 /λ s , a 3 := x 2 0 /λ s α s C 1 + γ f + C 4 γ sf , a 4 := α s C 1 + γ f + C 4 γ sf / 2γ f + C 4 γ ff , a 5 := α s |C 3 | / 2γ f + C 4 γ ff ,
Φ(z ) = ∞ i =0 Q i z i e (ln z )J ( 74
)
where J is the canonical form4 of R, if and only if R has characteristic roots which do not differ by positive integers, and

RQ 0 = Q 0 J, (75)
Q m+1 [J + (m + 1)I] = RQ m+1 + m k =0 A k Q m-k ,
where, in this context, I is the (4 × 4) unit matrix. So it follows that the fundamental matrix given as a function of ξ is Φ(ξ) = (Q 0 + exp (-ξ) Q 1 + exp (-2ξ) Q 2 + ...) e -ξ J .

(76)

When truncating this series at the first order term the non-dimensional solid and fluid densities are given by (see Fig. 2)

r s = k 2 + k 1 1 ξ 0 exp (-ξ) + 1 + 1 ξ 0 exp - ξ ξ 0 , r f = k 3 -k 1 a 4 a 6 1 ξ 0 exp (-ξ) + 1 + 1 ξ 0 exp - ξ ξ 0 , ( 77
)
where ξ 0 := (a 1a 3 a 4 ) -1/2 , k i are integration constants to be determined by imposing the boundary conditions implied by (61) at the first order; experiments must give information on their values. Considering that C 4 = ρ 0 and in a first gradient theory C 1 = ρ 0 s (compare with (69)) we may rewrite (78) on using (65) as

√ a 1 -a 3 a 4 = x 0 √ λ s det H H 22 , ( 79
)
where H is defined in (24). Thus

ξ 0 = √ λ s x 0 H 22 det H . ( 80
)
The eigenvalues of R are ± √ a 1a 3 a 4 = ±ξ -2 0 , 0, 0 . In the stability regime one has det H > 0, and H 22 > 0; therefore, ξ 0 is real valued. The fundamental solution of (73) corresponding to the negative eigenvalue of R grows with increasing ξ, so the boundary condition at ξ → ∞ inforcing regularity requires this solution to be absent in (77). from which the local statements (46) and (47) are now readily deduced.

Appendix B

In this Appendix we derive formulas (57) using ( 52) and (53). To this end, one needs

divB s = -div
t f = - ∂ψ ext ∂ρ f ρ f + ξ f ψ ext n, (87)
τ s = -ρ s ∂ψ ext ∂(∇ρ s ) • n n, τ f = 0,
which agrees with (57).

Fig. 1 .

 1 Fig. 1. Plot of P 2 (p * ) = 0 for a value of α = ±1. Two parabolas are obtained if csc f = 0. The parabolas have their vertices at p * = 1/2α and are open for negative values of p * . If p * > p * c then the static equilibrium of the basic state is unstable.For csc f = 0 the two parabolas reduce to the same pair of straight lines. The one with positive slope has an unbounded stability limit, while that with negative slope has a finite stability limit. These limits are denoted by p * c

2 sf 2ρ 3 0 |β| 2 ,

 22 is a condensed graphical representation in which the transformation reduces (26) toP 2 (p * ) := -αp 2 * ± p * + 1 > 0, α := ρ 2 0 csc f 2 ρ 6 0 γss γ ff -γthus collapsing the family of parabolas to two single graphs.

  ) is a positive function of c sc f , irrespective of whether β > 0 or β < 0; it is increasing (decreasing) for negative (positive) values of c sc f , assumes its maximum at c sc f = 0 and local minima at the boundaries |c sc f | = min(1/ρ 0 s , 1/ρ 0 f ). Denoting by c sc f B the point at which p ext 0 C attains its minimum, it follows from (31) that the minimum critical pressure p ext 0 C (cs-cf ) B depends on, see (27),

a 6 :

 6 = (C 4 -C 1 ) /C 1 . (72) Consider the following change of variable: z = exp (-ξ) ; the differential problem (70) becomesd Y dz = -1 z R + A 0 + z A 1 Y. (73)This change of variable maps the open set (0, ∞) onto the open set (0, 1); boundary conditions in ξ = 0 now are given at z = 1 and conditions at ξ → ∞ at z = 0.

Fig. 2 .

 2 Fig. 2. Scaled solid density rs plotted against dimensionless distance ξ parameterized for various values of the e-folding distance ξ 0 . All curves approach the asymptote as ξ → ∞ With the above expressions for a j (j = 1, 3, 4), we have √ a 1a 3 a 4 = x 0 √ λ s

We are grateful to P. Seppecher for drawing our attention to this fact.

The canonical form of a matrix is defined as its Jordan form. If the matrix admits linearly independent eigenvectors its canonical form is diagonal.

As beautifully summarized by[START_REF] De Boer | The development of the concept of effective stress[END_REF], Fillunger dismissed the fact that the pore pressure would affect the strength of the porous material.

Acknowledgements. The authors wish to thank Prof. Pierre Seppecher from Université de Toulon et du Var for his constructive criticism and long discussions about conservative boundary conditions in mixture theories. They also acknowledge the constructive reviews of two referees.

The parameter ξ 0 can be greater than unity, if det H = H 11 H 22 -H 2 12 < ρ 0 α s H 22 , (81) and it will eventually tend to infinity when H 2 12 → H 11 H 22 . This last condition occurs if the coupling γ sf is sufficiently large: the effect of the second-gradient-depending-deformation energy is an increasing widening of the boundary layer near the external surface of the body when the instability conditions are approached.

Conclusions

In this paper a binary mixture model was presented which possesses the ingredients that an external pressure may cause a dilatation of the pores, a phenomenon sometimes observed in heterogeneous porous materials. Terzaghi and in particular Fillunger were aware of this phenomenon and knew that their models could not predict this behaviour, and they dismissed the possibility after extensive search for evidence and own experimentation 5 . Perhaps they were too ambitious, for their arguments seem to suggest that they were looking for an explosion of the solid matrix of a pressurized solid-fluid mixture if only the external pressure would be sufficiently large.

The present paper showed within the context of a very simple mixture modeltoo simple to describe the deformation of the solid accurately, but sufficiently complex to isolate this detailhow Terzaghi's and Fillunger's search could be interpreted. To this end we assumed all constitutive quantities to depend on the apparent densities of the solid and the fluid (and eventually on their gradients) and no more. It turned out that a dependence of the internal free energy on the interaction term ρ s ρ f , i.e., on the product of the apparent solid and fluid densities, is important. Terzaghi's and Fillunger's explosion is interpreted here as the loss of stability of a pre-stressed reference state. The critical external pressure, which causes this reference state to become unstable, is dictated by two physically distinct properties: i) the coefficients of the quadratic terms of the free internal energy and in particular its interaction term and ii) the parameterization how the external pressure is distributed between the solid and the fluid normal boundary tractions. There are parameter sets for which instability never arises and others for which instability can, in principle, always occur, if only the external pressure is sufficiently large ( §4).

The second question addressed by Terzaghi and Fillunger is the mentioned dilatancy phenomenon viewed possibly by them to be the same as the explosive problem of the solid matrix. Our model also gives an answer to this question. If the pressure corresponding to our reference state is increased then the new stable state can exist only under the above mentioned stability conditions. This new state possesses a smaller or larger apparent solid density provided that the interaction term of the free energy obeys certain equalities involving the other coefficients of the free energy and the parametrization of the constituent boundary tractions ( §5). We find it most intriguing, that thermodynamics of bulk and boundary quantities provides the answer to this subtle behaviour of the mixture.

While we do now understand how the above mentioned dilatational effect can be predicted by the model equations, it does not have typical boundary layer structure in a first gradient theoretical setting. This boundary enhancement can be achieved by adding a density-gradient dependence of the solid phase to the free energy [START_REF] Dell'isola | A variational approach for the deformation of a saturated porous solid. A second gradient theory extending Terzaghi's effective stress principle[END_REF]. This dependence will lead to an enhancement of the pore space close to the interface between the mixture and the exterior world which dies out as one moves away from the interface, ( §6).

The one-dimensional problem which we solve in the last section ( §7) proves that second gradient regularization is necessary if one wants to describe the explosion phenomenon imagined by Fillunger and Terzaghi. Indeed if one interprets it as a loss of stability the transition from stable to unstable states can (in the framework of the present model) be parameterized by the coupling coefficient γ sf , ceteris paribus. When γ sf is increased, H 2 12 also increases and det H tends to zero: the boundary layer at the interface between the mixture body and the external world becomes wider and wider and eventually occupies the whole body before the instability conditions arise. As the second gradient boundary layer is characterized by a lower value of the apparent solid mass density one can state that instability is attained by a progressive dilatational process which is induced by the pore fluid pressure and initially arises at the boundary of the mixture body.

Appendix A

In this section we give the details of the calculation which permit us to obtain relation (45). Using the balance of mass for the solid gives

Consider the following identities for a second order tensor field A, a third order tensor field A 6 and a vector field v

using these identities in (82) for v = v s , A and A being the coefficients of ∇v s and ∇ ⊗ ∇v s in (82) respectively, the aforementioned equation takes the alternative form

Using the divergence theorem where appropriate yields formula (45). Combining this result with the "first gradient " expression of (40)using (42) and ( 43) yields now

(85) 6 We assume that a third order tensor A is a linear map defined as follows

V being a linear space, LIN V the collection of all linear endomorphisms on V and LIN V , LIN V the collection of all the linear morphisms mapping V into LIN V . The transpose of A is assumed to fulfill the following relation

for any u ∈V and any U ∈LIN V .