The partially saturated porous medium is described here adopting the classical approach of continuum poromechanics as presented by Coussy (2004), to which we refer for a detailed presentation; the notation used hereafter is consistent with that of the aforementioned reference. According to this classical approach to partial saturation, a porous solid is understood as a so-called ‘wetted’ porous skeleton with a thin layer of fluid attached to the pore walls, and thus a contribution to the stored energy of the porous solid is attributed to the fluid–fluid and solid–fluid interfaces. This interfacial energy density is obtained as the cumulative work carried out by a macroscopic capillary pressure difference between the wetting and the nonwetting fluids in modifying the relative fluid volume fractions within the pores.

Now, the formation of a fracture due to grain/pore reorganization, as hypothesized by Shin and Santamarina (2010, 2011), results in a discontinuity of the porous material. This “absence” of porous material within a fracture would mean, locally, a complete loss of retention property toward the wetting fluid and thus a loss of the ability of the porous medium to store interfacial energy. In the current study, a scalar damage variable, α, is introduced in order to describe the effects of grain/pore reorganization within the porous solid, in the regime of compressive effective stresses, as briefly discussed in the Introduction. Damage is thus assumed to induce a reduction in the air-entry pressure and in general to degrade the retention properties of the porous medium, which is equivalent to a reduction of the cumulative interfacial energy stored.

Accordingly, the degraded retention relation between the saturation degree of the wetting fluid, S _w , and the degraded macroscopic capillary pressure, p _c , reads, p c S w , α = a α p ̃ c S w ≈ − a α p ̃ w S w = − p w S w , α , (1) where a(α) is a damage law whose functional form is elaborated at a later point. p _w is the pore-pressure of the wetting fluid in a degraded porous solid, and p ̃ w its nondegraded counterpart. p ̃ c ( S w ) represents the standard retention relation, unaffected by damage. In the current study, the widely used van Genuchten form (van Genuchten, 1980), p ̃ c S w = ρ w g α v G S w n 1 − n − 1 1 n , (2) is assumed where, α _vG and n are the model parameters, g is the acceleration due to gravity, and ρ _w is the mass density of liquid water. It is to be noted that residual saturation is assumed to be vanishing for simplicity but its incorporation is trivial.

Also, the assumptions leading to Richards’ equation (Hilfer and Steinle, 2014) are adopted, in particular the hypothesis of a passive non-wetting phase is assumed, leading to p _c ≈ −p _w in Eq. 1. While the above development serves as a reasonable first assumption for investigation, a complete modeling approach should also involve the degradation of the resistance to fluid flow within the damaged zone and eventually along localized fracture planes. With the abovementioned assumptions, the transient hydraulic problem governing the evolution of pore-pressure, p _w , within the porous solid is written as, ∂ ϕ S w ∂ t − ϰ η w ∇ ⋅ K r w p w ∇ p w = 0 , (3) where the intrinsic permeability ϰ is assumed to be independent of α and η w is the dynamic viscosity of water. Since ϰ is expected to increase, in orders of magnitude, within a fracture, we make a rather simplifying assumption on the relative permeability function, K _rw (p _w ) = 1 to in part compensate for the lack of the aforementioned coupling. It is to be noted that since the primary unknown factor in the aforementioned equation is the pore-pressure, S _w can then be derived from it using the degraded inverse retention relation, S w p w , α = p ̃ c − 1 − p w a α . (4)

Appropriate boundary conditions compliment the aforementioned partial differential equation (PDE) which could be either an imposed flux representing a natural boundary condition or an imposed pressure representing a Dirichlet boundary condition.

We consider an n _d -dimensional domain, Ω ∈ R n d , representing the initial configuration of a damaging isotropic partially saturated porous solid whose boundary is denoted by ∂Ω. In line with standard gradient damage modeling (Marigo et al., 2016) and the theory of unsaturated poroelasticity, the following considerations are carried out.

1) The scalar internal variable, α ∈ [0, 1], is assigned the role of describing the extent of damage in the sense of loss of retention properties. α = 0 denotes an intact healthy skeleton whose retention properties are described by the standard retention curve, whereas α = 1 denotes a fully damaged skeleton whose retention properties are degraded and can only hold vanishing or residual amounts of wetting fluid.

The different terms of the aforementioned expression are explained below.

