J Duriez

M Eghbalian

R Wan

F Darve

The micromechanical nature of stresses in triphasic granular media with interfaces

Keywords: granular material, microstructure, stress, unsaturated conditions, interfaces

The total stress tensor as the average stress within a triphasic granular medium is formally derived from micromechanics where internal forces associated with the solid phase, the two immiscible fluid phases and the associated three interfaces are explicitly accounted for. It is demonstrated that for rigid solid particles, the contributions of all local solid-fluid surface tensions to the total stress are eventually zero. The present work gives the total stress expression as a function of a solid-phase specific stress tensor and a fluid mixture stress contribution that is related to the material's microstructure. A generally non-spherical fluid mixture stress is obtained in contrast to an averaged hydrostatic fluid pressure usually associated with standard thermodynamics. The tensorial nature of this fluid mixture stress contribution is highlighted through numerical experiments pertaining to an idealized granular material in the pendular regime at low wetting saturations. Numerical simulations providing full access to microstructural information are conducted using the Discrete Element Method (DEM), which describes internal forces using resultant forces that clearly deviate from the distributed nature of internal forces in triphasic granular media, e.g., fluid pressures. Nevertheless, this micro-scale representation is demonstrated to be indeed valid for macro-scale stress description in the pendular regime.

Introduction

Granular materials have many distinctive features, revealing among others, a dual nature: they can behave as either a fluid or a bulk solid as observed in an hourglass. In the bulk solid state, they possess mechanical strength and bearing capacity so that they are able to sustain foundation loads such as in civil engineering. However, this mechanical strength can drastically change in the presence of a mixture of fluids within the pore space due to the multiphasic condition that sets in, implicating various interactions between the individual solid and fluids phases. In the real world, this phenomenon is clearly evidenced by the change in consistency of sand in a beach as the water content gradually increases towards the shoreline.

As will become more evident later throughout this paper, challenges in the mechanical description of such a multiphasic porous material arise from the variety of stress variables that emerge. On the one hand, the equilibrium under external loads applied on the boundaries of a triphasic granular system involves an overall macroscopic stress, the so-called total stress. Probing into the microscopic scale on the other hand, shows that distinct stresses also exist in the fluids and the solid phase, with for instance distinct individual fluid pressures since the fluid phases are separated by singular surfaces showing interfacial surface tension.

A simplified interpretation of the behaviour through one single stress variable would be of great practical interest, which conceptually led to the Principle of Effective Stress [START_REF] Terzaghi | The shearing resistance of saturated soils and the angle between the planes of shear[END_REF][START_REF] Bishop | Some aspects of effective stress in saturated and partly saturated soils[END_REF][START_REF] Nuth | Effective stress concept in unsaturated soils: Clarification and validation of a unified framework[END_REF]. This concept delivers a very strong statement as it signifies that the mechanical response (strains) of a material is solely controlled by a so-called effective stress related to the solid phase, irrespective of the magnitudes of any of the stresses in the various phases. In other words, when comparing a dry (airsaturated) porous material with when it is invaded with various fluid phases, the Principle of Effective Stress implies that the same mechanical behaviour is to be expected as long as the effective stress loading path under multiphasic conditions coincides with the total stress path in the dry state.

The effective stress concept has been shown to be valid when the pore space is saturated by one single fluid, with the effective stress σ 1 being expressed under such conditions in terms of the total stresses Σ, the fluid pressure u, the solid matrix compressibility K s and the drained bulk compressibility (of the skeleton) K h (Biot, 1941):

σ 1 " Σ ´u ˆ1 ´Kh K s ˙δ « Σ ´u δ (1)
where δ is the identity tensor. Interestingly, the approximated version of Eq. ( 1), as proposed by [START_REF] Terzaghi | The shearing resistance of saturated soils and the angle between the planes of shear[END_REF], is the most venerable equation used in soil mechanics where soil particles are much stiffer than the skeleton as a whole.

