Analysis of stress partitioning in biphasic mixtures based on a variational purely-macroscopic theory of compressible porous media: recovery of Terzaghi’s law

IS - Computational Models and Methods for Multiphysics Processes in Multiphase Porous MediaAnalysis of stress partitioning in biphasic mixtures based on a variational purely-macroscopic theory of compressible porous media: Recovery of Terzaghi‘s law VI International Conference on Computational Methods for Coupled Problems in Science and Engineering COUPLED PROBLEMS 2015 B. Schrefler, E. Oñate and M. Papadrakakis(Eds) ANALYSIS OF STRESS PARTITIONING IN BIPHASIC MIXTURES BASED ON A VARIATIONAL PURELY-MACROSCOPIC THEORY OF COMPRESSIBLE POROUS MEDIA: RECOVERY OF TERZAGHI’S LAW ROBERTO SERPIERI∗, FRANCESCO TRAVASCIO†, SHIHAB ASFOUR†AND LUCIANO ROSATI+ ∗Dipartimento di Ingegneria, Università degli Studi del Sannio Piazza Roma, 21 - I, 82100, Benevento, Italy. e-mail: [email protected] †Biomechanics Research Laboratory, Department of Industrial Engineering, University of Miami 1251 Memorial Drive, MCA268, Coral Gables, FL, 33146. e-mail: f.trav[email protected], [email protected] +Dipartimento di Strutture per l’Ingegneria e l’Architettura (DIST), Università di Napoli Federico II Via Claudio, 21, 80125 Napoli, Italy. e-mail: [email protected] Key words: Variational poroelasticity, Effective stress, Terzaghi’s law, Consolidation, Porous media Abstract. The mechanics of stress partitioning in two-phase porous media is predicted on the basis of a variational purely-macroscopic theory of porous media (VMTPM) with compressible constituents. Attention is focused on applications in which undrained flow (UF) conditions are relevant, e.g., consolidation of clay soils and fast deformations in cartilagineous tissues. In a study of the linearized version of VMTPM we have recently shown that, as UF conditions are approached (low permeability or fast loading), Terzaghi’s effective stress law holds as a general property of rational continuum mechanics and is recovered as the characteristic stress partitioning law that a biphasic medium naturally complies with. The proof of this property is obtained under minimal constitutive hypotheses and no assumptions on internal microstructural features of a particular class of material. VMTPM predicts that such property is unrelated to compressibility moduli of phases and admits no deviations from Terzaghi’s expression of effective stress, in contrast with most of the currently available poroelastic theoretical frameworks. This result is presently illustrated and discussed. Simulations of compressive consolidation tests are also presented; they are obtained via a combined analytical-numerical integration technique, based on the employment of Laplace transforms inverted numerically via de Hoog et al.’s algorithm. The computed solutions consistently describe a transition from drained to undrained flow which confirms that Terzaghi’s law is recovered as the limit UF condition is approached and indicate a complex mechanical behavior. 1 35 R. Serpieri, F. Travascio, S. Asfour and L. Rosati 1 INTRODUCTION The analysis of the dynamic response of multiphase porous media in the limit of undrained flow conditions is relevant in all those applicative contexts where the nature of characteristic loadings, compared to the characteristic consolidation time of the medium, determines a flow regime close to the ideal conditions of complete prevention of fluid drainage. Canonical examples of applications in which this condition is met are the analysis of saturated clay soils subjected to seismic loading [1], or building-induced short time static loading, and the analysis of cartilaginous tissues subjected to physiological impulsive loading [2]. Among several approaches so far proposed for continuum modelling of multiphase porous media, variational approaches [3–5] provide a tool for addressing poroelastic multiphase problems by introducing the least possible number of mechanical assumptions and postulated balance laws. We have recently proposed a least-action based macroscopic continuum description of twophase poroelasticity to derive a general biphasic formulation at finite deformations based upon the inclusion of the intrinsic volume variations of the solid [6,7] among the kinematic descriptors. This theory is shortly referred to as Variational Macroscopic Theory of Porous Media (VMTPM). In VMTPM no Lagrange multipliers are employed to obtain the macroscopic local balance equations. The main consequence is that work-association between stress and strain measures is naturally preserved, and the local macroscopic Euler-Lagrange momentum balance equations include a momentum balance associated with intrinsic volumetric strains. For the subclass of undrained flow (UF) problems in which macroscopic fluid redistribution within the mixture is impeded and inertial forces can be neglected, the linearized version of VMTPM predicts that stress is partitioned in the two phases in strict compliance with Terzaghi’s law. This relation is found to hold irrespective of thermodynamic constraints, constitutive or microstructural features of the medium, and independently from intrinsic compressibility properties of phases. Furthermore, such a property, provided UF conditions are met, admits no deviations from Terzaghi’s expression of effective stress, in contrast with most of the currently available poroelastic theoretical frameworks. The objective of this study is to complement the general law of stress partitioning in UF, derived in [8], with solutions of the general 1D consolidation problem. The objective is the analysis of the transition from a drained to an undrained behavior so as to investigate the behavior of systems which admit superposition of solutions, and for which solution techniques based on Laplace transforms can be exploited. Accordingly, the governing equations herein considered are those pertaining to the linearized VMTPM in which elastic stiffness moduli are introduced in a standard purely-variational form as the second derivatives of strain potentials. The constitutive parameters employed in the numerical simulations are selected in compliance with the bounds reported in [9,10]. 2 36 R. Serpieri, F. Travascio, S. Asfour and L. Rosati 2 VMTPM IN QUASI STATIONARY FLOW The purely mechanical linearized version of VMTPM is considered herein. Accordingly, all forces are represented by a potential so that Hamilton’s principle, written for the solid porous phase, takes the form: δt2 t1T(s)−U(s)dt =0.(1) In the previous relation U(s)is the total potential energy of the porous solid body, T(s)is the kinetic energy and t1and t2are two generic times. In linearized kinematics, the infinitesimal deformation of the solid phase is defined, at the macroscale, by the solid macroscopic infinitesimal displacement field ¯ u(s), and by the macroscopic field of the intrinsic volumetric strain ˆe(s). The latter is related to the microscopic volumetric strain field of the pore-scale deformation, e(s), by the standard averaging relation: ˆe(s)(x)= 1 V(s)Ω(s) RV E(x) e(s)dV =dV (s) V(s).(2) where V(s)and dV (s)are the volumes of the subset of the representative volume element centered in x, and its first-order variation, respectively. The addition of ˆe(s)as a kinematic descriptor beside displacements is the unique essential kinematic/contitutive enhancement introduced with respect to standard Cauchy continuum theory. The reader is referred to [6–8] for a more extensive description of the VMTPM. Hereby, the kinematics and the linearized set of governing equations are briefly recalled. Specifically, employing the notation conventions and terminology used in [8], we recall from [7] the specialization of the linearized VMTPM, which stems by considering negligible inertia forces. This condition will be shortly termed henceforth Quasi Stationary Flow (QSF). The extrinsic volumetric strain, ¯e(s), is ordinarily related to the macroscopic solid displacement field as the trace of the solid strain tensor: ¯e(s)=tr¯ ε(s)=∇·¯ u(s),(3) where ¯ ε(s)=sym ¯ u(s) ∇, being ¯ u(s) ∇=¯ u(s)⊗∇the displacements gradient. The governing equations of linearized VMTPM at QSF are synoptically recalled in Table 1 where the notation employed in [8, 10] is used. Equations (4) and (5) are, respectively, the so-called extrinsic momentum balance of the solid phase and the intrinsic momentum balances, both stemming from (1). These two equations express stationarity of the Action functional with respect to variations of the solid displacements and variations of ˆe(s), respectively. Equation (7) represents the saturation constraint in linearized kinematics (see [8], and Bedford and Drumheller [3]), written in dimensionless form for a biphasic compressible medium. Accordingly, pis the fluid pressure while volumetric fractions φ(α) o(with α=s, f), due to the saturation hypothesis fulfill φ(f) o+φ(s) o=1. Moreover, ˇ σ(s)and ˆp(s), termed ’solid extrinsic stress tensor’ and ’solid intrinsic scalar pressure’, are the primary stress measures of the solid 3 37 R. Serpieri, F. Travascio, S. Asfour and L. Rosati Extrinsic linear momentum balance of solid phase ∇·ˇ σ(s)+πfs =0 (4) Intrinsic momentum balance of the solid phase ˆp(s)+ˆp(fs)=0 (5) Momentum balance of the fluid phase −φ(f) o∇p+πsf =0 (6) Combined fluid mass balance and saturation constraint (Dimensionless saturation constraint) φ(f) o¯e(f)+φ(s) o¯e(s)=φ(f) oˆe(f)+φ(s) oˆe(s)(7) Table 1: General set of equations that govern the mechanical behavior of the biphasic compressible poroelastic system in linearized VMTPM. phase workassociated with ¯ ε(s)and ˆe(s), respectively while ˆp(fs)=−∂¯ ψ(f) ∂ˆe(s)where ¯ ψ(f)denotes the strain energy density of the fluid phase. Mutual drag body volume forces acting over the solid and the fluid phase are respectively denoted by πfs and πsf , with πfs =−πsf . The linear theory in terms of elastic coefficients is derived introducing the elastic moduli in a standard form, as second order derivatives of strain potentials with respect to the primary strain measures [9, 10]. Compared to the general treatment, the QSF theory admits a simpler description in which only two additional moduli appear along with Lamé elastic moduli. These two moduli are a dimensionless coefficient ¯ krwhich is characteristic of VMTPM and couples extrinsic strains with intrinsic ones, plus an auxiliary intrinsic stiffness modulus, ˆ ks. They are defined as follows: ¯ kr=φ(s) o ∂2¯ ψ(s) ∂ˆe(s)∂¯e(s)∂2¯ ψ(s) ∂ˆe(s)∂ˆe(s)−1 ,ˆ ks=1 φ(s) o ∂2¯ ψ(s) ∂ˆe(s)∂ˆe(s).(8) Employing the coefficients defined in (8), a linear combination of the linear isotropic constitutive laws of the solid phase with the intrinsic momentum balance can be conveniently introduced in QSF, [10]: ˇ σ(s)= 2¯µ¯ ε(s)+¯ λ¯e(s)I−¯ krpI,φ(s) o ˆ ks p=−¯ kr¯e(s)−φ(s) oˆe(s)(9) Moreover, in QSF a combination of governing balances (4)-(7) yields a u-pform of the governing equations in which the primary variables are the solid macroscopic displacements and the fluid pressure, [9, 10]. This form turns out to be more practical from an engineering computational point of view, although limited to QSF problems. Upon introducing a simple linear 4 38 R. Serpieri, F. Travascio, S. Asfour and L. Rosati Darcy law for the drag body volume forces πfs: πfs =−πsf =K∂¯ u(f) ∂t −∂¯ u(s) ∂t (10) (with Kmeasured in Ns m4) the u-pform reads: ∇·σ(s) D−φ(f) o+¯ kr∇p=0 (11) −φ(f) o2 K(∇·∇)p+1+¯ kr∂ ∂t∇·¯ u(s)+1 ˆ ksf ∂p ∂t =0 (12) where σ(s) D=¯ λ¯e(s)I+ 2¯µ¯ ε(s)is the drained solid stress and ˆ ksf =φ(s) o ˆ ks+φ(f) o ˆ kf−1 , being ˆ kf the bulk modulus of the fluid phase. 3 STRESS PARTITITONING IN UNDRAINED FLOW The UF condition refers to the flow of a biphasic medium contained in a closed region of space ΩbU ⊆Ωbsuch that macroscopic fluid mass flow across the boundary ˜ Sof any closed region contained in ΩbU is zero. Conditions of UF are met when 1) the characteristic consolidation time is much higher than the timescale of observation of the problem; 2) when drainage is prevented at the boundaries and/or in presence of special symmetry conditions (such as those typically encountered in jacketed tests). The characterization of UF conditions from a kinematic point of view corresponds to the existence of a unique macroscopic continuum field of undrained displacements ¯ u(un), common to both solid and fluid phases, defining the macroscopic displacements of the undrained medium: ¯ u(un)(x)=¯ u(f)(x)=¯ u(s)(x),x∈ΩbU .(13) The dual static characterization of (13) has been determined, in linearized kinematics, in [8]. Upon deriving the characteristic forms, Uint bU ,Uext bU , achieved at UF by the internal and eternal energy potentials, the Principle of Virtual Deformations is applied in the framework of linearized VMTPM. The primary difference with respect to the derivation of the principle of virtual work for single phase continuum mechanics is that, in accordance with the presence in VMTPM of two macroscopic state fields ¯ u(un)and ˆe(s), two variational integral conditions (describing equilibrium in static UF configurations) are accordingly inferred. These equations are associated with virtual macroscopic displacements δ¯ u(un)and virtual intrinsic strain variations δˆe(s): ∂(¯ u(un),ˆe(s))Utot bU δ¯ u(un),δˆe(s)=0,∀δ¯ u(un),δˆe(s).(14) Derivation of the explicit strationarity integral equation stemming from (14) yields explicitly the following fundamental mechanical identification which holds within UF conditions: ∂¯ ψ(f) ∂¯e(un)=−p. (15) 5 39 R. Serpieri, F. Travascio, S. Asfour and L. Rosati As shown in [8], conversion from weak to strong form of the equilibrium equations yields a set of local equations which exactly match the classical tensorial statement of Terzaghi’s principle. Specifically, these equations state that on the boundary of a region undergoing undrained flow, tractions, t(ext), applied by the external environment are related to internal stress measures, ˇ σ(s)and p, by t(ext)=ˇ σ(s)−pIn,(16) where — it is worth being recalled — ˇ σ(s)is the stress tensor naturally work-associated with isochoric (volume preserving) strains of the solid phase, consistently with its role of ’effective’ stress. This relation holds irrespective of thermodynamic constraints, of the constitutive or microstructural features of the medium and independently from intrinsic compressibility of phases. 4 TRANSITION FROM DRAINED TO UNDRAINED FLOW The transition from a more general condition of drained flow to UF conditions has been investigated, as predicted by VMTPM, under conditions of quasi stationary flow for a simple 1D consolidation problem. A semi analytical numerical solution method based on Laplace transforms was applied to solve the QSF problem of a generic biphasic specimen subjected to uni-axial confined compression in a stress-relaxation test. The examined problem is the following. A biphasic sample is laterally and inferiorly confined in an impermeable chamber, and compressed by a porous plug allowing for fluid exudation. The setup considered is schematized in Figure 1. Figure 1: Simulated experimental setup for uni-axial confined compressive test: (a) biphasic mixture is confined in an impermeable chamber and axially compressed by a porous plug allowing for fluid exudation; (b) compressive displacement history applied to the plug for static tests (U) and associated boundary conditions. 6 40 R. Serpieri, F. Travascio, S. Asfour and L. Rosati A dimensionless treatment of this problem is exploited and tilde accents are used to indicate dimensionless quantities. Denoting by Land to, respectively, the characteristic sample length and time scale of observation, we introduce the following dimensionless quantities: ˜p=p (¯ λ+2¯µ)˜σ(s) Dxx =σ(s) Dxx (¯ λ+2¯µ),˜u(s) x=¯u(s) x L,˜x=x L,˜ t=t to (17) so that the 1D counterpart of u-pequations (11), (12) can be recast in the form: ∂2˜u(s) x ∂˜x2(˜x, ˜ t)−(φ(f) o+¯ kr)∂˜p ∂˜x(˜x, ˜ t)=0 (18) (1 + ¯ kr)∂2˜u(s) x ∂˜ t∂˜x(˜x, ˜ t)+Cr ∂˜p ∂˜ t(˜x, ˜ t)−1 De ∂2˜p ∂˜x2(˜x, ˜ t)=0 (19) Coefficients De and Crare dimensionless characteristic parameters of the system, defined by: De =τ to ,τ=L2(φ(f) o)2 K(¯ λ+2¯µ),C r=¯ λ+2¯µ ˆ ksf (20) where τis proportional to the consolidation time of the system. The physical meaning of coefficient De is closely related to the Deborah number [13], since τis a characteristic time proportional to the classic poromechanics definition of relaxation time (i.e. the time required for the biphasic medium to adapt to applied stresses or deformation). The range of this coefficient is 0≤De ≤∞. When De >> 1the system tends to behave as an elastic completely undrained medium. The opposite condition of completely drained behavior corresponds to De << 1while, for De ≃1, the fluid component of the mixture provides a significant contribution to time dependent stress/deformation phenomena. The dimensionless quantity Crrepresents the ratio of the aggregate modulus of the solid phase in the mixture (i.e., ¯ λand ¯µ) over the intrinsic stiffness (ˆ ksf ) of the medium. Accordingly, when Cr→0, the mixture can be assumed to be intrinsically incompressible. Space-time boundary conditions of the problem of Figure 1 are also directly examined in dimensionless form. At time t=0the mixture is in an undeformed state and in equilibrium with the environment so that the initial conditions are: ˜u(s) x=0,˜p=0 ˜ t=0,˜x∈[0,1] (21) At the bottom of the specimen (i.e., ˜x=0), solid displacement is null, and no fluid flows through the inferior wall of the chamber. This condition amounts to a null solid-fluid relative displacement and, recalling (6) and (10), implies a null pressure gradient. At the top (˜x= 1), fluid exudation occurs through the porous plug. In addition, for the case of displacement controlled compression, the displacement applied to the plug follows a ramp-and-hold time sequence, where the dimensionless ramp time is ˜ tR, and the final compressive displacement is 7 41 R. Serpieri, F. Travascio, S. Asfour and L. Rosati −Uo, (with Uo>0) see figure 1b. Based on these assumptions, the boundary conditions of the problem are collected below: ˜u(s) x(0,˜ t)=0,∂˜p ∂˜x(0,˜ t)=0,˜u(s) x(1,˜ t)= ˜ U(˜ t)=−Uo L˜ t−H(˜ t−˜ tR),˜p(1,˜ t)=0 (22) where His the Heaviside function ˜ U=U/L and ˜ tR=tR/to. The absolute value of the displacement applied to the porous plug, finally attained in the ramp-and-hold time sequence, is Uo=0.01 L(i.e., 1 percent of the total length of the sample) while the ramp time is ˜ tR=0.1 (see figure 1b). Solutions to the boundary value problem composed of equations (18), (19), (21), (22) are first determined in the complex Laplace space. Expressing (18), (19) in terms of Laplace transforms, denoted by star superscripts, yields: ∂2(˜u(s) x)∗ ∂˜x2(˜x, s)−(φ(f) o+¯ kr)∂˜p∗ ∂˜x(˜x, s)=0 (23) s(1 + ¯ kr)∂(˜u(s) x)∗ ∂˜x(˜x, s)+sCr˜p∗(˜x, s)−1 De ∂2˜p∗ ∂˜x2(˜x, s)=0 (24) Denoting by ˜ U∗the Laplace transform of ˜ U(˜ t), the transforms of the dimensionless solid displacement, of the solid displacements gradient and of the fluid pressure turn out to be: (˜u(s) x)∗(˜x, s)= ˜ U∗β1 α√s sinh(α√s˜x) cosh(α√s)+ (1 −β)˜x β1 α√stanh(α√s)+1−β(25) ∂¯u∗ x ∂˜x(˜x, s)= ˜ U∗βcosh(α√s˜x) cosh(α√s)+1−β β1 α√stanh(α√s)+1−β(26) ¯p∗(˜x, s)= ˜ U∗cosh(α√s˜x) cosh(α√s)−1 (φ(f) o+¯ kr)1 α√stanh(α√s)+ Cr (1+¯ kr)(27) where parameters α,βare defined as follows: α=De(Cr+γ),β=γ Cr+γ,γ=(φ(f) o+¯ kr)(1 + ¯ kr).(28) Further mathematical details for obtaining the solution of equations (23) and (24) are reported in [11]. Except for special cases, Laplace anti-transforms of (25)-(27) cannot be performed analytically retaining a closed form. Thus, they were carried out by numerical computation via de Hoog et al’s algorithm [15]. A special case is represented by the UF limit at De →∞. In such a limit the antitrasform of (27) to the time domain is easily obtained analytically. This solution 8 42 R. Serpieri, F. Travascio, S. Asfour and L. Rosati consistently corresponds to strain and pressure fields that are uniform in space. In particular the constant value attained by the pressure field, written in dimensional form, is: p=−1+¯ krˆ ksf Uo(t) L(29) and corresponds to the general undrained solution, compliant with Terzaghi’s law (16), whose full expression is reported in [10]. For solutions outside of the UF limit, De =∞, the observation of the system is carried out over a value of the dimensionless time equal to 1 and is obtained numerically. Concerning the selection of the values of the dimensionless coefficients employed in the simulations (φ(f) o,¯ λ+2¯µ,¯ kr,ˆ ksfor the solid, and ˆ kffor the fluid phase), these were set so as to describe the behavior of a generic biphasic medium, without referring to a particular medium of a specific application. In absence of an experimental characterisation of these parameters, special attention has been payed in assigning to φ(f) o,¯ λ+2¯µ,¯ kr,ˆ ksconsistent values preserving physical admissibility and meaningfulness. In this respect, it is important to observe that, while ˆ kfis unrelated to the parameters pertaining to the solid phase (φ(f) o,¯ λ+ 2¯µ,¯ kr,ˆ ks) the latter four parameters cannot be independently assigned since they all belong to the solid phase. To recognize this it can be easily verified that, if these values are improperly selected so that the condition Cr+γ>0is violated, the solution achieves a singularity with unbounded values for βand imaginary values for α. It can be shown that, by considering standard Mori-Tanaka estimates for the macroscopic shear modulus and Hashin’s Composites Spheres Assemblage (CSA) estimates for ¯ kr,ˆ ksand ¯ λ, the range of variation of ¯ kris −1≤¯ kr≤0[9, 10] and that the quantity (Cr+γ)is strictly positive [11]. Conversely, as φ(f) oincreases, function βturns from negative to positive values. In the example of Figure 3 below it is shown that this change of sign determines significant qualitative and quantitative differences in the strain profiles. For De =1and different values of φ(f) o,¯ krand Cr(details on the employed data are reported in [11]) the shape of the curves of stress history at the plug is the same. Stress increases throughout the ramp phase of the displacement and reaches a peak value at ˜ tR. Subsequently, while plug displacement is held fixed, stress relaxes to an equilibrium value corresponding to the condition in which no relative motion between solid and fluid phases occurs (see Figure 2). To appreciate the role played by ¯ krin determining the character of the solution, sensitivity analyses of the space-time solutions have been carried out numerically [11] by employing dimensionless coefficients provided by CSA and Mori-Tanaka estimates [9–11]. Herein two solutions are presented from those considered in [11]; they have been obtained with the two sets of dimensionless parameters of Table 2, which in turn correspond to two different values of the porosity φ(f) o=0.25 (A) and φ(f) o=0.5(B). Figure 3 shows, at several dimensionless times, the computed strain profiles. While the final stationary states of cases Aand Bare coincident, significant qualitative and quantitative differences are detected in the transient response. Details of these analyses and an extensive discussion of the related results are reported in [11]. 9 43