a) The energy density of the porous solid, Ψ _s (ɛ, ϕ, S _w , α), encompasses the contributions due to the deformation of the porous skeleton, resulting in solid strains and changes in porosity and, also, the interfacial energy contribution due to formation/annihilation of interfaces. Adopting the concept of energy separation (Coussy, 2004) between the free energy of the porous skeleton and the interfacial energy we obtain,

Here, a specific choice is made that would result in the aforementioned postulate that damage means the degradation of only the retention properties and not the poroelastic properties. The latter would have been the choice for modeling tensile mode fracture as carried out, for instance, by Cajuhi et al. (2018). The choice made here allows to isolate the investigation of fracture formation just due to the ‘release’ of interfacial energy. Specifically, we assume a degradation of the form, U S w , α = a α U ̃ S w , (7) consistent with the definition of the degraded capillary pressure in Eq. 1. In Figure 1, the two functions p _c (S _w , α) and U(S _w , α) are depicted with the assumed degradation behavior. U ̃ ( S w ) represents the undamaged interfacial energy given by, U ̃ S w = ∫ S w 1 p ̃ c s d s . (8)

While various choices of the functional form of a(α) are possible, we make a specific choice following the standard gradient damage modeling, a α = 1 − α 2 + k ℓ , (9) where k _ℓ is a small positive constant used solely for numerical purposes to govern the fully damaged state. Accordingly, when α grows from 0 to 1, the interfacial energy density degrades to a vanishing value. The consequence on the damage evolution due to the abovementioned choices is discussed further.

Now concerning the free energy density of the skeleton, a possible form that retrieves back the classical Biot’s linear constitutive relations as its state equations reads as, ψ s ε , ϕ = 1 2 C 0 ε − ε 0 ⋅ ε − ε 0 + 1 2 b 2 N ϵ − ϵ 0 2 − Δ ϕ b N ϵ − ϵ 0 − p 0 ∗ + 1 2 N Δ ϕ 2 , (10) where the notations are borrowed from Coussy (2004). C 0 is the elastic stiffness tensor, b the Biot coefficient, and N the tangent Biot modulus. Δϕ represents the increment of Lagrangian porosity, (ϕ − ϕ ₀). ϕ ₀, ɛ ₀, and p 0 ∗ denote the initial reference states of the Lagrangian porosity, the linearized strain tensor, and the average pore-pressure, respectively.

b) The second term in Eq. 5 is the local part of the dissipated energy,

Now, we pose the evolution problem for the solution triplet (u, ϕ, α) of the poromechanical damage problem using a variational approach similar to the developments in the studies by Pham and Marigo (2010a, 2010b). The principles governing the evolution are elaborated for the current problem concerning the porous solid.

The variational approach is based on the definition of a suitable total energy of the body under consideration (i.e., the porous solid), that is associated with the triplet of admissible states ( v , φ , β ) ∈ C × P × D . These functional spaces are defined as, C = v ∈ H 1 Ω n d : v = 0 on ∂ Ω D , P = φ ∈ L 2 Ω : 0 < φ < 1 in Ω , D = β ∈ H 1 Ω : 0 ≤ β < 1 in Ω , (13) where L ²(Ω) is the Lebesgue space of square-integrable functions equipped with an L ²-norm and H ¹(Ω) is the Sobolev space of square-integrable functions whose weak gradients are also also square-integrable. ∂ Ω D is a subset of ∂ Ω where displacement, u , is imposed.

This definition of total energy at each time t is given as the integral over the whole domain of the bulk energy density minus the potential of external forces, E t v , φ , β = ∫ Ω W t v , φ , β , ∇ β ; S t d x − W t e . (14)

Since the focus here is on the porous solid, the external loading concerns not only bulk forces and tractions exerted on the solid skeleton, which in this case are assumed to vanish for the sake of simplicity, but also zeroth-order bulk actions tending to modify from inside the porosity of the solid. This time-dependent bulk action is obviously related to the presence of fluids in the porous network and is, therefore, coupled with the solution of the hydraulic problem, governed by Eq. 3. In a similar way as in the case of fully saturated porous materials, working this action of the change of Lagrangian porosity, it will coincide with the average pore-pressure, p t ∗ = p t S t + p n w S n w , for all t > 0. Pore pressure, p _w , and saturation degree, S _w , of the wetting fluid at time t are further denoted p _t and S _t , respectively, in order not to overload the notation. The assumption of a passive nonwetting phase leads to p t ∗ ( S t , α t ) ≈ p t S t . It is to be noted that the dependence of p t ∗ on damage is due to the coupling of the fluid problem to the damage evolution within the domain. This contribution to the potential due to external efforts defined on the set of admissible porosity fields, φ ∈ P , reads, W t e φ = ∫ Ω p t ∗ φ d x . (15)