By construction, Eq. ( 1) is restricted to the biphasic case, which led Bishop (1959;[START_REF] Bishop | Some aspects of effective stress in saturated and partly saturated soils[END_REF] to address the more commonly encountered triphasic condition, proposing a variant of Terzaghi's effective stress expression, i.e.

σ 1 " Σ ´un p1 ´χq δ ´uw χ δ (2) with χ P r0; 1s the so-called Bishop parameter, and u n , u w the pressures of the non-wetting (n) and wetting (w) immiscible fluid phases, respectively. Water and air in so-called unsaturated soils, or water and oil in petroleum reservoirs

are typical examples of wetting and non-wetting fluids-with the wetting fluid type being dependent on the host soil or rock for oil-water mixtures.

Unfortunately, the application of Eq. (2) to experimental observations on triphasic media faced major difficulties, as reviewed by [START_REF] Nuth | Effective stress concept in unsaturated soils: Clarification and validation of a unified framework[END_REF].

In particular, no stress variable was found by [START_REF] Nuth | Effective stress concept in unsaturated soils: Clarification and validation of a unified framework[END_REF] that would relate to strain in a consistent manner between unsaturated and dry conditions, thus challenging the principle of effective stress. The same conclusion has been reached by [START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF] based on the numerical modelling of dry and unsaturated granular media. On the other hand, considering a solid-phase related stress variable, [START_REF] Khalili | Effective stress in unsaturated soils: Review with new evidence[END_REF]; [START_REF] Lu | Suction stress characteristic curve for unsaturated soil[END_REF]; [START_REF] Nuth | Effective stress concept in unsaturated soils: Clarification and validation of a unified framework[END_REF] obtained a unique stress limit criterion from experimental data at failure in both triphasic and biphasic conditions, which has also been obtained

numerically by [START_REF] Hicher | A microstructural elastoplastic model for unsaturated granular materials[END_REF]; [START_REF] Scholtès | Micromechanics of granular materials with capillary effects[END_REF]; [START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF]. In the end, a consensus is that some aspects of the mechanical behaviour in triphasic conditions, such as failure description, can be inferred from biphasic conditions using a single effective stress variable, whereas a full stress-strain description requires including additionnal stress variables, e.g. the capillary pressure u c " u n ´uw , in constitutive models [START_REF] Nuth | Effective stress concept in unsaturated soils: Clarification and validation of a unified framework[END_REF]Alonso et al., 2010). Basically, this means that the Principle of Effective Stress as stated previously does not hold in triphasic conditions. Nevertheless, it appears that the effective stress denomination is often still associated to the solid-phase related variable mentionned in the above [START_REF] Nuth | Effective stress concept in unsaturated soils: Clarification and validation of a unified framework[END_REF][START_REF] Borja | On the effective stress in unsaturated porous continua with double porosity[END_REF]Alonso et al., 2010). Connected with the possible breakdown of the Principle of Effective Stress is the "two stress variables approach" that does not rely on any effective stress [START_REF] Fredlund | Unsaturated soil mechanics in engineering practice[END_REF]Alonso et al., 1990), as well as the suction hardening concept according to which the capillary pressure enters the hardening parameters [START_REF] Borja | On the effective stress in unsaturated porous continua with double porosity[END_REF], leading to discrepancies between biphasic and triphasic conditions.

The search for a single effective stress variable in the realm of geomechanics may be formally compared to the search for a back stress variable, as coined by [START_REF] Nowick | Dislocation theory as applied by N.A.C.A. to the creep of metals[END_REF], in Solid Mechanics. Tensorial back stress variables enter elasto-viscoplastic models as the difference between the total stress tensor and a dissipative stress variable which governs the material behaviour [START_REF] Prager | A new method of analysing stresses and strains in workhardening plastic solids[END_REF][START_REF] Lemaitre | Mechanics of solid materials[END_REF][START_REF] Frederick | A mathematical representation of the multiaxial Bauschinger effect[END_REF][START_REF] Henann | A large deformation theory for rate-dependent elastic-plastic materials with combined isotropic and kinematic hardening[END_REF]. Evolving yield surface centers in kinematic hardening models are typical back stress examples.

Modelling triphasic porous media according to the Principle of Effective

Stress would also require subtracting from the total stress a stress contribution accounting for internal forces due to the fluid mixture, leaving the solid-phase related effective stress to govern material deformation and failure. The fluid mixture stress contribution is denoted as suction stress in the general case [START_REF] Lu | Suction stress characteristic curve for unsaturated soil[END_REF][START_REF] Nikooee | Effective stress in unsaturated soils: A thermodynamic approach based on the interfacial energy and hydromechanical coupling[END_REF] and reduces to the so-called capillary stresses [START_REF] Hicher | A microstructural elastoplastic model for unsaturated granular materials[END_REF][START_REF] Scholtès | Micromechanics of granular materials with capillary effects[END_REF][START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF] for granular media with non-viscous fluids. Indeed, the pore fluid mixture is exclu-sively governed by capillarity in such media, with other physico-chemical forces and adsorbed liquid layers being negligible. The present manuscript focuses on capillary stresses only, hoping that a correct capillary stress determination as proposed herein will subsequently bring some insights into a more general effective stress discussion. At the moment, the following results have been obtained

for the nature of capillary stresses.

In line with Bishop's (1959) expression given in Eq. ( 2), capillary stresses are often expressed in terms of the fluid pressures only (e.g., [START_REF] Zienkiewicz | Computational geomechanics with special reference to earthquake engineering[END_REF][START_REF] Nuth | Effective stress concept in unsaturated soils: Clarification and validation of a unified framework[END_REF][START_REF] Borja | On the effective stress in unsaturated porous continua with double porosity[END_REF][START_REF] Lu | A closed-form equation for effective stress in unsaturated soil[END_REF]. However, triphasic conditions imply singular discontinuity surfaces as fluid-fluid interfaces which are absent in the biphasic case. This motivated other approaches to include interfacial surface tension forces in the stress analysis. Among those, [START_REF] Lunati | Young's law and the effects of interfacial energy on the pressure at the solid-fluid interface[END_REF][START_REF] Shikhmurzaev | On Young's (1805) equation and Finn's (2006) 'counterexample[END_REF]. Nevertheless, solid-fluid surface tensions, being described as early as in Young's equation ( 1805), are evidenced by the contraction of very small metal particles (e.g. [START_REF] Vermaak | On surface stress and surface tension: I. theoretical considerations[END_REF][START_REF] Jiang | Lattice contraction and surface stress of fcc nanocrystals[END_REF] or may let sand adopt a hydrophobic behaviour [START_REF] Vitz | Magic sand: Modeling the hydrophobic effect and reversed-phase liquid chromatography[END_REF].

A first objective of the present work is to clarify under which conditions the solid-fluid interfaces contribute or not to the stress of triphasic granular media. In fact, it will be rigorously demonstrated that the stress contribution of solid-fluid surface tension is zero on average for negligible solid-solid surfaces, an important result absent from previous works bearing similarities with Section 2 such as [START_REF] Chateau | Micromechanics of unsaturated granular media[END_REF]; Dormieux et al. (2006a); [START_REF] Gray | The solid phase stress tensor in porous media mechanics and the Hill-Mandel condition[END_REF].

Also, most approaches represent the capillary stresses as a scalar or a spherical, i.e. isotropic, tensor [START_REF] Hassanizadeh | Mechanics and thermodynamics of multiphase flow in porous media including interphase boundaries[END_REF][START_REF] Nuth | Effective stress concept in unsaturated soils: Clarification and validation of a unified framework[END_REF][START_REF] Borja | On the effective stress in unsaturated porous continua with double porosity[END_REF][START_REF] Gray | The solid phase stress tensor in porous media mechanics and the Hill-Mandel condition[END_REF][START_REF] Lu | A closed-form equation for effective stress in unsaturated soil[END_REF][START_REF] Nikooee | Effective stress in unsaturated soils: A thermodynamic approach based on the interfacial energy and hydromechanical coupling[END_REF].

This feature arises in these works from the usual consideration of mixture theory and thermodynamics framework, often adopting a restrictive scalar or spherical description of stresses and strains. Indeed, other micromechanical approaches suggest that capillary stresses show a generic tensorial nature [START_REF] Chateau | Micromechanics of unsaturated granular media[END_REF][START_REF] Hicher | A microstructural elastoplastic model for unsaturated granular materials[END_REF][START_REF] Scholtès | Micromechanics of granular materials with capillary effects[END_REF][START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF][START_REF] Duriez | Stress in wet granular media with interfaces via homogenization and discrete element approaches[END_REF][START_REF] Li | Meso-hydromechanically informed effective stresses and effective pressures for saturated and unsaturated porous media[END_REF]. It is tempting to consider micromechanics to bring more pertinent insights to the mechanical behaviour, such as done by [START_REF] Sibille | Granular plasticity, a contribution from discrete mechanics[END_REF], because granular materials show a clear microstructure with also a fluid 'fabric' in triphasic conditions [START_REF] Manahiloh | Characterizing liquid phase fabric of unsaturated specimens from X-Ray computed tomography images[END_REF].

A second objective of the manuscript is thus to clearly establish and illustrate how the microstructure enters into the stress state of triphasic granular media, leading to deviatoric capillary stresses as shown in Section 4. It is noteworthy that [START_REF] Madeo | A continuum model for deformable, second gradient porous media partially saturated with compressible fluids[END_REF] evidenced the need of a second gradient theory to describe the mechanical behaviour of unsaturated porous media. However, the scope of the present paper is restricted to the stress tensor and the first displacement gradient being work-conjugate.

The present paper begins with the micromechanical derivation of stresses in a triphasic granular medium in Section 2, isolating the sought after capillary stresses. A clear microstructure dependency is proposed for the capillary stresses, and because experimental microstructure measurements in triphasic granular media, e.g from X-ray computed tomography [START_REF] Culligan | Interfacial area measurements for unsaturated flow through a porous medium[END_REF][START_REF] Bruchon | X-ray microtomography characterisation of the changes in statistical homogeneity of an unsaturated sand during imbibition[END_REF], are not yet precise enough to be used in conjunction with the analytical developments presented herein, Section 3 recalls a numerical model of ideal triphasic granular media at low wetting saturation that provides full access to pertinent microstructural information, including interfaces (Duriez and Wan, 2016a). Based on the Discrete Element Method [START_REF] Cundall | A discrete numerical model for granular assemblies[END_REF], the model describes by design the micro-scale stress field from resultant vector forces, that are distinct in nature from distributed fluid pressure or surface tension internal forces. [START_REF] Chateau | Micromechanics of unsaturated granular media[END_REF]; [START_REF] Duriez | Stress in wet granular media with interfaces via homogenization and discrete element approaches[END_REF] nevertheless suggested that such DEM models are still able to describe stresses in triphasic granular media, but in the absence of any formal analytical demonstration.

Hence, a third objective of the paper is to formally provide such a demonstration as touched upon in Section 3.2. Finally, Section 4 explores in details important characteristics and properties of the analytically derived capillary stresses with the aid of DEM.

The notation adopted in this paper is such that a superimposed arrow designates first order tensors (vectors) such as x, ν, and bold symbols are employed for second order tensors like δ, Σ. Also, the geomechanics sign convention is used, i.e. σ n corresponds to an applied loading on a given body when the normal vector n points inwards to the boundary, leading to a positive sign for compressive stresses σ.

Stress homogenization of triphasic granular media with interfaces

The representative elementary volume (REV) V of a triphasic granular material features, at the micro-scale, different volume phases V α , α P tn; s; wu for the solid (s) and fluids (n, w) phases, as well as interface surfaces S αβ , α and β ‰ α P tn; s; wu, and contact lines Γ where all three phases meet. Extending [START_REF] Duriez | Stress in wet granular media with interfaces via homogenization and discrete element approaches[END_REF], surface tension is here described along all three interfaces, and we first underscore the corresponding consequences that ensue for the stress description at the micro-scale.

Micro-scale stress description along interfaces and contact lines

An interface S αβ considered here as a distinct medium is subject to internal surface tension forces proportional to the surfacic energy γ αβ . The following stress-like tensor π αβ describes these internal surface tension forces that are tensile in nature [START_REF] Chateau | Micromechanics of unsaturated granular media[END_REF]Dormieux et al., 2006b;[START_REF] Gray | The solid phase stress tensor in porous media mechanics and the Hill-Mandel condition[END_REF][START_REF] Duriez | Stress in wet granular media with interfaces via homogenization and discrete element approaches[END_REF]:

π αβ " γ αβ p n b n ´δq (3)
Following the geomechanics sign convention, Eq. (3) basically expresses the surface stress-like tensor π αβ in terms of γ αβ and the projection tensor onto S αβ , δ αβ " pδ ´ n b nq, with δ the identity tensor, and n the normal to S αβ .

n The interface material points experience along the contour C the surface tension forces described by π αβ , as well as internal forces within α (described by Cauchy stress σ α ), and internal forces within β (described by σ β ). Hence, the equilibrium of forces at the interface material points leads to:

ż S αβ pσ α ´σβ q n dS ´żC γ αβ ν αβ dl " 0 (4)
with n pointing from α to β and ν αβ the inward conormal tangent to S αβ and normal to C, and dl ą 0 the arc length (Fig. 1). Next, the contour integral in Eq. ( 4) is classically rewritten using the following equality (see e.g., [START_REF] Gurtin | A continuum theory of elastic material surfaces[END_REF][START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF]:

ż C ν αβ dl " ż S αβ div S p nq n dS (5)
with div S pq the surface divergence operator as defined by [START_REF] Gurtin | A continuum theory of elastic material surfaces[END_REF]: roughly speaking, div S pq refers to the projection of the classical differential operator ∇. " divpq onto the tangent plane of the surface S of interest.

Eq. ( 5) finally gives the normal stress discontinuity across S αβ due to capillarity; in case S αβ is curved:

pσ α ´σβ q n " γ αβ div S p nq n (6)
Considering the fluid-fluid interface S nw , σ α " u α δ, α P tn; wu with u α the fluid pressures, together with Eq. ( 6) give the Laplace-Young equation:

u c " u n ´uw " γ nw div S p nq (7)
with n towards the wetting fluid. The capillary pressure u c can be identified as the matric suction in the present case of granular media.

Considering the solid-fluid interface S sβ , Eq. ( 6) serves as a boundary condition for the solid stresses that accounts for solid-fluid surface tension [START_REF] Chateau | Micromechanics of unsaturated granular media[END_REF]Dormieux et al., 2006a;[START_REF] Lunati | Young's law and the effects of interfacial energy on the pressure at the solid-fluid interface[END_REF][START_REF] Gray | The solid phase stress tensor in porous media mechanics and the Hill-Mandel condition[END_REF], with n being the outward solid normal:

σ s n " pu β `γsβ div S nq n (8)
A different boundary condition for the solid stresses holds along the contact line Γ where the solid phase meets the two fluid phases (Fig. 2). Here, equilib- 9), dl is a portion of Γ, and dS is an elemental surface area constructed from dl and an infinitesimally small length (magnified in the Figure) along t rium consideration for a peripheral solid facet dS located along Γ and sustaining solid internal forces and three surface tension forces leads to:

σ s n dS`p γ nw ` γ ns ` γ sw q dl " 0 ô σ s n dS " `´ γ nw `pγ sw ´γns q t ˘dl (9)
where n is the outward solid normal and t a tangent vector. In connection with

Young's equation ( 1805) given as

γ nw cos θ " γ ns ´γsw (10)
with θ being the contact or wetting angle (Fig. 2), a local zero tangential solid stress along Γ can be deduced from Eq. (9):

pσ s nq. t " γ nw cos θ `γsw ´γns " 0 (11)
Now that the micro-scale stress description has been recalled, the macroscale stress tensor is next derived.

Macro-scale stresses of triphasic granular media

The total stress Σ for the REV is volume-averaged using the microscopic stress tensors within all three volume phases and three interfaces [START_REF] Chateau | Micromechanics of unsaturated granular media[END_REF]Dormieux et al., 2006a;[START_REF] Gray | The solid phase stress tensor in porous media mechanics and the Hill-Mandel condition[END_REF]. Classically, no stress tensor or other mechanical properties is considered for the contact line itself [START_REF] Hassanizadeh | Mechanics and thermodynamics of multiphase flow in porous media including interphase boundaries[END_REF][START_REF] Gray | The solid phase stress tensor in porous media mechanics and the Hill-Mandel condition[END_REF][START_REF] Nikooee | Effective stress in unsaturated soils: A thermodynamic approach based on the interfacial energy and hydromechanical coupling[END_REF]:

Σ " 1 V « ÿ α"n,s,w ż Vα σ α dV `ÿ α,β"n,s,w ż S αβ π αβ dS ff (12)
The fluid stresses are the isotropic pressures u n δ, u w δ. Considering thermodynamic equilibrium, the fluid pressures are assumed to be uniform, leading to a uniform capillary pressure, whatever the saturation regime and the connected or disconnected feature of the fluid phases.

Regarding the solid phase, the latter is made of distinct particles p, V s " Ť p V p , and the volume integral for each particle is classically rewritten as a surface integral from the divergence theorem, considering equilibrium of particles without body forces:

ż Vp σ s dV " ż Sp pσ s nq b x dS (13)
with n the outward normal and the positions x being defined from an adequate centroid for each particle thanks to particles equilibrium.

Paying attention to the microstructure of triphasic granular media, the particle surfaces S p show solid-fluid interfaces S p,sα , α P tn; wu, contact lines Γ p , and contact surfaces S p,ss . Eq. ( 8), resp. ( 9), applies along S p,sα , resp. Γ p .

Considering quasi-rigid particles, the contact surfaces S p,ss are negligible and thus, ż

Sp,ss

pσ s nq b x dS " ´ÿ cp f c b x c (14)
with f c being the contact forces, as sustained by p, due do the different contacts

tc p u.
Accounting for the above, algebraic manipulations of Eq. ( 12) lead to the following expression for the total stresses:

Σ ´un δ " 1 V " ř c f c b l ´uc ´Vw δ `şSsw n b x dS γnw ´şSnw pδ ´ n b nq dS `şΓ ν nw b x dl ř α"n;w γ sα ´şSsα pdiv S n n b x ` n b n ´δq dS ´şΓ ν sα b x dl ¯ı (15)
The above stress expression together with the physical interpretation of each term will be amply discussed in the following sections.

Solid-fluid surface tension zero stress contribution

We herein focus on the last stress contributional term on the r.h.s of Eq. ( 15), arising from solid-fluid surface tensions γ sα . First, we observe that such stress contributions are due directly to the internal forces within the interfaces S sα by integrating surface tension force field π sα over the membrane surface. Other γ sα -stress contributions come from the stress application by S sα onto the solid phase, through the two boundary conditions for the solid stresses, i.e. Eq. ( 8)

along the surfaces S sα excluding the contact lines, and Eq. ( 9) along the contact lines. However, these γ sα -stress contributions cancel out for each interface since the following equality holds for any surface S of contour C " Γ, conormal ν, and normal n (Fig. 1):

ż S pdiv S p nq n b x ` n b n ´δq dS " ż Γ ν b x dl (16)
Eq. ( 16) is formally demonstrated using Stokes theorem in Appendix A. In fact, this important result is a generalization of the following equation proposed by [START_REF] Rosenkilde | Surface-energy tensors[END_REF] that is restricted to closed surfaces S c :

ż Sc div S p nq n b x dS " ż Sc pδ ´ n b nq dS (17)
Obviously, Eq. ( 16) reduces to (17) for closed surfaces showing no contour Γ.

The conclusion of a zero average stress contribution from solid-fluid surface tensions is valid for any porous media including one or two immiscible fluids and showing negligible solid-solid surfaces S ss such that the solid-fluid surfaces S sα are always bound by contact lines Γ. In this sense, the present result obtained from quasi-static micromechanics is consistent with the thermodynamical approach followed by [START_REF] Nikooee | Effective stress in unsaturated soils: A thermodynamic approach based on the interfacial energy and hydromechanical coupling[END_REF]. Indeed, Nikooee et al.

(2013) also showed, from alternative free energy considerations, that solid-fluid surface tensions would affect the behaviour only for deformable solid particles.

It should be noted that [START_REF] Chateau | Micromechanics of unsaturated granular media[END_REF] previously mentioned this zero stress contribution from solid-fluid surface tensions without a clear justification.

Also, it is worth mentioning that the herein proposed Eq. ( 16) suggests possible simplications in the final result of [START_REF] Gray | The solid phase stress tensor in porous media mechanics and the Hill-Mandel condition[END_REF].

Capillary stresses and microstructure

As a consequence of the result presented in the previous section, Eq. ( 15)

together with ( 16) lead to the total stress expression for the triphasic REV as:

Σ ´un δ "σ cont ´1 V " u c ˆVw δ `żSsw n b x dS γnw ˆżSnw pδ ´ n b nq dS `żΓ ν nw b x dl ˙ (18) σ cont " 1 V ÿ c f c b l (19)
The first stress contribution on the r.h.s of Eq. ( 18) is denoted herein as the contact stress tensor σ cont since it accounts for the stress interaction at any contact c between a solid particle pair 1-2 through the contact force f c as sustained by 2, l being the branch vector from the centre of 1 to the one of 2. This first term describes the part of the solid phase internal forces that is due to the solid phase itself, excluding stress applications by the other phases.

In this sense, it has thus been denoted as effective stress by Lu and Likos Finally, the third stress contribution in Eq. ( 18) arises from the γ nw surface tension forces related to the fluid-fluid interface S nw . In particular, the internal forces within the S nw interface are described in terms of γ nw and the local interface orientation pδ ´ n b nq. The stress application by S nw onto the solid particles along the contact lines, depending on the wettability ν nw , is also accounted for.

These two former stress contributions coming from the fluid phases and interfaces correspond to the capillary stresses σ cap mentioned in the introduction and defined as:

σ cap " Σ ´un δ ´σcont (20)
We show herein the capillary stresses depend not only on the fluid pressures u n , u w , such as proposed by Bishop (1959); [START_REF] Zienkiewicz | Computational geomechanics with special reference to earthquake engineering[END_REF]; [START_REF] Nuth | Effective stress concept in unsaturated soils: Clarification and validation of a unified framework[END_REF]; [START_REF] Borja | On the effective stress in unsaturated porous continua with double porosity[END_REF]; [START_REF] Lu | A closed-form equation for effective stress in unsaturated soil[END_REF], but also on the surface tension γ nw and more importantly fluid-fluid and solid-fluid interface topologies. Furthermore, we propose a clear microstructure dependency through the following microstructure tensors defined as:

µ V w " V w δ (21) µ Ssw " ż Ssw n b x dS (22) µ Snw " ż Snw pδ ´ n b nq dS (23) µ Γ " ż Γ ν nw b x dl (24)
such that the capillary stresses express as:

σ cap " ´1 V " u c pµ V w `µSsw q `γnw pµ Snw `µΓ q ı (25)
It is important to recall that a general tensorial (non-spherical) form of the capillary stress is possible due to the nature of all microstructure tensors other than µ V w in Eq. ( 25). This general non-spherical form for the capillary stress could not be obtained from a thermodynamics approach as in Hassanizadeh and [START_REF] Hassanizadeh | Mechanics and thermodynamics of multiphase flow in porous media including interphase boundaries[END_REF]; [START_REF] Lu | A closed-form equation for effective stress in unsaturated soil[END_REF]; [START_REF] Nikooee | Effective stress in unsaturated soils: A thermodynamic approach based on the interfacial energy and hydromechanical coupling[END_REF]. The shortcoming arises from the constitutive assumptions considered in those works where the free energies are systematically expressed according to scalars such as saturation, whereas tensorial dependencies would be necessary to account for a possibly oriented microstructure.

Being obtained from static considerations without any kinematics, the stress decomposition in Eq. ( 18) however does not attach any constitutive meaning to the contact and capillary stress tensors. Further investigations are thus led in the remaining developments, where the focus is on the capillary stresses for two extreme saturation regimes.

Capillary stresses for high wetting saturation

In biphasic conditions, the capillary stresses are spherical, with the stress decomposition, Eq. ( 18), reducing to a form equivalent to the effective stress expression for such conditions, Eq. (1):

Σ " σ cont `u δ (26)
with the contact stress tensor σ cont that operates in the solid phase serving as the effective stress for granular materials with negligible solid matrix compressibility.

The same stress decomposition as Eq. ( 26), together with a spherical nature for the capillary stress, also holds under triphasic conditions in the so-called capillary regime showing high wetting saturations. In this regime, the nonwetting fluid forms isolated bubbles which are of spherical shape to minimize fluid-fluid specific energy. This leads in particular to the following expression for the microstructure tensor µ Snw , obtained from the algebraic equality (17) applied to the closed interfaces S nw together with Laplace-Young Eq. ( 7):

µ Snw " u c γ nw V n δ (27)
such that it is easy to get:

σ cap " ´uc δ (28)
and Eq. ( 26) with u " u w is recovered, as suggested by [START_REF] Chateau | Micromechanics of unsaturated granular media[END_REF] using a different justification.

The sphericity of capillary stresses in such a saturation regime corresponds to the absence of any preferred orientation at the micro-scale when solid surfaces are completely wetted and nw interfaces are spherical. On the other hand, a preferred orientation is expected for low saturation, in the so-called pendular regime involving oriented capillary bridges.

3. Multiscale simulation of ideal granular media in the pendular regime

Model formulation

The DEM model considered herein (Duriez and Wan, 2016a) is built within the code Yade ( Šmilauer et al., 2015), being inspired by the one by [START_REF] Scholtès | Micromechanics of granular materials with capillary effects[END_REF]. It applies to low wetting saturations-the so-called pendular regime for which the wetting fluid forms isolated capillary bridges, i.e. menisci between particle pairs (Fig. 3). This saturation regime is of limited extent, up to 5-10% wetting saturation ratio roughly, but nevertheless provides granular media with very distinct features (e.g. apparent cohesion) compared with biphasic conditions, making this regime worthy of investigation [START_REF] Chateau | Micromechanics of unsaturated granular media[END_REF][START_REF] Scholtès | Micromechanics of granular materials with capillary effects[END_REF][START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF][START_REF] Wang | Anisotropy of a tensorial Bishop's coefficient for wetted granular materials[END_REF].

R 2 ≥ R 1 R 1 θ θ δ 2 δ 1 d s s w n Figure 3: A capillary bridge
The DEM model considers a collection of spherical discrete elements (DE), possibly of different radii R p , that correspond in shape and size to the solid particles of an idealized granular medium. For such media, the capillary stresses follow directly from Eq. (25) as:

σ cap " ´1 V « u c ˜Vw δ `ř p R p ş Sp,sw n b n dS γnw ˜şSnw pδ ´ n b nq dS `ř p R p ş Γp ν nw b n dl ¸ff (29)
Since actual solid particles and DE are quasi-rigid, strains can only arise from particle relative displacements under pair-wise mutual force interactions between the particles.

As a starting point, contacting DE sustain point forces approximating contact interaction following simplified contact laws that express the normal (resp. tangential) contact force as a function of the normal (resp. tangential) relative displacement. In the normal direction, the contact law is linearly elastic in compression and disregards tension. For the tangential direction, the behaviour is governed by Coulombian friction leading to a linear elastic-plastic contact law. Such contact description involves three model parameters, relevant to the contact scale (see Wan et al., 2015, for details).

The next step is to incorporate capillary forces to account for the pore fluid mixture in terms of menisci topology between pairs of particles. A rigorous capillary force determination starts with solving numerically the Laplace-Young

Eq. ( 7) in order to compute the axisymmetric geometry of capillary bridges, gravity being neglected (Fig. 3). The existence and properties of the menisci depend on the particle radii R 1 , R 2 , their separating distance d, the fluid-fluid surface tension, and assumed uniform contact angle and capillary pressure values. The model predicts by design a fluid distribution obeying uniform capillary pressure conditions, as done by [START_REF] Scholtès | Micromechanics of granular materials with capillary effects[END_REF]; [START_REF] Wang | Anisotropy of a tensorial Bishop's coefficient for wetted granular materials[END_REF].

These uniform suction conditions are specific to thermodynamic equilibrium and the model does not include any pore flow consideration that could address transient phenomena, contrary to [START_REF] Gili | Microstructural deformation mechanisms of unsaturated granular soils[END_REF]; [START_REF] Mani | Role of contact-angle hysteresis for fluid transport in wet granular matter[END_REF]; [START_REF] Sivanesapillai | A CSF-SPH method for simulating drainage and imbibition at pore-scale resolution while tracking interfacial areas[END_REF]. Also, simulations are currently conducted using a constant contact angle irrespective of advancing or receding movements of the contact line along the solid surface, thereby neglecting contact angle hysteresis phenomena. However, the model has previously been compared in terms of water retention behavior with experiments [START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF][START_REF] Duriez | Stress in wet granular media with interfaces via homogenization and discrete element approaches[END_REF].

Once the capillary bridge distribution is determined, attractive capillary forces f cap between capillary bonded particle pairs are finally obtained from the capillary pressure and fluid-fluid surface tension loadings on either spherical particle:

|| f cap || " π R 1 sin δ 1 pu c R 1 sin δ 1 `2 γ nw sinpθ `δ1 qq " π R 2 sin δ 2 pu c R 2 sin δ 2 `2 γ nw sinpθ `δ2 qq (30)
with δ 1 , δ 2 the filling angles on each particle 1, 2 (Fig. 3). Considering either particle to compute the capillary force norm-using index 1 or 2 in Eq. ( 30)-is immaterial due to the meniscus equilibrium.

In line with the previous focus on solid-fluid surface tensions, it is worth to clarify that these do not contribute to the capillary force. A first possible justification to this assertion is to consider the solid-fluid surface tensions as an internal force to the solid particles, thus leading to no resultant external forces on the DE. Alternatively, in case the DE were assimilated to solid particles separately from the solid-fluid interfaces, two external loadings proportional to γ sα , α " n; w would apply on the DE:

1. a first one along S sα : ş Ssα ´γsα div S p nq n dS, see Eq. ( 8) 2. a second one along the contour C of S sα : ş C γ sα ν sα dl, see Eq. ( 9)

and the resulting loading would be zero, see Eq. ( 5) and also [START_REF] Lunati | Young's law and the effects of interfacial energy on the pressure at the solid-fluid interface[END_REF]; [START_REF] Shikhmurzaev | On Young's (1805) equation and Finn's (2006) 'counterexample[END_REF]. The absence of solid-fluid surface tension γ sα from the resultant forces used in DEM to describe internal forces is consistent-though not equivalent-with their zero contribution to the total stresses of a triphasic material, discussed in Section 2.2.

Through the capillary force determination, the particle size distribution, the fluid-fluid surface tension and the contact angle constitute additional model parameters (Duriez and Wan, 2016a).

Numerical description of stress

According to their calculation paradigm, DEM models of multiphasic conditions [START_REF] Gili | Microstructural deformation mechanisms of unsaturated granular soils[END_REF][START_REF] Scholtès | Micromechanics of granular materials with capillary effects[END_REF][START_REF] Richefeu | Stress transmission in wet granular materials[END_REF] including the present one describe all internal forces of a triphasic granular medium exclusively with interaction forces between DE. In particular, non-solid phases with their internal forces are not directly included, but rather reduced to the resultant forces induced on the solid phase. A total stress tensor is then defined for the DEM sample from contact and capillary forces for all contacts c and menisci m, applying the classical homogenization techniques to the DE (Duriez and Wan, 2016c):

Σ " 1 V ÿ DE ż VDE σ dV " 1 V ÿ c f c b l `1 V ÿ m f cap b l " σ cont `σcap DEM (31) σ cap DEM " 1 V ÿ m f cap b l (32)
By contrast with the computation of particles relative displacements and strains in granular media, the sole consideration of resultant forces may be inadequate to describe stresses in a general manner, as evidenced by any solid sustaining non-zero stresses while in equilibrium. For instance, the DEM is unable to capture the stress state of a biphasic granular media showing a uniform fluid pressure because the fluid-induced resultant force on solid particles is zero in such case. For this reason, the DEM total stress in Eq. ( 31) is to be compared with the net stress pΣ´u n δq in Eq. ( 18) derived from analytical homogenization.

The current stress analysis of the DEM model involves the same contact stress tensor σ cont as the one derived previously from analytical homogenization, Eq. ( 19). As for the fluid mixture stress contribution, the DEM framework leads to a capillary stress expression from resultant capillary forces, clearly distinct in nature from the previous expression obtained from analytical homogenization, Eq. ( 29). As such, we refer to the capillary stress tensor defined by Eq. ( 32) as σ cap DEM , as opposed to the previous capillary stress tensor σ cap expressed by Eq. ( 29). Nevertheless, [START_REF] Chateau | Micromechanics of unsaturated granular media[END_REF] suggested with little justification that a stress description from resultant forces is still possible, such that σ cap " σ cap DEM . Duriez and Wan (2016b,c) presented numerical comparisons confirming it, provided the fluid-fluid interface is accounted for in addition to bulk phases [START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF]. In order to further assess the equivalence of the stress state of DEM models with the stress state of simulated triphasic materials, as derived from homogenization, the equality of the first invariants of σ cap and σ cap DEM is herein demonstrated from the Laplace-Young equation as follows.

The mean capillary stress as predicted by the homogenization approach, denoted as p cap hom , is obtained taking the trace of Eq. ( 29):

p cap hom " ´1 3 V « u c ˜3 V w `ÿ p R p S p,sw ¸`γ nw ˜2 S nw `ÿ p R p Γ p sin θ ¸ff (33)
By construction of the DEM model, the comparison is restricted to the pendular regime for which the wetting fluid forms menisci between particle pairs pi, jq " p1, 2q. Then, Eq. ( 33) is rewritten accounting for menisci microstructure:

p cap hom " ´1 3 V ÿ m « u c ˜3 V m `ÿ i"1,2 R i S i,sw ¸`γ nw ˜2 S m `ÿ i"1,2 R i Γ i sin θ ¸ff (34)
For a given meniscus of volume V m and fluid-fluid surface S m , the wetted surface S i,sw on each spherical particle i is S i,sw " 2 π R i 2 p1 ´cos δ i q and the meniscus contour Γ i is equal to 2 π R i sin δ i , see Fig. 3.

On the other hand, the mean stress p cap DEM of σ cap DEM involves the scalar products f cap . l through the trace of f cap b l, see Eq. ( 32), such that:

p cap DEM " ´1 3 V ÿ m `uc πR 1 2 sin 2 δ 1 `γnw Γ 1 sinpθ `δ1 q ˘pR 1 `R2 `dq (35)
Comparing the two expressions of the mean capillary stress, qualitative differences clearly appear with the menisci surface S m and volume V m entering Eq. ( 34) and not (35) for instance. However, using the Laplace-Young equation and applying Eq. ( 16) to the bridge surface S m , Appendix B demonstrates that for any meniscus:

u c ˜3 V m `ř i"1,2 R i S i,sw ¸`γ nw ˜2 S m `ř i"1,2 R i Γ i sin θ " `uc πR 1 2 sin 2 δ 1 `γnw Γ 1 sinpθ `δ1 q ˘pR 1 `R2 `dq (36)
Eq. ( 36) is a general property of capillary bridges between spherical particles, under mechanical equilibrium conditions i.e. complying with the Laplace-Young equation. It could be used as such to validate triphasic microstructure experimental measurements at the capillary bridge scale. Used here at the REV scale, it remarkably proves that the mean capillary stresses provided by DEM modelling or analytical homogenization are finally equivalent, validating the stress description by DEM in spite of the above mentionned possible reasons for discrepancy.

From a numerical point of view, the complete equality between σ cap DEM and σ cap , considering other invariants than p cap as well as principal directions, has been recently obtained by [START_REF] Duriez | Stress in wet granular media with interfaces via homogenization and discrete element approaches[END_REF].

Microstructure and anisotropy of capillary stresses in the pendular regime

The DEM model is finally applied to a slightly polysized numerical packing of spherical particles. Table 1 presents the retained parameters. [START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF]. Their values are taken from comparable DEM simulations on glass beads [START_REF] Hanes | Simulations and physical measurements of glass spheres flowing down a bumpy incline[END_REF][START_REF] Favier | Predicting the drag coefficient of a granular flow using the discrete element method[END_REF]. Regarding the capillary interaction, the particle size distribution is chosen to be uniform in number between two extreme diameters D min , D max .

The void ratio e of the considered packing is such that the sample shows a dense, i.e. dilatant behaviour under shearing for the considered confining pressures.

Finally, an air-water fluid mixture is considered with the relevant surface tension value at ambient temperature and a small, though not zero, contact angle value θ, in accordance with the wetting of glass beads by water [START_REF] Richefeu | Stress transmission in wet granular materials[END_REF][START_REF] Scheel | Morphological clues to wet granular pile stability[END_REF].

Isotropic solid packing

An isotropic packing is first considered with the numerical sample being µ " m δ (37)

with m defined as the first invariant of µ:

m " trpµq 3 (38)
and being expressed for the various µ tensors as:

m V w "V w (39) m Ssw " 1 3 ÿ p R p S p,sw (40)
m Snw " 2 3 S nw (41) m Γ " 1 3 ÿ p R p Γ p sin θ (42)
Thus, the capillary stresses σ cap reduce to their spherical part p cap δ, with p cap previously expressed according to the microstructure in Eq. ( 33):

σ cap " ´1 V " u c pm V w `mSsw q `γnw pm Snw `mΓ q ı δ " p cap δ (43)
The capillary stress together with its distinct contributions are depicted according to the wetting saturation in Fig. 5, based on microstructural information obtained from the DEM model as given in Fig. 4. The capillary stress shows a non-zero limit when the saturation asymptotically tends to zero, while it is strictly zero for the dry case. This non-zero limit is related to the non-zero limit of the capillary force for a capillary bridge (between contacting particles) whose volume asymptotically vanishes [START_REF] Gili | Microstructural deformation mechanisms of unsaturated granular soils[END_REF][START_REF] Herminghaus | Dynamics of wet granular matter[END_REF]. For such infinitesimal saturations, the capillary stress mostly depends on the capillary pressure and the minutely wetted solid surfaces, the other stress contributions being negligible. For higher saturations close to the upper limit of the pendular regime, a significant contribution from surface tension arises through fluid-fluid interfaces and contact lines, see Fig. 5(a). Generally speaking, and as it is obvious for the contact line term, the relative weight of the distinct contributions depends on the wetting angle Even though an increasing saturation magnifies the relevant micro-features, e.g.

the wetted surfaces (Fig. 4), it is accompanied by a decrease in capillary pressure which directly enters the capillary stress expression, Eq. ( 43).

Also, the mean capillary stress as obtained here from micromechanics is compared with the classical expression p cap " ´uc S r . The latter expression corresponds to the early proposal by [START_REF] Bishop | Some aspects of effective stress in saturated and partly saturated soils[END_REF], taking the corresponding χ parameter equal to saturation S r . This expression has more recently also been proposed from thermodynamical considerations [START_REF] Hassanizadeh | Mechanics and thermodynamics of multiphase flow in porous media including interphase boundaries[END_REF][START_REF] Lu | A closed-form equation for effective stress in unsaturated soil[END_REF]. However, such an expression significantly underestimates the capillary stresses in the pendular regime, as compared with the present micromechanical approach, see Fig. 5(b). The discrepancy is believed to arise because of the missing microstructure dependencies in the constitutive assumptions formulated in such thermodynamic approaches. It should be noted that [START_REF] Nikooee | Effective stress in unsaturated soils: A thermodynamic approach based on the interfacial energy and hydromechanical coupling[END_REF] later on introduced new dependencies in the constitutive assumptions made by [START_REF] Hassanizadeh | Mechanics and thermodynamics of multiphase flow in porous media including interphase boundaries[END_REF], modifying the initial result.

Anisotropic solid packing

A second configuration is herein considered such as the classical critical state which is reached after sustained shearing of a granular material. Here, the numerical sample is brought to critical state by imposing an axisymmetric, socalled "triaxial", compression up to 45% axial strain ε 1 " ε yy , under 10 kPa constant lateral pressure and dry conditions.

At this stage, an anisotropic contact network has been induced by the deviatoric loading, starting from an initially isotropic configuration. Imposing subsequently triphasic conditions to the sample with various capillary pressure values, and solving for capillary bridges between all contacting and separated particles, an anisotropic liquid distribution arises since capillary bridges preferably form between close particles. The microstructure tensors then may show both a spherical and deviatoric component, µ V w excepted. Taking advantage of axisymmetry and by analogy to the classical stress variables p " pΣ xx `Σyy `Σzz q{3

and q " Σ yy ´Σxx " Σ yy ´Σzz , the spherical and deviatoric part of any microstructure tensor µ are respectively quantified by the scalars m, Eq. ( 38), and t defined as:

t " c 3 2 ||µ ´mδ|| " µ yy ´µxx " µ yy ´µzz (44)
where, obviously:

t V w " 0 (45)
For this critical state configuration, deviatoric components are measured as expected for some microstructure tensors (Fig. 6), namely µ Ssw and µ Γ that refer to the wetted solid surfaces and contact lines, respectively.

The deviatoric nature of µ Ssw , as measured by t Ssw {m Ssw , can be compared to the fabric tensor F that describes the contact network:

F " 1 n c ÿ cont. n b n (46)
where n c is the total number of contacts between spherical particles, each of them being oriented by n. An exact equality between t F {m F " 3 t F and t Ssw {m Ssw would hold for a monosized solid packing under very high capillary pressure for which capillary bridges exist at contacts only and show negligible wetted surfaces. Here, the discrepancy between the two ratios narrows to 11% for u c " 40 kPa and S r « 0.063 % (Fig. 7). For higher saturations when cap- illary bridges between separated particles exist, the deviatoric nature of µ Ssw decreases while the anisotropy of the solid packing i.e. t F is barely affected by the hydraulic loading. In fact, deviatoric strains and changes in t F arise transitioning from dry to unsaturated conditions but become negligible upon subsequent hydraulic loading in triphasic conditions (see also [START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF].

Focusing now on the the fluid-fluid interface, the microstructure tensor µ Snw interestingly turns out to be spherical with a very good approximation. At the capillary bridge scale, a deviatoric nature exists for µ Snw (Fig. 8 and 9(a)), which is in relationship with an oriented meniscus shape. However, these devi- and 9(b)). The sphericity of µ Snw at the macro-scale has also been observed by [START_REF] Duriez | Stress in wet granular media with interfaces via homogenization and discrete element approaches[END_REF] for other numerical samples, still in the pendular regime. This sphericity of the nw interface was previously proposed as an assumption by [START_REF] Gray | The solid phase stress tensor in porous media mechanics and the Hill-Mandel condition[END_REF], without any possible measures to support it.

Finally, the capillary stresses depending on µ V w , µ Ssw , µ Snw , µ Γ deviate from an averaged fluid pressure (Fig. 10) showing a non-zero deviatoric component q cap defined as:

q cap " σ cap yy ´σcap xx " σ cap yy ´σcap zz " ´c 3 2 ||σ cap ´pcap δ|| (47)
For very low saturations, deviatoric capillary stresses arise mostly from the capillary pressure internal forces along the wetted solid surfaces, in relationship with the anisotropic solid packing as discussed in the above. For higher saturations, surface tension internal forces along the menisci contours contribute

Evolving solid packing

Finally, the case of a solid packing subject to important deformations is tackled considering the following mechanical loadings in triphasic conditions: two triaxial compressions under 4 and 40 kPa of capillary pressure and 10 kPa confining pressure. Such mechanical loading is considered to occur within smaller time scales than the hydraulic loadings previously considered under thermodynamic equilibrium. As such, the possibility to solve the Laplace-Young equation is complemented by additional conditions before considering capillary bridges.

Namely, capillary bridges are now assumed to form between contacting particles only, in accordance with observations by [START_REF] Herminghaus | Dynamics of wet granular matter[END_REF]. In case particles subsequently separate because of the mechanical loading, the menisci continue to be considered as long as the interparticle distance allows for a solution of the Laplace-Young equation. These processes of bridge creation and rupture occur in a fully drained manner, and during the mechanical loading, the saturation slightly evolves in r2.8%; 3.6%s and r0.06%; 0.09%s for a 4 and 40 kPa capillary pressure, respectively. Similar macroscopic behaviours are obtained during the two loading paths, see Fig. 11(a), the capillary stresses being a small portion of the total stresses (Fig. 11). This is due to the particle size distribution that confers the numerical sample a limited sensitivity to triphasic conditions: the latter would be greater for smaller particles.

Starting from an isotropic initial state, the solid phase turns anisotropic because of the deviatoric mechanical loading, as evidenced by the fabric tensor F , see Fig. 12(a). As observed during hydraulic loadings in the above, the wetted solid surfaces tend to follow the same anisotropy as the solid phase, especially for higher capillary pressure, see Fig. 12(a when the microstructure anisotropy is the most developed (Fig. 12).

Let us finally recall that a failure criterion encompassing stress limit states for granular materials in both dry and wet conditions can be expressed according to the contact stress tensor σ cont " Σ ´un δ ´σcap [START_REF] Hicher | A microstructural elastoplastic model for unsaturated granular materials[END_REF][START_REF] Scholtès | Micromechanics of granular materials with capillary effects[END_REF][START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF]Duriez and Wan, 2016a). For axisymmetric loading paths such as triaxial compressions, it is sufficient to consider the stress ratio q cont {p cont that obeys:

q cont p cont " q ´qcap p ´un ´pcap (48)
Obviously, predicting failure for a given external loading in terms of the total stresses p, q requires a precise estimation of both p cap and q cap . Neglecting in particular the deviatoric capillary stress q cap could lead to underestimating failure.

Conclusion

The present micromechanical analysis elucidates the nature of stresses in triphasic granular materials by accurately taking into account details of internal forces within the three bulk phases and the three associated interfaces. It is in particular demonstrated analytically that solid-fluid surface tensions do not contribute to the average stress whenever solid contact surfaces are negligible.

The proposed tensorial stress decomposition involves a capillary stress that is derived in terms of capillary pressure, fluid-fluid surface tension, as well as other tensorial quantities that depend on microstructure. An important implication of this expression is that the mixture of the two non-viscous fluids is able to sustain deviatoric stresses in triphasic media as demonstrated in DEM numerical experiments that precisely describe the micro-scale. Numerical experiments considered idealized spherical granular materials. If the particle shape of real granular materials were to be considered, it would be expected that the anisotropy of the capillary stress due to the fluid mixture would be even greater than in the present analysis.

Another major finding of the present work is that the DEM stress description from resultant forces acting on solid particles is equivalent in the pendular regime to the stress decomposition herein proposed from analytical homogenization of all phases and interfaces. The demonstration of the DEM stress description for granular media under low-saturation fully validates its use under such wet conditions. This also implies that possible simplifications in the realm of microstructure experimental measurements can be made by using Eq. ( 36) that expresses complex relationships between micro-scale parameters.

Contrary to alternate approaches based on thermodynamics, the proposed stress decomposition does not attach any constitutive meaning to the distinct stress contributions since only statical considerations are made without any kinematics. On the other hand, and with respect to [START_REF] Hassanizadeh | Mechanics and thermodynamics of multiphase flow in porous media including interphase boundaries[END_REF]; thermodynamics account for the microstructure of triphasic granular materials, which is believed to be of utmost importance. In fact, it is believed that tensorial free energy dependencies, such as the microstructure tensors proposed here, should be introduced in addition to scalar dependencies to capture the deviatoric nature of capillary stresses from thermodynamics.

Ongoing work is actually being pursued in order to assess the effective nature of σ cont " Σ ´un δ ´σcap from a constitutive point of view, and not restricted to failure description. In the meantime, the proposed stress decomposition

s ˜3 V m `ř i"1,2 R i S i,sw ¸`γ nw ˜2 S m `ř i"1,2 Γ i R i sin θ " ř i"1,2
s π R i 2 sin 2 δ i HO i `γnw Γ i HO i sin pδ i `θq (B.9)

Since s π R 1 2 sin 2 δ 1 `γnw Γ 1 sin pδ 1 `θq " s π R 2 2 sin 2 δ 2 `γnw Γ 2 sin pδ 2 `θq from Eq. ( 30), the r.h.s of (B.9) is equal to:

`s π R 1 2 sin 2 δ 1 `γnw Γ 1 sin pδ 1 `θq ˘ÿ i"1,2 HO i (B.10) Considering finally that ř HO i " R 1 `R2 `d, this readily demonstrates Eq. ( 36).

Fredlund ( 2006

 2006 );[START_REF] Madeo | A continuum model for deformable, second gradient porous media partially saturated with compressible fluids[END_REF];[START_REF] Duriez | Stress in wet granular media with interfaces via homogenization and discrete element approaches[END_REF] account for the fluid-fluid interface only, while Gray et al. (2009); Nikooee et al. (2013) include solid-fluid interfaces as well. The lack of consensus may arise from the fact that capillarity is best evidenced and understood for fluid-fluid interfaces e.g., in capillary tubes, leading even to controversies about solid-fluid interfaces in the fluid mechanics community

Figure 1 :

 1 Figure 1: Interface S αβ with n the normal and ν αβ the inwards conormal

Figure 2 :

 2 Figure 2: Contact line Γ: 3D (left) and 2D (right) views. With respect to Eq. (9), dl is a portion of Γ, and dS is an elemental surface area constructed from dl and an infinitesimally small length (magnified in the Figure) along t

(

  2006);[START_REF] Gray | The solid phase stress tensor in porous media mechanics and the Hill-Mandel condition[END_REF];[START_REF] Li | Meso-hydromechanically informed effective stresses and effective pressures for saturated and unsaturated porous media[END_REF]. Also, as a partial confirmation of this denomination,[START_REF] Hicher | A microstructural elastoplastic model for unsaturated granular materials[END_REF];[START_REF] Scholtès | Micromechanics of granular materials with capillary effects[END_REF];[START_REF] Wan | A tensorial description of stresses in triphasic granular materials with interfaces[END_REF];Duriez and Wan (2016a) considered the above-mentioned stress term to unify the failure description between triphasic and biphasic conditions.The second stress contribution in Eq. (18) describes the internal forces of a triphasic material due to the capillary pressure u c . A spherical stress p´u c V w {V δq is due to the fluid volumes and corresponding isotropic pressures.Another contribution encompasses the stress application onto the solid phase by the fluid phases through the suction surfacic loadings along the wetted solid surfaces: (´u c {V ş Ssw n b x dS). The latter contribution directly depends on the wetted solid surfaces orientations.

  under hydrostatic total stresses (with Σ xx " Σ yy " Σ zz " 10 kPa) and various capillary pressure, i.e. wetting saturation S r values. According to the model design presented back in Section 3.1, wetting saturation (capillary pressure) is uniformly imposed to the sample: capillary bridges are introduced solving the Laplace-Young equation between particles pairs for a given capillary pressure value. Doing so, capillary bridges are introduced between both contacting and separated particles, up to some interparticle distance beyond which the Laplace-Young equation does not give any solution. Due to the isotropic configuration of the solid phase, the resulting fluid phases and interfaces distributions are also isotropic and any microstructure tensor µ among µ V w , µ Ssw , µ Snw , µ Γ is spherical, i.e.

  Figure 4: Scalar representation of spherical microstructure tensors during the hydraulic loading of an isotropic sample

  Figure 5: Capillary stresses during the hydraulic loading of an isotropic sample

  Figure 6: Scalar representations of axisymmetric microstructure tensors during the hydraulic loading of a sample at critical state: spherical (m) and deviatoric (t) parts

Figure 7 :

 7 Figure 7: Deviatoric natures of the fabric tensor F and the wetted solid surfaces microstructure tensor µ Ssw

Figure 8 :

 8 Figure 8: Deviatoric nature of µ Snw for capillary bridges between two particles showing R 2 {R 1 " 1.1, considering various dimensionless capillary pressures u c and interparticle distance d ˚" d{R 2

Figure 9 :

 9 Figure 9: Orientation of the nw interface at different scales and for different capillary pressures: polar representation of 3{Snw vpαq.p ş Snw nb n dSq vpαq, v being an arbitrary unit vector oriented by the polar angle α

Figure 10 :

 10 Figure 10: Capillary stresses during hydraulic loading of a sample at critical state

Figure 12 :

 12 Figure 11: Triaxial compression loading paths in triphasic conditions

  [START_REF] Gray | The solid phase stress tensor in porous media mechanics and the Hill-Mandel condition[END_REF];[START_REF] Madeo | A continuum model for deformable, second gradient porous media partially saturated with compressible fluids[END_REF];[START_REF] Nikooee | Effective stress in unsaturated soils: A thermodynamic approach based on the interfacial energy and hydromechanical coupling[END_REF] who considered an elastic behaviour for the solid phase, the present result is not restricted to any particular mechanical behaviour of the different phases-except that fluid viscosity was neglected. Comparing the present work with thermodynamics approaches allows us to assess how well stress decompositions obtained from

Table 1 :

 1 DEM model parameters

	Contact k n { D k t {k n ϕ D min D max Packing	e	Capillarity γ nw θ
	(MPa)	(-)	( 0 ) (mm) (mm)	(-)	(N/m) ( 0 )
	10	0.21	20	0.25	0.3	0.61 0.073	10
	Three mechanical parameters, k

n { D, k t {k n , and ϕ, are relevant to the two constant local stiffnesses k n , k t and the local friction angle ϕ of the contact interaction, as detailed by

  already sheds further light on the old debate regarding the validity of Bishop's effective stress equation in partially saturated soils in that it identifies missing microstructural dependencies in Bishop's capillary stress equation.For each wetted surface S i,sw (i " 1, 2):´ÝÝÑ HO i . n `Ri n. n ¯dS " ´HO i π R i 2 sin 2 δ i `Ri S i,sw `Ri S i,sw ´HO i π R i 2 sin 2 δ i ˘¸´2 S m (B.5)Focusing now on the r.h.s of Eq. (B.1), we have: ÝÝÑ HO i `Ri nq dl " Γ i p´HO i sin pδ i `θq `Ri sin θq (B.7)Finally, the r.h.s of Eq. (B.1) is:

									S m	R 2
	R 1 O 1 δ 1 Equating (B.5) and (B.8) according to (B.1) gives: θ θ δ 2 O 2 S 2,sw S 1,sw H
								Figure 13: Meniscus half-geometry
	theorem with the classical divergence operator div gives:
			ż
	´uc γ nw	Sm	x. n dS
	"	´uc γ nw ˜żVm	div x dV	´ÿ i"1,2	ż Si,sw	x. n dS	¸;	n inwards to the particles
	"	´uc γ nw	˜3 V m	`ÿ i"1,2	ż Si,sw	x. n dS	¸;	n outwards to the particles (B.3)
	ż				ż
		x. n dS "
	Si,sw					Si,sw
									(B.4)
	The l.h.s of Eq. (B.1) is finally equal to:
		´uc γ nw	˜3 V m	`ÿ i"1,2
								ż	ż
									ν. x dl "	ν. x dl	`żΓ2	ν. x dl	(B.6)
								Γ	Γ1
	For each wetted contour Γ i (i " 1, 2):
	ż				ż
		ν. x dl "		ν.p
	Γi				Γi
					ż
						ν. x dl "	ÿ	Γ i p´HO i sin pδ i `θq `Ri sin θq	(B.8)
					Γ				i"1,2

Acknowledgements

Fruitful discussions with Rakulan Sivanesapillai (Institute of Mechanics, Ruhr-Universität Bochum, Germany) and Dr. Mehdi Pouragha (University of Calgary) are gratefully acknowledged. This work has been funded by the Natural Science and Engineering Research Council of Canada and Foundation Computer Modelling Group. Alonso, E., Pereira, J.-M., Vaunat, J., Olivella, S., 2010. A microstructurally based effective stress for unsaturated soils. Géotechnique 60, 913-925(12). Alonso, E. E., Gens, A., Josa, A., 1990. A constitutive model for partially saturated soils. Géotechnique 40 (3), 405-430. Biot, M. A., 1941. General theory of three-dimensional consolidation. Journal of Applied Physics 12 (2), 155-164. Bishop, A. W., 1959. The principle of effective stress. Teknisk Ukeblad 106, 859-863.

Appendix A: Demonstration of Eq. ( 16)

Stokes' theorem states the surface integral of the curl of an arbitrary vector field F over a surface S equals to the line integral of that vector field over the surface's boundary Γ. Adopting dl " h dl and the index notation with Einstein convention:

In Eq. (A.1), ǫ ijk is the Levi-Civita symbol. It is possible to replace the vector field F with a second-order tensor field W to get:

It is convenient here to define W from the Levi-Civita symbol, a second order tensor field w and a constant vector b: W ij " ǫ jnl w in b l . Then, the l.h.s of Eq. (A.2) is:

Using the identity:

we get finally for the l.h.s of (A.2):

On the other hand, for the r.h.s of Eq. (A.2) we have:

Because b k are arbitrary, they can be eliminated from Equations (A.2),(A.5),(A.6); which gives the following second order tensorial equality:

Now, choosing w " x b n and knowing that ǫ kni n n h i " ν k , Eq. (A.7) takes the form:

Also knowing that:

Equation (A.8) can finally be stated in the form:

which demonstrates Eq. ( 16). The lateral surface S m does not completely enclose the meniscus volume V m , however S m Y S 1,sw Y S 2,sw does (Fig. 13). Then, applying the divergence


Gerard Chalhoub
email: gerard.chalhoub@udamail.fr

Eric Perrier

La Bâthie

Michel Misson

Overhead caused by WiFi on ZigBee networks using slotted CSMA/CA

Keywords: coexistence, 802.11, 802.15.4, CSMA/CA, Zig-Bee, WiFi, MAC

This paper presents an in-depth analysis on the effect of the coexistence in the 2.4 GHz band between IEEE 802.15.4 and IEEE 802.11 on the medium access control. We focus on the slotted CSMA/CA medium access algorithm that is used by IEEE 802.15.4 and we analyse the delay that frames undergo when suffering from interference from a WiFi activity. We measure the overhead caused by the additional delay spent in the MAC sublayer frame queue before accessing the medium due to the presence of WiFi interference. We experimented different scenarios with overlapping channels and non-overlapping channels. We show that the two wireless protocols can coexist if we take into considerations the relative positions of the nodes to avoid very high interference and if we avoid monopolizing the channel with a very high rate WiFi traffic even under overlapping channels.

I. INTRODUCTION

Wireless networks are more and more deployed in different areas of applications and in our everyday life. The WiFi standard, which is based on IEEE 802.11 standard [START_REF]Part 11: Wireless LAN medium access control (MAC) and physical layer (PHY) specifications[END_REF], has been widely deployed in the past 25 years. The most used version of this standard uses the 2.4 GHz ISM band and is based on a random medium access control that allows different WiFi networks to coexist even when they are working on the same channel or overlapping channels. Nevertheless, it is highly recommended to use 5-channel separation to avoid interference because WiFi channels are overlapping as it is shown in figure 1. For example, channels 1, 6 and 11 in figure 1 are highlighted as possible 3 non-overlapping channels.

Another standard which has recently gained in popularity is ZigBee [START_REF] Zigbee | Zigbee Specification[END_REF]. It is based on IEEE 802.15.4 standard [START_REF]Part 15.4: Wireless medium access control (MAC) and physical layer (PHY) specifications for low-rate wireless personal area networks (WPANs)[END_REF]. ZigBee also uses the 2.4 GHz ISM band and is often deployed in areas where WiFi is already installed. Other more recently approved wireless standards such as WirelessHART [START_REF]HART field communication protocol specifications[END_REF] and ISA100.11a [START_REF]ISA100.11a: 2009 wireless systems for industrial automation: Process control and related applications[END_REF], that are proposed for industrial wireless sensor networks, are also based on the the physical layer of IEEE 802. [START_REF] Sikora | Coexistence of ieee 802.15.4 with other systems in the 2.4 ghz-ism-band[END_REF].4 and work in the 2.4 GHz ISM band. Figure 1 also shows how IEEE 802.15.4 channels overlap with those of the IEEE 802.11. WirelessHART and ISA100.11a are supposed to be deployed in more controlled environments such as indoor industrial monitoring and control applications and they are mostly based on TDMA (Time Division Multiple Access) and frequency hopping. Unlike ZigBee which is essentially based on CSMA/CA (Carrier Sense Multiple Access with Collision Avoidance) for accessing the medium.

In most cases, these networks need to coexist in the same physical space and share the same bandwidth in the 2.4 GHz ISM band. The increasing number of applications using wireless technology is creating congestion in the 2. [START_REF]HART field communication protocol specifications[END_REF] GHz frequency band and making the coexistence issue even more critical. Many studies have been done to evaluate the effect of the coexistence in this band as discussed in the state of the art of this paper.

In this paper, we study in great details the behaviour of slotted CSMA/CA of IEEE 802.15.4 in the presence on a WiFi network. Indeed, we investigate how WiFi affects the performance of slotted CSMA/CA algorithm in terms of additional delays caused on the backoff algorithm. According to our knowledge, this level of analysis has not been done before. We give an in-depth analysis of the behaviour of CSMA/CA according to the relative positions of IEEE 802.15.4 transmitters and receivers and WiFi nodes. We show that the effect of WiFi on the transmitters is not the same as on the receivers. In addition, we will show how the effect of WiFi on IEEE 802.15.4 depends on the traffic load of WiFi. We also measured the interference in overlapping channels and non-overlapping channels. We also propose a new overhead estimation caused by the coexistence. In what follows we will refer to IEEE 802.15.4 as ZigBee and IEEE 802.11 as WiFi.

The remainder of the paper is organized as follow. In section II we summarize some of the related work concerning coexistence in the 2.4 GHz band. In section III we describe how slotted CSMA/CA operates. In section IV we analyse the measurements from the different experiments that we did to evaluate to effect of WiFi on ZigBee and we present our overhead estimation. We conclude the paper in section V.

II. STATE OF THE ART

In the literature, many works have been done in this domain. These studies were limited to the Received Signal Strength Indicator (RSSI), the Packet Error Rate (PER), and the Link Quality Indicator (LQI).

In [START_REF] Yang | Wireless coexistence between ieee 802[END_REF], authors made a survey on the coexistence problem of these 2 wireless standards and discussed the state of the art of the different evaluations that exist. The paper discusses This work has been done with the collaboration of EDF R&D. We would like to thank Kevin Verny for his important contribution to the realization of the experiments and Andre Misson for providing the experimental site. the solutions and good practices to ensure none destructive coexistence based on frequency domain and time domain. These solutions for ZigBee standard, published in technical reports [START_REF] Company | CC2420 coexistence[END_REF], [START_REF]ZigBee and wireless frequency coexistence[END_REF], [START_REF]Co-existence and security whitepaper[END_REF], [START_REF]Avoiding RF interference between WiFi and Zigbee[END_REF], remain very costly such as reporting interference on the current channel and asking the ZigBee coordinator to order a channel switch for the whole network. The paper also points out the lack in the literature of evaluations that differ between transmitters and receivers of IEEE 802.15.4, others studied the effect of the relative angles between the links of both technologies.

In [START_REF] Rihan | On zigbee coexistence in the ism band: Measurements and simulations[END_REF], authors evaluated the effect of WiFi activity on ZigBee communications. Only the three following parameters RSSI, LQI and PER were evaluated. Results showed how the link quality degrade when the distance between the ZigBee communicating nodes increases. This was done for different channels of IEEE 802.15.4, different packet sizes for WiFi, and for both IEEE 802.11g and IEEE 802.11b. In addition, authors presented some simulation results using OPNET to investigate the end-to-end delay and the medium access delay. Authors concluded that cognitive radio should be enabled in ZigBee to cope with the coexistence.

Other studies have also been made by simulation such as [START_REF] Arif | Simulation and analysis of zigbee -wifi interference[END_REF] where authors concluded that throughput of ZigBee network decreases in the presence on a WiFi network. Evaluations where made on application level evaluations.

Authors in [START_REF] Wagh | Performance evaluation of IEEE 802.15.4 protocol under coexistence of wifi 802.11b[END_REF] conducted simulations and showed based on PER and BER results that using frequency agility can efficiently mitigate the effect of WiFi interference and enhance the performance of ZigBee networks.

In [START_REF] Shin | Packet error rate analysis of ieee 802[END_REF], authors studied the interference issue by analysis and simulation using OPNET. The only factor that was studied in this paper was the PER based on the collision timing and the Bit Error Rate (BER). The authors concluded that the activity of WiFi has very little effect on the the communication of IEEE 802.15.4 in two cases: if the WiFi interferer is more than 8 meters away from the IEEE 802.15.4 nodes, and if the frequency offset between the two technologies is bigger than 7 MHz.

In [START_REF] Sikora | Coexistence of ieee 802.15.4 with other systems in the 2.4 ghz-ism-band[END_REF], authors also evaluated the interference of WiFi on IEEE 802.15.4 activity in terms of PER. Results are based on real measurements using CC2420 Chipcon transceivers. In this paper authors tested a different type of Clear Channel Assessment (CCA) based only on the energy level but results were not compared to the normal CCA mode which takes into account the modulation conformity. They also evaluated the effect of Bluetooth and microwave ovens. Bluetooth caused bursty packet loss and microwave ovens have an impact only in very close range less than 1 meter away from the nodes.

In [START_REF] Angrisani | Experimental Study of Coexistence Issues Between IEEE 802.11b and IEEE 802.15.4 Wireless Networks[END_REF], authors also measured the coexistence effect on the PER by means of measurements. Consequences on both sides, that is on WiFi and on IEEE 802.15.4, were given. Results show that IEEE 802.15.4 is more damaged by the activity of WiFi than the other way around.

Other researchers such as [START_REF] Khanna | Eliminating the problem of coexistence between zigbee and wifi using EDCA synchronization and multi header transmission[END_REF], [START_REF] Valck | Exploiting programmable architectures for wifi/zigbee inter-technology cooperation[END_REF] and [START_REF] Jacob | Enabling coexistence of ZigBee and WiFi[END_REF] focused on new techniques or equipment for eliminating the interference effect but these techniques or equipment are not currently supported by the standards or provided in standard hardware architecture.

According to our knowledge, the study presented in this paper has never been done before.

III. IEEE 802.15.4 SLOTTED CSMA/CA

The IEEE 802.15.4 standard [START_REF]Part 15.4: Wireless medium access control (MAC) and physical layer (PHY) specifications for low-rate wireless personal area networks (WPANs)[END_REF] supports two operational modes: (i) the beacon-enabled mode in which periodic beacon frames are transmitted by the coordinator to synchronize nodes according to a superframe structure depicted in Fig. 2, and (ii) the non-beacon-enabled mode in which unslotted CSMA/CA is used. In this paper, we focus on the beaconenabled mode. Although the study could have been made on unslotted CSMA/CA, similar behaviour should be observed when unslotted CSMA/CA is used. In the beacon-enabled mode, the superframe is divided into two periods: the CAP (Contention Access Period) during which slotted CSMA/CA is used to avoid collisions, and the CFP (Contention Free Period) where the coordinator allocates time slots for nodes to ensure that nodes transmit frames in a TDMA manner in order to avoid collisions. The interval that separates two consecutive beacons (BI) and the superframe duration (SD) are determined by two parameters: the superframe order (SO) and the beacon order (BO). BI and SD are defined as follows:

BI = aBaseSuperframeDuration.2 BO , SD = aBaseSuperframeDuration.2 SO ,
where 0 ≤ SO ≤ BO ≤ 14 and aBaseSuperframeDuration = 15.36 ms.

Two topologies are proposed by the standard: star and peer-to-peer. The peer-to-peer topology suffers from beacon frame collisions when the beacon-enabled mode is used (see proposals of Task Group 4b of IEEE 802.15.4 or [START_REF] Koubaa | A time division beacon scheduling mechanism for IEEE 802.15.4/Zigbee cluster-tree wireless sensor networks[END_REF], [START_REF] Chalhoub | Cluster-tree based energy efficient protocol for wireless sensor networks[END_REF]). Thus, in what follows we only consider the star topology which covers the majority of deployed topologies.

Slotted CSMA/CA algorithm is applied before the transmission of a data frame during the Contention Access Period of the superframe. This algorithm is based on a 320 µs time unit called backoff period. The boundaries of the backoff periods are aligned with the boundaries of the superframe slots. Any activity of the MAC sublayer (such as the channel sensing and the transmissions) starts at the boundary of a backoff period.

Slotted CSMA/CA algorithm uses three variables: NB (Number of Backoffs), which is the number of times a backoff has been drawn for this transmission attempt, CW (Contention Window), which is the number of consecutive backoff periods during which a device senses the channel, and BE (Backoff Exponent) which defines the range of possible backoff periods a device waits for until it assesses the channel.

Figure 3 represents the steps of slotted CSMA/CA. In Step 1, the MAC sublayer initializes the three variables and locates the boundary of the next backoff period: NB is set to 0, CW is set to 2 and BE is set to 3. In Step 2, a random number of backoff periods is chosen from [0; 2 BE -1]. If the number of remaining backoff periods in the CAP is less than the chosen random number, the MAC sublayer draws a new backoff at the start of the CAP of the next superframe.

Once the backoff has been decremented, the MAC sublayer makes sure that the remaining number of backoff periods in the CAP is enough to perform two CCAs (since CW = 2), to send the frame and to receive the acknowledgment (if required). If this is the case, it proceeds to the next step. Otherwise, it waits for the start of the CAP of the next superframe and goes back to Step 2. This is called the CCA deference.

In Step 3, the MAC sublayer asks the physical layer to perform a CCA at the next backoff boundary. The next step depends on the result of the CCA. If the channel is assessed to be busy, the MAC sublayer goes to

IV. RESULTS OF COEXISTENCE

We conducted different experiments in order to measure the effect of WiFi on ZigBee. All experiments were done using Cisco Aironet access points with certified WiFi clients, and Texas Instruments (TI) wireless sensor nodes that implement IEEE 802.15.4 slotted CSMA/CA with a CC2420 transceiver [START_REF] Company | CC2420 coexistence[END_REF]. WiFi traffic was generated without nearby IEEE 802.11b stations, RTS/CTS were not used, application layer acknowledgements were not used, and ERP-OFDM modulation with a short preamble was used. We made modifications on the TI code in order to obtain the different timing elements. TI nodes were connected to a computer using a serial cable and timing indicators for the node activity were logged on the computer.

We studied the effect of WiFi duty cycle on the packet loss ratio of ZigBee in IV-A. In IV-B, we argue that a good link quality has better chances to resist to interferences, and in IV-C and IV-D we evaluate the overhead in terms of delay while accessing the channel due to interferences. Details about each experiment as well as metrics that are evaluated are provided in each section.

A. Effect of WiFi duty cycle

In what follows, we evaluated the effect of the WiFi traffic profile on the packet loss ratio of ZigBee. The goal behind this series of measurements is to find the threshold at which the WiFi traffic has no or little effect on the ZigBee traffic. We considered that 1% of packet loss is enough to judge that WiFi has an effect. Thus, we plotted the RSSI values of interference for each duty cycle when an effect is detected.

Indeed, low duty cycles represent low channel occupancy by WiFi. A duty cycle of 35% for example means that the WiFi channel is occupied only 35% of the time. We made measurements with the following WiFi traffic profiles:

• a duty cycle of 2,25% with UDP traffic at a 54 Mbps modulation,

• a duty cycle 5% with UDP traffic at a 24 Mbps modulation,

• a duty cycle 17,5% with UDP traffic at a 6 Mbps modulation,

• a duty cycle 23% with UDP traffic at a 24 Mbps modulation,

• a duty cycle 35% with UDP traffic at a 9 Mbps modulation,

• a duty cycle 44% with UDP traffic at a 12 Mbps modulation,

• a duty cycle 59% with TCP traffic at a 54 Mbps modulation,

• a duty cycle 88% with TCP traffic at a 6 Mbps modulation, With these different WiFi traffic profiles we tested three scenarios: overlapping channels (WiFi on channel 11 and ZigBee on channel 23), adjacent channels (WiFi on channel 11 and ZigBee on channel 25), and non-overlapping channels (WiFi on channel 11 and ZigBee on channel 26).

Note that, when we a WiFi duty cycle of 88% has an effect on ZigBee at a given level of interference (-70 dBm for example), we conclude that there is no need to do additional tests for higher levels of interference. We generated only one ZigBee packet per second, and an acknowledgement is requested for each packet. Each message is 30-byte long.

Figures 5, 6 and7 show the threshold that should not be exceeded in terms of received signal power observed at the receiver and the transmitter sides to avoid packet loss on the ZigBee traffic when WiFi and ZigBee channels are overlapping, adjacent and non-overlapping respectively. Figure 5. Threshold of RSSI of interference level that defines the limit not to be exceeded in order for WiFi to have no effect on ZigBee in terms of packet loss: overlapping channels. Figure 6. Threshold of RSSI of interference level that defines the limit not to be exceeded in order for WiFi to have no effect of ZigBee in terms of packet loss: adjacent channels.

Results show that when the interference is applied on the receiver, the effects on packet loss is more significant and start at lower RSSI values compared to the case where the interference is applied on the transmitter. Indeed, this is an expected result because on the transmitter side the node applies CSMA/CA and hope to get a chance to find the medium idle to send its packet. Once it starts transmitting the packet, it will not stop even if the channel becomes occupied. This is not the case for the receiver which needs to have lager windows of idle channel (channel non occupied by the WiFi traffic) in order Figure 7. Threshold of RSSI of interference level that defines the limit not to be exceeded in order for WiFi to have no effect of ZigBee in terms of packet loss: non-overlapping channels.

to successfully receive the packet during the whole period of reception. Thus, the risk of losing packets is higher on the receiver side. For example, in figure 5, if we consider a WiFi duty cycle of 23%, ZigBee receivers should not experience an interference level above -75 dBm in order to avoid packet loss.

B. Link quality effect

In the following scenario, we show that with a high RSSI level, the ZigBee link resists much better to WiFi perturbation. In order to do so, we maintained the same level of interference coming from WiFi traffic and we varied the RSSI level on the ZigBee link. Small RSSI values represent a low quality link and high RSSI values represent a high quality link.

We put a WiFi access point at a 10-meter distance from the ZigBee receiver. The access point is transmitting a TCP traffic at 100 mW of transmit power using a 6 Mbps bit rate. The received signal power of WiFi at the ZigBee receiver is around -62 dBm. The WiFi channel is channel 11 and ZigBee channel is channel 23. We varied the distance between the ZigBee nodes between 50 cm and 13 m. The Zigbee transmitter is sending 1 message every 9 ms, each message is 110-byte long. Results of figure 8 show that the ZigBee link can resist to perturbations from WiFi as long as the RSSI of the WiFi signal is not more than 10 dB above the RSSI value of the ZigBee signal. This phenomenon is related to what we call the capture effect in a wireless transmission [START_REF] Whitehouse | Exploiting the capture effect for collision detection and recovery[END_REF].

C. Overhead of the coexistence

In what follows, we study the overhead of the coexistence in terms of additional delay that the ZigBee packets endure before accessing the medium. Indeed, we calculated the average delay for accessing the medium without the presence of interference and we considered that all additional average delays when a WiFi interference is present is an overhead generated by the coexistence. These delays would come from longer backoffs for successful transmissions, time spent on unsuccessful transmissions and repetitions. Thus, we computed the overhead according to the following formula in order to evaluate each of the additional delay sources: Overhead = (ADwP -ADwoP) + AAF*(NbF/NbS) + AAS*((NbS-NbR)/NbS), where

• ADwP (Average Delay with Perturbation) is the average medium access delay with the presence of WiFi interference,

• ADwoP (Average Delay without Perturbation) is the average medium access delay without WiFi interference,

• AAF (Average Access Failure) is the average time for an access failure to occur,

• NbF (Number of Failures) is the number of frames that were not transmitted due to failure to access the medium,

• NbS (Number of Successes) is the number of frames that were successfully transmitted,

• AAS (Average Access Success) is the average time for accessing the medium without errors,

• NbR (Number of Received) is the number of received frames. Where,

• (ADwP -ADwoP) will give us the average additional delay due to interference from WiFi.

• AAF*(NbF/NbS) will give us the additional average time spent on trying to send frames that were not sent due to interference.

• AAS*((NbS-NbR)/NbS) will give us the additional average time spent on repeating frames due to interference. In these scenarios, the ZigBee node that is being the object of interference will be suffering from an interference level at -65 dBm and the other node at -80 dBm in the first case. We have made the tests with 2 ZigBee traffic profiles. On one hand, a relatively low ZigBee traffic with one message of 30 bytes every second, each message is acknowledged and repeated 3 times in case the acknowledgement is not received. This traffic profile represents typical notification messages sent by sensors to a control center. And on the other hand, a very high traffic rate with one message of 30 bytes every 9 ms, messages are not repeated, but acknowledgements are kept for calculating the statistics. This profile represents the maximum transmission data rate with acknowledgements that we were able to generate using our motes without causing synchronization errors.

Results in figure 9 show that with the presence of WiFi interference, each ZigBee packet suffers between 30 ms and 55 ms of additional delay due to access delay, repetitions and packet loss. In addition, it is clear to notice that the overhead is higher when the interference is more significant on the transmitter. This is essentially due to the fact that when the transmitter is under high interference, it will struggle to send a frame and will keep waiting for an idle channel to be able to send it. This will cause significant delay before accessing the medium.

In addition, a high rate ZigBee traffic is more likely to suffer from interference because it is more likely to find the medium idle to receive or send packets when the number of packets is very high.

We also tested the high data rate profile with UDP WiFi traffic with 1 Mbps and 5 Mbps over a 6 Mbps modulation. As figure 10 shows, UDP generates less interference than TCP and the 1 Mbps generates less interference than the 5 Mbps. The overhead in additional delay reflects the consequences of these interferences.

D. Overhead caused by interference with non-overlapping channels

In figure 11 we show the overhead on the ZigBee activity when we use channels 20 to 26 while WiFi is active on channel 6 with an interference level at around -30 dBm. Notice that interference on the transmitter has significantly more effect on the overhead than interference on the receiver, this is due to delays caused by extended backoffs and unsuccessful transmissions due to busy channel detections. On the other hand, when interference is applied to the receiver, overhead is mainly caused by collisions on the receiver that produce more repetitions.

We did a similar experiment with WiFi on channel 11 and ZigBee on channels 25 and 26. We started with an interference level at -55 dBm and then -65 dBm. Results show significant overhead on channel 25 and on channel 26 when the interference level is at -55 dBm as shown in figure 12. In complementary results, with an interference level at -65 dBm there was no overhead on neither channels.

V. CONCLUSION AND PERSPECTIVES

In this paper, we evaluated the effect of the coexistence between two wireless network standards IEEE 802.11 and IEEE 802.15.4. We showed that ZigBee is effected by WiFi on many levels: frame loss, frame repetitions, access delay, and consequently, energy consumption. But this is very related to the channel used in both standards. The effect is more significant when dealing with overlapping or adjacent channels.

We made a very detailed analysis on the overhead caused by the presence of interference. We measured the additional delay and actions that a node might undergo before being able to correctly access the channel. The results of overhead can directly be converted to energy consumption overhead. Indeed, these operations are causing the nodes to consume more energy under the presence of interference.

Results showed that using ZigBee and WiFi in the same environment is possible if special precautions are taken into account. For example, placing the ZigBee nodes in places where the RSSI level of WiFi is below -65 dBm when possible will significantly protect ZigBee traffic from interference. Choosing non-overlapping channels is always encouraged. And if possible, reduce WiFi activity in such a way to avoid occupying the channel with very high duty cycles; duty cycles below 33% are very much encouraged.

Results in this paper only presented evaluations on link levels. In future works, experimentations should be made on network scale taking into account application scenarios. According to our results, application level results will very much depend on the level of exposure of the network to the interference and the importance of exposed links in relation to the application.

Figure 1 .

 1 Figure 1. The overlap between WiFi channels from 1 to 13 and the ZigBee channels from 11 to 26 in the 2.4 GHz spectrum.

Figure 2 .

 2 Figure 2. Beacon interval in IEEE 802.15.4 beacon-enabled mode.

CW = 2 , 1 NB

 21 Figure 3. IEEE 802.15.4 slotted CSMA/CA (without the battery life extension).

  exceeds macMaxCSMABackoffs (denoted by macMaxCB. on the figure), the algorithm terminates with a channel access failure, otherwise it goes back to Step 2. In Step 5, CW is decremented. If CW reaches zero, the transmission begins at the boundary of the next backoff period. Otherwise, the MAC sublayer returns to Step 3.

Figure 4 .

 4 Figure 4. Equipment used for the experiments. One Cisco Aironet access point and two TI ZigBee nodes. WiFi clients were used in order to generate traffic according to the different scenarios used in the experiments. Relative positions of the access point and ZigBee nodes change according to each experiment.

Figure 8 .

 8 Figure 8. Percentage of packet loss at the ZigBee receiver with different RSSI values of the ZigBee signal.

Figure 9 .Figure 10 .

 910 Figure 9. Overhead: average additional delay due to interference with TCP traffic.

Figure 11 .Figure 12 .

 1112 Figure 11. Overhead with very high interference and non-overlapping channels.