Remark: We assume that all the fields are sufficiently smooth in time, allowing for the following developments. Also, we consider only the states of damage where α < 1 and the saturation degree S _w > 0. This is done so that in the following analysis, the total energy remains finite.

Now, we are in a position to setup the principles of the variational approach for the evolution in Ω of ( u t , ϕ t , α t ) ∈ C × P × D for all t ≥ 0 which read:

1) Irreversibility of damage: t↦α _t must be nondecreasing. Consequently, the admissible states accessible from α _t belong to the space D t defined as,

Equilibrium equation and damage criterion can be obtained as first-order stability conditions starting from Eq. 17. In addition, since the solution now involves ϕ _t , an associated zeroth-order balance equation is expected to appear.

We start again by expanding the perturbed energy up to the second-order, h E t ′ u t , ϕ t , α t v − u t , φ − ϕ t , β − α t + h 2 2 E t ″ u t , ϕ t , α t v − u t , φ − ϕ t , β − α t 2 + o h 2 ≥ 0 . (19)

E t ′ and E t ′ ′ represent, respectively, the first and second directional derivatives of E t further referred to as FDD and SDD, respectively, for compactness. The representation E t ′ ′ ( u t , ϕ t , α t ) ( ⋅ ) 2 is to be understood as a shorthand for the quadratic form, whereas the associated symmetric bilinear form is represented E t ′ ′ ( u t , ϕ t , α t ) ⟨ ⋅ , ⋅ ⟩ that is, the application of E t ′ ′ ( u t , ϕ t , α t ) to the pair of directions ⟨⋅, ⋅⟩. The directional derivatives in Eq. 19 have the following forms in the general direction ( v ̂ , φ ̂ , β ̂ ) , E t ′ u t , ϕ t , α t v ̂ , φ ̂ , β ̂ = ∫ Ω ∂ W t ∂ ε ⋅ ε v ̂ + ∂ W t ∂ ϕ − p t ∗ φ ̂ + ∂ W t ∂ α − p t ∗ ′ ϕ t β ̂ + ∂ W t ∂ ∇ α ⋅ ∇ β ̂ } d x , (20) E t ″ u t , ϕ t , α t v ̂ , φ ̂ , β ̂ 2 = ∫ Ω C 0 ε v ̂ ⋅ ε v ̂ + N b ϵ v ̂ − φ ̂ 2 − ϕ t π t ″ β ̂ 2 − 2 π t ′ φ ̂ β ̂ + w 1 ℓ 2 ∇ β ̂ ⋅ ∇ β ̂ } d x , (21) where for the sake of compactness of the notation, we denoted the directions ( v − u t , φ − ϕ t , β − α t ) with ( v ̂ , φ ̂ , β ̂ ) and define their associated function spaces as C × P ̂ × D with, P ̂ = φ ̂ ∈ H 1 Ω : − 1 < φ ̂ < 1 in Ω . (22)

In Eq. 20, the partial derivatives of bulk energy density are functions of state (u _t , ϕ _t , α _t ) given by Eq. 12.

Dividing Eq. 19 by h and passing to the limit h → 0 gives, E t ′ u t , ϕ t , α t v − u t , φ − ϕ t , β − α t ≥ 0 , (23) which is the so-called first-order stability condition and can be viewed as characterizing stationarity of the state (u _t , ϕ _t , α _t ). In Eq. 23, testing with φ = ϕ _t and β = α _t and noting that C is a linear space, one can obtain the variational (weak) form of the classical equilibrium condition: ∫ Ω ∂ W t ∂ ε ⋅ ε v − u t d x = 0 , ∀ v − u t ∈ C , (24) where a stress tensor can obviously be identified as, σ t = ∂ W t ∂ ε = C 0 ε t − ε 0 + b N b ϵ t − ϵ 0 − Δ ϕ t I . (25)

Similarly, testing with v = u _t and β = α _t in Eq. 23, one obtains the zeroth-order balance equation associated with small variations in ϕ _t , ∫ Ω ∂ W t ∂ ϕ − p t ∗ φ − ϕ t d x = 0 , ∀ φ ∈ P , (26)

Finally, using Eq. 24 and Eq. 26 in Eq. 23 gives the variational form of the non-local damage criterion ∫ Ω ∂ W t ∂ α − p t ∗ ′ ϕ t β − α t + ∂ W t ∂ ∇ α ⋅ ∇ β − α t d x ≥ 0 , ∀ β ∈ D t (27)

Using classical localization arguments of calculus of variations, one can obtain from Eq. 24 and Eq. 26 and Eq. 27 the following local forms, respectively, with corresponding boundary conditions, ∇ ⋅ σ t = 0 in  Ω , σ t ⋅ n ̂ = 0 on  ∂ Ω \ ∂ Ω D , (28) − b N ϵ t − ϵ 0 + p 0 ∗ + N Δ ϕ t + U α t ; S t = p t ∗ in  Ω , (29) − π t ′ ϕ t + w ′ α t − w 1 ℓ 2 Δ α t ≥ 0  in  Ω , ∂ α t ∂ n ̂ ≥ 0    on    ∂ Ω , (30) where n ̂ denotes the outward unit normal vector associated with part of the boundary wherever invoked. Note the absence of surface and volume forces in the equilibrium equation according to earlier assumption.

Equation 29 is the local form of the zeroth-order balance law associated with variations in porosity. Rearranging it, one can identify the relation which in classical poromechanics is characterized as the constitutive relation for porosity, ϕ t = b ϵ t − ϵ 0 + 1 N p t ∗ − p 0 ∗ − U α t ; S t + ϕ 0 . (31)

Remark: The aforementioned formulation results in an equivalent equilibrium equation to be resolved in tandem with the zeroth-order balance law for variations in porosity. This can be seen by substituting Eq. 31 into the equilibrium equation, Eq. 28, giving, ∇ ⋅ C 0 ε − ε 0 − b p t ∗ − p 0 ∗ − U α t ; S t I = 0 in  Ω . (32)

In the abovementioned formula, the classical stress tensor for unsaturated poroelasticity can be identified as, σ t = C 0 ε t − ε 0 − b p t ∗ − p 0 ∗ − U α t ; S t I . (33)

This equivalence between the formulations is exploited further for the purpose of numerical approximation.

Since we have assumed that the evolution is smooth in time, taking a time derivative of the Energy balance leads to, 0 = d d t E t u t , ϕ t , α t + ∫ Ω ϕ t ∂ π t ∂ S t S ̇ t d x = ∫ Ω ∂ W t ∂ ε ⋅ ε u t ̇ + ∂ W t ∂ ϕ − p t ∗ ϕ t ̇ + ∂ W t ∂ α − p t ∗ ′ ϕ t α t ̇ + ∂ W t ∂ ∇ α ⋅ ∇ α t ̇ d x . (34)

Integrating by parts the gradient terms and further using the local form of the equilibrium equations Eq. 28 and the zeroth-order balance law, Eq. 29 gives, 0 = ∫ Ω ∂ W t ∂ α − p 0 ∗ ′ ϕ t − ∇ ⋅ ∂ W t ∂ ∇ α α t ̇ d x + ∫ ∂ Ω ∂ W t ∂ ∇ α ⋅ n ̂ α t ̇ d x . (35)

Owing to the irreversibility of damage everywhere in Ω and the local inequalities Eq. 30, the two integrals on the right-hand side of the abovementioned equation are positive or equal to zero. So, both of them should vanish. Further using classical localization arguments, we obtain the so-called consistency conditions or the Karush–Kuhn–Tucker (KKT) conditions applicable everywhere in Ω and on the boundary ∂Ω respectively, − π t ′ ϕ t + w ′ α t − w 1 ℓ 2 Δ α t α t ̇ = 0  in  Ω , ∂ α t ∂ n ̂ α t ̇ = 0  on  ∂ Ω . (36)

These conditions can be read in tandem with the Irreversibility of damage. The first condition states that everywhere in Ω, the damage increases only if the local form of the damage criterion Eq. 30 is an equality, and if it is a strict inequality, then damage does not increase. The second condition states that everywhere on the boundary ∂Ω, if damage increases then the spatial derivative normal to the boundary vanishes.

The reference initial configuration is a rectangular domain Ω = (0, L) × (0, H) as shown in the Figure 2 with the boundary ∂Ω = {(x ₁ = 0) ∪ (x ₁ = L) ∪ (x ₂ = 0) ∪ (x ₂ = H)}. At time t = 0, it is assumed that the porous solid is completely intact, stress-free, and fully saturated giving, α 0 x = 0 , u 0 x = 0 , ε 0 x = 0 , S 0 x = 1 , p 0 x = 0 . (39) x ∈ Ω is the position vector given by x = x ₁ e ₁ + x ₂ e ₂. The horizontal and vertical displacements are set to vanish, respectively, on the lateral (x ₁ = 0, x ₁ = L) and the bottom (x ₂ = H) faces along with the shear stresses at those boundaries. The top face at x ₂ = 0 is a free surface. These set of mechanical boundary conditions ∀t ≥ 0 read, u t ⋅ − e 1 | x 1 = 0 = 0 , u t ⋅ e 1 | x 1 = L = 0 , u t ⋅ e 2 | x 2 = H = 0 , e 2 ⋅ σ t ⋅ − e 1 | x 1 = 0 = 0 , e 2 ⋅ σ t ⋅ e 1 | x 1 = L = 0 , σ t ⋅ − e 2 | x 2 = 0 = 0 , e 1 ⋅ σ t ⋅ e 2 | x 2 = H = 0 . (40)

Damage is not prescribed and is free to evolve at all the boundaries whenever the damage criterion is met. Accordingly, the natural boundary condition on damage reads, ∂ α t ∂ n ̂ | x ∈ ∂ Ω = 0 ∀ t ≥ 0 . (41)

For the hydraulic problem, the lateral and bottom faces are set to be impermeable. The loading is carried out by a constant imposed outward flux, − q _f e ₂, on the top face, thus inducing a drying effect. While the drying flux measured using the discharge rates in Stirling (2014) is found to be a function of time, we choose to make a simplifying assumption of constant flux. These boundary conditions ∀t ≥ 0 read, ∇ p t ⋅ − e 1 | x 1 = 0 = 0 , ∇ p t ⋅ e 1 | x 1 = L = 0 , − ϰ η w K r w p t ∇ p t ⋅ − e 2 | x 2 = 0 = q f , ∇ p t ⋅ e 2 | x 2 = H = 0 . (42)

It can be seen from the problem setup that the geometry, material properties, and initial conditions are invariant along the x ₁-direction. The loading and boundary conditions on the top and bottom faces are as well-invariant along the x ₁-direction. The boundary conditions on the lateral faces demand vanishing gradients of the solution along the x ₁-direction. So, even though the problem does not render itself for obtaining an easy exact solution due to its nonlinearity, one can a priori anticipate the existence of a particular class of solutions that are homogeneous in the x ₁-direction and are dependent only on x ₂. Such solutions are termed here as base solutions. In fact, to obtain the base solutions, one just needs to solve the problem in one-dimension along the x ₂-direction with appropriate boundary conditions, instead of the full two-dimensional problem posed earlier. This is the purpose that is explained in this section. Rigorous mathematical proof of the existence and uniqueness of such base solutions for all times is out of the scope of the current work. Instead, we exhibit numerically these solutions for t > 0.

The domain of the one-dimensional problem is defined Ω ̃ = ( 0 , H ) along the x ₂-coordinate, with boundary ∂ Ω ̃ = { ( x 2 = 0 ) ∪ ( x 2 = H ) } . The initial and boundary conditions from Section 4.1 are adapted to this domain. While studying the influence of depth H could be interesting, we choose a depth such that H ≫ ℓ.

In view of the damage criterion, Eq. 30, with the definition of average pore-pressure, p t ∗ ( S t , α t ) = p t ( S t , α t ) S t , the solution states can be classified into two: undamaged ( α t = 0 ∀ x 2 ∈ Ω ̃ ) and damaged ( ∃ x 2 ∈ Ω ̃ ∣ α t > 0 ) . Since the initial state at t = 0 is assumed to be that of a uniformly intact solid, damage would initiate if the following local damage criterion is an equality. ϕ t U ′ α t ; S t − p t ′ S t ϕ t + w ′ α t ≥ 0  in  Ω ̃ , (43) where Eq. 1 gives − p t ′ = p c ′ ( S t , α t ) = a ′ ( α t ) p ̃ c ( S t ) and Eq. 7 gives U ′ ( α t ; S t ) = a ′ ( α t ) U ̃ ( S t ) . In Eq. 43, the terms that drive the damage evolution are first and second. Now for α _t = 0, a′(α _t ) = −2, and for 0 < S _t < 1, U ̃ ( S t ) > 0 and p ̃ c ( S t ) > 0 . Since ϕ t ∈ P , the first driving term is negative, whose magnitude increases with decreasing S _t and conversely with increasing α _t . With a similar reasoning also, the second term is negative. The third term is a positive constant, w′(α _t ) = w ₁, acting as the threshold for damage initiation. Owing to the drying flux applied at x ₂ = 0, the saturation degree is the lowest here. Thus, one can anticipate that starting from a uniformly undamaged state, the damage initiates at x ₂ = 0 within a finite depth and propagates into the domain driven by desaturation. However, at initial times, since a minimum amount of capillary energy needs to be stored before damage initiation according to Eq. 43, the saturation degree everywhere in the domain could be such that the damage criterion is a strict inequality and damage does not initiate at all. So, one can envisage a finite time t _d > 0 beyond which the damage criterion is an equality at the drying boundary and so damage initiates. All the solution states for t ≥ t _d are called the damaged states, and this t _d depends on the intensity of the drying flux for a given material.

Figures 3, 4 show the time evolution of the solution and functions of the solution, with H = 1.0m, q _f = 10m.s⁻¹, and w ₁ = 1000N.m⁻³. The material properties correspond to sand and the fluid combination to air–water. It can be seen that initially, damage does not evolve anywhere in the domain even though flux-induced desaturation initiates at x ₂ = 0 and progresses into the domain. During this stage of evolution, the strain remains compressive and the porosity reduces from the initial value as expected within the unsaturated zone. However, once the damage criterion becomes equality, damage initiates at x ₂ = 0. With damage growing from 0, the saturation degree degrades close to 0, indicating invasion of the air phase while the pore-water pressure only reduces slightly in magnitude. The paths followed by the solutions for different thresholds, w ₁, in the space (S _t , p _t ) are shown in Figure 5. It is clear that the path followed during the initial damage evolution is that of a ‘softening’ with respect to pore-water pressure. It is also interesting to observe that even though the drying flux at the boundary continues to drive the desaturation process, the damage stagnates at a maximum value at larger times, and this is replicated also within the domain, thus creating a zone of uniform damage. This is due to the simultaneous reduction in the magnitude of the first driving term in Eq. 43 with increasing damage so that after a certain value of damage, the criterion does not become an equality anymore, even though desaturation progresses.

Figure 6 shows the time evolution of damage for different thresholds, w ₁. In all cases, the damage initiates at x ₂ = 0 as anticipated and then propagates into the domain supported by a finite depth, d _t , that increases with time. Lower thresholds are characterized by earlier initiation of damage and a damage value closer to 1 at the drying boundary.

The one-dimensional base solutions resolved in Section 4.2 are the x ₁-homogeneous solutions or the so-called fundamental states of the corresponding two-dimensional problem. These solutions are denoted by ( u ̄ t , ϕ ̄ t , α ̄ t ) at time t and, for consistency, the solution of the hydraulic problem by ( S ̄ t , p ̄ t ) . The branch of the evolution along which all the states remain x ₁-homogeneous is called the fundamental branch. Following the works of Benallal and Marigo (2006); Pham et al. (2011b); Sicsic et al. (2014) for standard gradient damage models, we investigate the possibility of bifurcation from the fundamental branch of evolution. A bifurcation in this sense is understood as the availability of admissible solution state/s other than the fundamental one. Whether or not the solution shifts from the fundamental branch depends on the loss of stability of the fundamental state itself and the characteristics of the perturbation that can cause such a shift. In the earlier mentioned works concerning standard gradient damage models, the possibility of bifurcations was associated with a loss of uniqueness criterion of the fundamental state. We follow a similar approach by first analyzing the loss of stability of the fundamental state and then the possibility of bifurcation.