2.1. Model problem

We consider a one-dimensional sample of soft porous material of relaxed length $L$ and relaxed porosity (fluid fraction) $\phi _{f,0}$. The left boundary of the material is permeable and located at $x=a(t)$ (moving); the right boundary is impermeable and located at $x=L$ (fixed in place) (figure 1). We impose a periodic, displacement-driven loading via the position of the left boundary,

(2.1)\begin{equation} a(t)= \frac{A}{2} \left[1-\cos\left(\frac{2{\rm \pi} t}{T}\right)\right], \end{equation}

where $A$ and $T$ are the amplitude and period of the loading, respectively. Macroscopically, the deformation is strictly compressive in the sense that $a(t)\geq {}0$.

Figure 1. We consider a one-dimensional sample of soft porous material of relaxed length $L$, subject to a periodic, displacement-driven loading at its left boundary (white arrows). The left boundary is permeable, thus allowing fluid flow in or out (blue squiggles) to accommodate the loading. The right boundary is impermeable and fixed in place.

We consider deformations ranging from small to large macroscopic nominal strain ($-0.4\,\%$ to $-20\,\%$ or $0.004\leq {}A/L\leq {}0.2$), with commensurate macroscopic changes in bulk (total) volume. We assume that the fluid and the solid phases are individually incompressible, so that the total volume of solid is constant and any change in bulk volume must correspond to a change in total pore (fluid) volume via a rearrangement of the pore structure and an influx or an efflux of fluid at the left boundary.

2.2. Kinematics

We work in an Eulerian reference frame, in which the solid displacement field is $\boldsymbol {u_s}=\boldsymbol {x}-\boldsymbol {X}(\boldsymbol {x},t)$, with $\boldsymbol {X}(\boldsymbol {x},t)$ the reference position of the material point that at time $t$ occupies the position $\boldsymbol {x}$. We choose $\boldsymbol {X}(\boldsymbol {x},0)\equiv \boldsymbol {x}$ so that $\boldsymbol {u_{s}}(\boldsymbol {x},0)=0$, in which case the reference configuration is the relaxed configuration. We denote the true volume fractions of fluid and solid by $\phi _f$ and $\phi _{s}$, respectively, where $\phi _{f}+\phi _{s}=1$. The flow and deformation are uniaxial, such that

(2.2a–c)\begin{equation} \boldsymbol{u_s} = u_s(x,t)\boldsymbol{\hat{e}_x}, \quad \boldsymbol{v_s} = v_s(x,t)\boldsymbol{\hat{e}_x}, \quad \boldsymbol{v_f} = v_f(x,t)\boldsymbol{\hat{e}_x}, \end{equation}

where $u_s$, $v_s$ and $v_f$ are the $x$-components of the solid displacement field and the solid and fluid velocity fields, respectively, and $\boldsymbol {\hat {e}_x}$ is the unit vector in the $x$-direction.

In one dimension, the local state of deformation is fully characterised by the Jacobian determinant $J= (1-\partial {u_s}/\partial {x})^{-1}$, which measures the local current volume per unit reference volume. For incompressible constituents and uniform initial porosity $\phi _{f,0}$, the local change in volume is linked to the change in porosity according to

(2.3)\begin{equation} J(x,t) = \frac{1-\phi_{f,0}}{1-\phi_f} \quad\to\quad \frac{\partial{u_s}}{\partial{x}} =\frac{\phi_f-\phi_{f,0}}{1-\phi_{f,0}}. \end{equation}

Continuity for this one-dimensional system can be written

(2.4a,b)\begin{equation} \frac{\partial{\phi_f}}{\partial{t}} +\frac{\partial}{\partial{x}}(\phi_f v_f) = 0 \quad\mathrm{and}\quad \frac{\partial{\phi_f}}{\partial{t}} -\frac{\partial}{\partial{x}}{[(1-\phi_f)v_s]} = 0, \end{equation}

which together imply that the total flux $q=\phi _fv_f+(1-\phi _f)v_s$ is uniform in space, $\partial {q}/\partial {x}=0$.

2.3. Darcy's law

We assume that the movement of the fluid relative to the solid skeleton is described by Darcy's law,

(2.5)\begin{equation} \phi_f(v_f-v_s) ={-}\frac{k(\phi_f)}{\mu}\frac{\partial{p}}{\partial{x}}, \end{equation}

where $k(\phi _f)$ is the permeability of the solid skeleton, $\mu$ is the dynamic viscosity of the fluid and $p$ is the fluid (pore) pressure, and where we have neglected gravity. We have taken the permeability $k(\phi _f)$ to be a function of porosity only. For simplicity, we use a normalised Kozeny–Carman relation for deformation-dependent permeability (MacMinn et al. Reference MacMinn, Dufresne and Wettlaufer2016):

(2.6)\begin{equation} k(\phi_f) = k_0\frac{(1-\phi_{f,0})^2}{\phi_{f,0}^3} \frac{\phi_f^3}{(1-\phi_f)^2}, \end{equation}

where $k_0\equiv {}k(\phi _{f,0})$ is the reference permeability. Kozeny–Carman permeability is a common choice for gels and soft tissues (Malandrino et al. Reference Malandrino, Lacroix, Hellmich, Ito, Ferguson and Noailly2014; Sacco et al. Reference Sacco, Causin, Zunino and Raimondi2014; Rahbari et al. Reference Rahbari, Montazerian, Davoodi and Homayoonfar2017; Gao & Cho Reference Gao and Cho2022). We compare it with a simpler power-law formulation in Appendix A, showing that the two are qualitatively and quantitatively similar and are expected to produce similar behaviour.

2.4. Fluid flow

Combining (2.4a,b) and (2.5), we arrive at the nonlinear flow equations:

(2.7a,b)\begin{equation} \frac{\partial{\phi_f}}{\partial{t}} +\frac{\partial}{\partial{x}}\left[{\phi_f q}-(1-\phi_f) \frac{k(\phi_f)}{\mu}\frac{\partial{p}}{\partial{x}}\right]=0 \quad\mathrm{and}\quad \frac{\partial{q}}{\partial{x}}=0, \end{equation}

where the total flux $q$ is again

(2.8)\begin{equation} q\equiv\phi_f v_f + (1-\phi_f)v_s. \end{equation}

The Darcy-scale fluid velocity and solid velocity are then given by

(2.9a,b)\begin{equation} v_f =q-\frac{1-\phi_f}{\phi_f}\,\frac{k(\phi_f)}{\mu}\frac{\partial{p}}{\partial{x}} \quad\mathrm{and}\quad v_s=q+\frac{k(\phi_f)}{\mu}\frac{\partial{p}}{\partial{x}}. \end{equation}

and the local fluid flux is

(2.10)\begin{equation} q_f=\phi_f v_f. \end{equation}

2.5. Mechanical equilibrium

The true Cauchy total stress $\boldsymbol {\sigma }$ is supported jointly by the fluid phase and the solid phase. The total stress can be decomposed into a contribution from the fluid pressure $p$ and a contribution from Terzaghi's effective stress $\boldsymbol {\sigma '}$,

(2.11)\begin{equation} \boldsymbol{\sigma}=\boldsymbol{\sigma'}-p\boldsymbol{I}, \end{equation}

where we adopt the sign convention of tension being positive. Neglecting inertia and body forces, mechanical equilibrium can be written

(2.12)\begin{equation} \boldsymbol{\nabla}\boldsymbol{\cdot}\boldsymbol{\sigma} = \boldsymbol{\nabla} \boldsymbol{\cdot}\boldsymbol{\sigma'}-\boldsymbol{\nabla}p=0. \end{equation}

In one dimension, (2.12) implies that

(2.13)\begin{equation} \frac{\partial{\sigma'}}{\partial x}=\frac{\partial{p}}{\partial{x}}, \end{equation}

where $\sigma ^\prime$ is the $xx$ component of $\boldsymbol {\sigma }^\prime$.

2.6. Elasticity law

The effective stress is the portion of the total stress that contributes to deformation of the solid skeleton. We take the solid skeleton to be elastic, with no viscous or dissipative behaviours. For confined compression in one dimension, as considered here, any elasticity law can be written in the form $\sigma ^\prime =\sigma ^\prime (\phi _f)$. Thus, (2.7a,b) can be rewritten as a nonlinear advection–diffusion equation:

(2.14a,b)\begin{equation} \frac{\partial{\phi_f}}{\partial{t}} +\frac{\partial}{\partial{x}}\left[{\phi_f q}-D_f(\phi_f)\frac{\partial{\phi_f}}{\partial{x}}\right]=0 \quad\mathrm{and}\quad \frac{\partial{q}}{\partial{x}}=0, \end{equation}

where the nonlinear composite constitutive function

(2.15)\begin{equation} D_f(\phi_f)=(1-\phi_f)\frac{k(\phi_f)}{\mu}\frac{\mathrm{d}\sigma^\prime}{\mathrm{d}\phi_f} \end{equation}

is the poroelastic diffusivity.

Hencky elasticity is a simple nonlinear hyperelasticity model that considers logarithmic strains (‘true strains’ or ‘Hencky strains’) to capture kinematic nonlinearity (Hencky Reference Hencky1931), and which is commonly used to model soft rubbers and foams (Hencky Reference Hencky1933; Anand Reference Anand1979; Xiao & Chen Reference Xiao and Chen2002) and sometimes for soft biological tissues (Marchesseau et al. Reference Marchesseau, Heimann, Chatelin, Willinger and Delingette2010; Fraldi et al. Reference Fraldi, Palumbo, Carotenuto, Cutolo, Deseri and Pugno2019). Hencky elasticity is convenient for our present purposes because (i) it reduces to linear elasticity for small strains and (ii) it uses the same two elastic parameters as linear elasticity. It is straightforward to replace Hencky elasticity in the present formulation with a different elastic behaviour, such as a neo-Hookean model, as appropriate for the problem/material of interest. Neo-Hookean elasticity is commonly used as a simple model for soft tissues (e.g. Sengers et al. Reference Sengers, Oomens and Baaijens2004; Ehlers, Karajan & Markert Reference Ehlers, Karajan and Markert2009); in the present context, we expect Hencky elasticity to provide a qualitatively similar mechanical response (see Appendix A).

For a uniaxial deformation, the relevant component of the effective stress tensor for Hencky elasticity is (MacMinn et al. Reference MacMinn, Dufresne and Wettlaufer2016; Auton & MacMinn Reference Auton and MacMinn2018)

(2.16)\begin{equation} \sigma^\prime= \mathcal{M}\frac{\ln(J)}{J}=\mathcal{M} \left(\frac{1-\phi_f}{1-\phi_{f,0}}\right) \ln\left(\frac{1-\phi_{f,0}}{1-\phi_f}\right), \end{equation}

where $\mathcal {M}$ is the $p$-wave or oedometric modulus. With appropriate initial and boundary conditions, (2.14a,b), (2.15) and (2.16) form a closed model for the evolution of the porosity.

2.7. Initial and boundary conditions

Finally, we specify appropriate initial and boundary conditions for the solid skeleton and for the fluid. As noted previously, we locate the left and right boundaries of the solid at $x=a(t)$ and $x=L$, respectively.

2.7.1. Initial conditions

Equation (2.1) suggests that $a(0)=0$. The solid is therefore initially relaxed and the initial porosity is uniform and equal to the relaxed porosity,

(2.17a,b)\begin{equation} u_s(x,0)=0 \quad\mathrm{and}\quad \phi_f(x,0)=\phi_{f,0}. \end{equation}

2.7.2. Left boundary

For $t>0$, we apply a displacement-controlled mechanical loading at the left boundary, which is therefore a moving boundary (see (2.1)). The associated boundary conditions are

(2.18a)\begin{equation} u_s(a,t)=a(t) \quad\mathrm{and}\quad v_s(a,t)=\frac{\mathrm{d}a}{\mathrm{d}t}. \end{equation}

The left boundary is also permeable, so we take

(2.18b)\begin{equation} p(a,t)=0. \end{equation}

2.7.3. Right boundary

The right boundary is impermeable and fixed in place, such that

(2.19)\begin{equation} u_s(L,t)=v_s(L,t)=v_f(L,t)=0. \end{equation}

This condition and the requirement that $q$ be uniform in space (see the end of § 2.2 and (2.7a,b)) together imply that $q\equiv 0$, meaning that there is no net flow through any cross-section. Equation (2.8) then requires that

(2.20)\begin{equation} v_f={-}\frac{(1-\phi_f)}{\phi_f} v_s, \end{equation}

meaning that the fluid and the solid always locally move in opposite directions.

2.8. Linear poroelasticity

For comparison with the fully nonlinear model, we linearise the relations given previously to arrive at linear poroelasticity, which is valid for infinitesimal deformations, $(\phi _f-\phi _{f,0})/(1-\phi _{f,0})=\partial {u_s}/\partial {x}\ll {}1$. In this limit, (2.7a,b) reduces to the linear-poroelastic diffusion equation,

(2.21)\begin{equation} \frac{\partial{\phi_f}}{\partial{t}} -\frac{\partial}{\partial{x}}\left(D_{f,0}\frac{\partial{\phi_f}}{\partial{x}}\right)\approx 0, \end{equation}

where Hencky elasticity reduces to linear elasticity,

(2.22)\begin{equation} \sigma^\prime\approx \mathcal{M}\,\frac{\partial{u_s}}{\partial{x}}=\mathcal{M}\left(\frac{\phi_f-\phi_{f,0}}{1-\phi_{f,0}}\right), \end{equation}

and $D_{f,0}=k_0\mathcal {M}/\mu$ is the constant linear-poroelastic diffusivity. With appropriate initial and boundary conditions, (2.21) is a closed linear model for the evolution of the porosity.

In the linear poroelastic model, the initial conditions and the boundary conditions for the right boundary are again (2.17a,b) and (2.19), respectively. The linearised boundary conditions for the left boundary are

(2.23a,b)\begin{equation} u_s(0,t)\approx a(t) \quad\mathrm{and}\quad v_s(0,t) \approx \frac{\mathrm{d}a}{\mathrm{d}t}, \end{equation}

where $a(t)$ is again given by (2.1). Note that these conditions are linearised relative to (2.18) by virtue of being applied at $x=0$ rather than at $x=a(t)$.

2.9. Scaling and summary

We make the above problem dimensionless via the scaling

(2.24a–h)\begin{align} \tilde{x}=\frac{x}{L},\quad \tilde{u}_s=\frac{u_s}{L},\quad \tilde{t}=\frac{t}{T_{{pe}}},\quad \tilde{\sigma}^\prime=\frac{\sigma'}{\mathcal{M}},\quad \tilde{p}=\frac{p}{\mathcal{M}},\quad \tilde{k}=\frac{k(\phi)}{k_0},\; \tilde{v}_f=\frac{v_f}{L/T_{{pe}}}, \quad \tilde{v}_s=\frac{v_s}{L/T_{{pe}}}, \end{align}

where $T_{{pe}}=L^2/D_{f,0}=\mu {}L^2/(k_0\mathcal {M})$ is the classical poroelastic timescale, which is the characteristic diffusion time for the relaxation of pressure over a distance $L$. Now taking $q\equiv {}0$, as required by the boundary conditions (see § 2.7.3), the nonlinear flow equation can be rewritten in dimensionless form as

(2.25)\begin{equation} \frac{\partial{\phi_f}}{\partial{\tilde{t}}} -\frac{\partial}{\partial{\tilde{x}}}\left[\tilde{D}_f(\phi_f)\frac{\partial{\phi_f}}{\partial{\tilde{x}}}\right]=0 \end{equation}

with nonlinear-poroelastic diffusivity

(2.26)\begin{equation} \tilde{D}_f=\frac{D_f}{D_{f,0}}=(1-\phi_f)\tilde{k}(\phi_f)\frac{\mathrm{d}\tilde{\sigma}^\prime}{\mathrm{d}\phi_f}, \end{equation}

elasticity law

(2.27)\begin{equation} \tilde{\sigma}^\prime=\left(\frac{1-\phi_f}{1-\phi_{f,0}}\right) \ln\left(\frac{1-\phi_{f,0}}{1-\phi_f}\right), \end{equation}

initial conditions

(2.28a–c)\begin{equation} \tilde{a}(0)=0, \quad \phi_f(\tilde{x},0)= \phi_{f,0}, \quad \tilde{v}_f(\tilde{x},0)=\tilde{v}_s(\tilde{x},0)=0, \end{equation}

left boundary conditions

(2.29a,b)\begin{equation} \tilde{u}_s(\tilde{a},\tilde{t})=\tilde{a}(\tilde{t})= \frac{\tilde{A}}{2} \left[1-\cos\left(\frac{2{\rm \pi}\tilde{t}}{\tilde{T}}\right)\right] ,\quad \tilde{v}_s(\tilde{a},\tilde{t})=\frac{\mathrm{d} \tilde{a}}{\mathrm{d}\tilde{t}}\quad \mathrm{and}\quad \tilde{p}(\tilde{a},\tilde{t})=0, \end{equation}

and right boundary conditions

(2.30)\begin{equation} \tilde{u}_s(1,\tilde{t}) =\tilde{v}_s(1,\tilde{t}) =\tilde{v}_f(1,\tilde{t})=0, \end{equation}

where $\tilde {A}=A/L$ and $\tilde {T}=T/T_{{pe}}$. We consider only dimensionless quantities in the following, dropping the tildes for convenience.

This one-dimensional model describes flow and mechanics in a poroelastic material subject to periodic loading. The kinematics are rigorous and nonlinear, the elasticity law is Hencky elasticity, and the permeability law is the Kozeny–Carman relation. The full and linearised problems share the same three dimensionless control parameters: the dimensionless amplitude and period of the loading, $A$ and $T$, and the relaxed porosity, $\phi _{f,0}$.

3.1. Average porosity

We begin by deriving some basic kinematic results for the average porosity. The macroscopic total volume at any instant is $1-a(t)$, whereas the total volume of solid is constant and equal to $1-\phi _{f,0}$. As a result, the total volume of fluid is $\phi _{f,0}-a(t)$ and the spatially averaged porosity is

(3.1)\begin{equation} \langle \phi_f \rangle (t) = \frac{\phi_{f,0}-a(t)}{1-a(t)}. \end{equation}

The average of $\langle \phi _f\rangle (t)$ in time over any integer number of loading cycles is then

(3.2)\begin{equation} \overline{\langle\phi_f\rangle} = \frac{1}{mT}\int_{nT}^{(n+m)T} \langle{\phi_f}\rangle \,\mathrm{d}t = 1 - \frac{1-\phi_{f,0}}{\sqrt{1-A}} \end{equation}

for any $n\geq 0$ and integer $m\geq 1$. Note that both the spatial and overall averages are smaller than $\phi _{f,0}$ for $a>0$ and $A>0$ because the loading has a non-zero mean (i.e. $\bar {a}=A/2>0$), so the material is on average compressed.

3.2. Linear poroelasticity: early time solution

The linear problem posed in § 2.8 can be rewritten as a bounded linear diffusion problem for the displacement. When the loading begins, information about the motion of the left boundary propagates into the domain via poroelastic diffusion. At early times, before this information has had time to reach the right boundary, the response is the same as if the material were semi-infinite in the $x$ direction. The corresponding semi-infinite diffusion problem involves applying the right boundary conditions at $x\to \infty$. This early time (‘$\mathrm {et}$’) solution can be derived via Laplace transform and written as a convolution integral,

(3.3)\begin{equation} u_{s,{et}}(x,t)=\frac{{\rm \pi}{}A}{T}\int_0^t\,\mathrm{erfc}\left(\frac{x}{2\sqrt{\tau}}\right)\sin\left[\frac{2{\rm \pi}(t-\tau)}{T}\right]\,\mathrm{d}\tau. \end{equation}

The corresponding porosity field can be derived from (3.3) via (2.3), and is given by

(3.4)\begin{equation} {\phi_{f,{et}}}(x,t)=\phi_{f,0}-(1-\phi_{f,0})\frac{{\rm \pi}{}A}{T}\int_0^t\,\frac{1}{\sqrt{{\rm \pi}\tau}}\exp\left(-\frac{x^2}{4\tau}\right)\sin\left[\frac{2{\rm \pi}(t-\tau)}{T}\right]\,\mathrm{d}\tau. \end{equation}

This solution is valid until the deformation spans the domain. Introducing the penetration length $\delta _{et}(t)$ as the distance from the left boundary over which the change in porosity has exceeded a threshold, the nature of the problem and the structure of the solution suggest that $\delta _{et}(t)\sim {}2\sqrt {t}$. We confirm this reasoning in Appendix B. Thus, the above solution is valid for $\delta _{et}\lesssim {}1\,\to \,t\lesssim {}1/4$.

Having originated from linear poroelasticity, the above solution is limited to small deformations ($A\ll 1$). Naïvely, this solution is valid for any loading period $T$; however, sufficiently small values of $T$ also violate the assumption of small deformations because deformation of size $\sim {}A$ are localised to a region of size $\sim \delta _{et}(t)$ at early times. As a result, we expect the maximum strain near the left boundary to be roughly of size $\max [a(t)/\delta _{et}(t)]\sim {}A/\sqrt {T}$. More precisely, (3.4) can be used to show that the extreme values of ${\phi _{f,{et}}}$ will always occur at the left boundary and that the evolution of ${\phi _{f,{et}}}$ at the left boundary is given by

(3.5)\begin{equation} {\phi_{f,{et}}}(0,t)=\phi_{f,0}-(1-\phi_{f,0})A\sqrt{\frac{\rm \pi}{T}}\, \mathcal{I}(t/T), \end{equation}

where

(3.6)\begin{equation} \mathcal{I}(y)=C(2\sqrt{y})\sin(2{\rm \pi}{}y)-S(2\sqrt{y})\cos(2{\rm \pi}{}y) \end{equation}

and $C$ and $S$ are the Fresnel cosine and sine integrals, respectively. The extreme values of ${\phi _{f,{et}}}$ are then given by

(3.7)\begin{equation} \min({\phi_{f,{et}}}) =\phi_{f,0}-(1-\phi_{f,0})A\sqrt{\frac{\rm \pi}{T}}\,\mathcal{I}_{max}, \end{equation}

and

(3.8)\begin{equation} \max({\phi_{f,{et}}}) =\phi_{f,0}-(1-\phi_{f,0})A\sqrt{\frac{\rm \pi}{T}}\,\mathcal{I}_{min}, \end{equation}

where $\mathcal {I}_{max}=\mathcal {I}(698/1909)\approx {}0.9491$ and $\mathcal {I}_{min}=\mathcal {I}(310/353)\approx {}-\!0.5406$ are the maximum and minimum values of $\mathcal {I}$, respectively. Strictly, the above solution is non-physical for parameter combinations for which the porosity decreases to 0 or increases to 1, corresponding to $\min ({\phi _{f,{et}}})=0$ and $\max ({\phi _{f,{et}}})=1$, respectively. In practice, these solutions become inaccurate for much less extreme parameter combinations as kinematic and constitutive nonlinearities become increasingly important. Hewitt et al. (Reference Hewitt, Paterson, Balmforth and Martinez2016) showed that very fast monotonic compression of a soft porous material can lead to extreme localisation near the piston in the form of a ‘bloated’ low-porosity boundary layer, the formation and evolution of which depends sensitively on the particular constitutive functions $k(\phi _f)$ and $\sigma ^\prime (\phi _f)$. We avoid these extreme loading conditions in the present study, focusing instead on scenarios that are likely to have more biological relevance. The above results confirm that the maximum local strain is indeed of size $A/\sqrt {T}$, suggesting that the above solution and the linear model in general are valid for $A\ll {}\sqrt {T}$. A similar condition can be derived by noting that, in the linear-poroelastic case, the boundary conditions at the left are applied at $x=0$ rather than at $x=a(t)$ (see (2.23a,b)). Thus, the validity of the linear-poroelastic model requires that the error in this linearisation, which is $\sim {}A$, must be negligible relative to the poroelastic diffusion length associated with the deformation, which is $\sim \sqrt {T}$.

3.3. Linear poroelasticity: full solution

The original bounded linear diffusion problem can be solved analytically via separation of variables. The resulting linear-poroelastic (‘$\mathrm {lpe}$’) displacement field is

(3.9)\begin{align} u_{s,{lpe}}(x,t) &= a(t)(1-x) \nonumber\\ &\quad - \sum_{n=1}^{\infty} 2A \sin{(n{\rm \pi} x)} \dfrac{\left[2 e^{{-}n^2{\rm \pi}^2t}-2\cos{\left(\dfrac{2{\rm \pi} t}{T}\right)}+n^2{\rm \pi} T \sin{\left(\dfrac{2{\rm \pi} t}{T}\right)}\right]}{n {\rm \pi}(4+ n^4 {\rm \pi}^2 T^2)}. \end{align}

The corresponding porosity field can be derived from (3.9) via (2.3) and is given by

(3.10)\begin{align} &\phi_{f,{lpe}}(x,t)=\phi_{f,0}- (1-\phi_{f,0})\nonumber\\ &\qquad \times \left\{ a(t) + \sum_{n=1}^{\infty} 2A \cos{(n{\rm \pi} x)} \dfrac{\left[2 e^{{-}n^2{\rm \pi}^2t}-2\cos{\left(\dfrac{2{\rm \pi} t}{T}\right)}+n^2{\rm \pi} T \sin{\left(\dfrac{2{\rm \pi} t}{T}\right)}\right]}{ (4+ n^4 {\rm \pi}^2 T^2)}\right\} . \end{align}

The corresponding fluid velocity field can be derived by taking the time derivative of (3.9) to obtain the solid velocity and then using (2.20) to arrive at

(3.11)\begin{align} &v_{f,{lpe}}(x,t)={-}\dfrac{(1-\phi_{f,0})}{\phi_{f,0}} \nonumber\\ &\qquad \times \left\{ \dot a(t)(1-x) - \sum_{n=1}^{\infty} 4 A \sin{(n{\rm \pi} x)} \dfrac{\left[ - n^2{\rm \pi} e^{{-}n^2{\rm \pi}^2t}+\dfrac{2}{T}\sin{\left(\dfrac{2{\rm \pi} t}{T}\right)}+n^2{\rm \pi} \cos{\left(\dfrac{2{\rm \pi} t}{T}\right)}\right]}{n \left(4+ n^4 {\rm \pi}^2 T^2\right)}\right\}. \end{align}

Like the early time solution, this solution provides a good approximation for $A\ll 1$ and for $A\ll {}\sqrt {T}$. Also like the early time solution, this solution provides general insight into the poromechanical response of the system to periodic loading. The expression for ${\phi _{f,{lpe}}}$ can be divided into three parts:

1. (i) a uniform quasi-static part proportional to $a(t)$, which is the linearised form of the nonlinear quasi-static solution derived in the following;

2. (ii) an early time transient that decays exponentially in time at a rate that is independent of $A$ and $T$; and

3. (iii) a periodic forced response with period $T$.

The early time transient captures the early time solution derived previously, which spans the domain after a time $t\approx {}1/4$ and then decays exponentially relative to the periodic forced response that dominates the solution thereafter. Going forward, we focus on this periodic regime. We show in Appendix C that the same reasoning also applies for large deformations.

3.4. Linear poroelasticity: response to very fast loading (Stokes’ second problem)

As noted previously, the response at early times will be confined to a region of size $\sim \sqrt {t}$, spreading diffusively until the entire domain is engaged, at which point the response will evolve exponentially toward its periodic regime. However, the oscillations in the periodic regime will be confined to a region of size $\sim \sqrt {T}$ (or $\sim {}A$ if larger, but recall that linear poroelasticity requires that $A\ll {}\sqrt {T}$). Thus, if the period is sufficiently small (i.e. $\sqrt {T}\ll 1$), the material near the piston will oscillate while the far field exists in a state of static compression. This response to very fast loading is well known from Stokes's classical ‘second problem’, in which oscillations diffuse into a semi-infinite domain with an amplitude that decays exponentially in space, much like an evanescent wave. The corresponding analytical solution to linear poroelasticity for the periodic response to very fast loading (‘$\mathrm {vf}$’) is

(3.12)\begin{equation} u_{s,{vf}}(x,t) = \frac{A}{2}\left[ 1 -x -\exp\left({-}x\sqrt\frac{\rm \pi}{T}\right)\cos\left(\frac{2{\rm \pi}{}t}{T}-x\sqrt\frac{\rm \pi}{T}\right)\right] \end{equation}

and

(3.13)\begin{align} \phi_{f,{vf}}(x,t) &= \phi_{f,0} - (1-\phi_{f,0})\frac{A}{2}\left\{ 1-\sqrt{\frac{\rm \pi}{T}}\exp\left({-}x\sqrt\frac{\rm \pi}{T}\right)\right.\nonumber\\ &\left.\quad \times\left[\cos\left(\frac{2{\rm \pi}{}t}{T}-x\sqrt\frac{\rm \pi}{T}\right)-\sin\left(\frac{2{\rm \pi}{}t}{T}-x\sqrt\frac{\rm \pi}{T}\right)\right]\right\}, \end{align}

where the first two terms in $u_{s,{vf}}$ and in $\phi _{f,{vf}}$ give the static far-field compression, which is also the (linearised) overall mean compression. This solution confirms that the oscillations will be increasingly localised near the piston as $T$ decreases, featuring near-piston oscillations with an amplitude proportional to $1/\sqrt {T}$ that decay exponentially in space over a characteristic distance $\sqrt {T}$. This solution is illustrated and discussed further in Appendix B.

3.5. Response to very slow loading (quasi-static solution)

In the full linear-poroelastic solution given previously, the porosity and displacement fields ((3.9) and (3.10)) converge to the quasi-static limits ${\phi _{f,{lpe}}}(x,t) \to \phi _{f,0}-(1-\phi _{f,0})a(t)$ and $u_{s,{lpe}}\to {}a(t)(1-x)$ as $T\to \infty$. This limit is a uniform state of strain in which poroelastic transients are fast relative to the loading period, and are therefore negligible. The fully nonlinear problem can be solved analytically in the same limit by taking $\partial {\phi _f}/\partial {t}\to {}0$. The resulting quasi-static (‘$\mathrm {qs}$’) solution is

(3.14)\begin{gather} u_{s,{qs}}(x,t)= \frac{a(t)}{1-a(t)} (1-x), \end{gather}

(3.15)\begin{gather} {\phi_{f,{qs}}}(t)= \frac{\phi_{f,0}-a(t)}{1-a(t)}, \end{gather}

and

(3.16)\begin{equation} v_{f,{qs}}(x,t)={-}\left[\frac{1-\phi_{f,0}}{\phi_{f,0}-a(t)}\right]\left\{\frac{1-x}{[1-a(t)]^2}\right\}\dot{a}(t). \end{equation}

This solution is kinematically exact for arbitrarily large values of $A$, but is only valid for $T\gg {}1$. Note that this expression for $\phi _{f,{qs}}$ is the same as that in (3.1) for $\langle {\phi _f}\rangle (t)$ because the quasi-static porosity is spatially uniform and must therefore be equal to the average porosity.

3.6. Scaling quantities

In the absence of a net flow, fluid motion is directly related to changes in porosity. Hence, we present our solutions and results below in terms of the change in porosity with respect to the overall average porosity,

(3.17)\begin{equation} \Delta{\phi_f} \equiv \phi_f - \overline{\langle{\phi_f}\rangle}. \end{equation}

This definition accounts for the non-zero mean compression.

The various results given previously suggest simple scaling relationships for the magnitudes of $\Delta {\phi _f}$ and $v_f$ in terms of the dimensionless control parameters $A$ and $T$. In particular, the magnitude of $\Delta {\phi _f}$ is captured by the spatially averaged change in porosity at mid-cycle,

(3.18)\begin{equation} \Delta\phi_f^M= \left|\langle{\Delta \phi_f}\rangle(T/2)\right| =\frac{\phi_{f,0}-A}{1-A}-\overline{\langle\phi_f\rangle}. \end{equation}

The quasi-static solution suggests that appropriate scales for the magnitude of the solid and fluid velocity are

(3.19)\begin{equation} v_s^{{\ast}}=\frac{2A}{T}, \end{equation}

and

(3.20)\begin{equation} v_f^{{\ast}}= \left( \frac{1-\phi_{f,0}}{\phi_{f,0}} \right) \frac{2A}{T}. \end{equation}

An appropriate scale for the fluid flux is therefore

(3.21)\begin{equation} q_f^{{\ast}}= \phi_{f,0} v_f^{{\ast}}. \end{equation}

4.1. Comparison with analytical solutions

We next compare the numerical solution to the linear-poroelastic and nonlinear quasi-static analytical solutions described in § 3, each of which is appropriate for a specific range of $A$ and $T$. The aim of this comparison is to quantify these ranges of validity and to examine the convergence of the numerical results to each of these special cases. To do so, we calculate all three solutions over a wide range of $A$ and $T$ and then calculate the root-mean-square (r.m.s.) relative difference between the numerical and linear-poroelastic solutions (figure 3), and between the numerical and nonlinear quasi-static solutions (figure 4). We calculate these differences using $\phi _f(a,t)$ during one cycle in the periodic regime. In both figures, we plot these differences against $T$ for several values of $A$ (left panels) and then against $A$ for several values of $T$ (right panels).

Figure 3. The r.m.s. relative difference between the full numerical solution and the linear-poroelastic analytical solution (from (3.10)) based on $\phi _f(a,t)$ during one cycle in the periodic regime. On (a), we plot the difference against $T$ for fixed values of $A$ ranging from $0.02$ to $0.2$ (dark to light). On (b), we plot the same difference against $A$ for fixed values of $T$ ranging from $0.001{\rm \pi}$ to $10{\rm \pi}$ (dark to light).

Figure 4. The r.m.s. relative difference between the full numerical solution and the nonlinear quasi-static solution (from (3.15)) based on $\phi _f(a,t)$ during one cycle in the periodic regime. On (a), we plot the difference against $T$ for fixed values of $A$ ranging from $0.02$ to $0.2$ (dark to light). On (b), we plot the same difference against $A$ for fixed values of $T$ ranging from $0.001{\rm \pi}$ to $10{\rm \pi}$ (dark to light).

As expected, figure 3 shows good agreement between the numerical and linear-poroelastic solutions for small deformations, worsening as $A$ increases. The difference scales with $A^2$ for a given value of $T$ (right panel), as expected from linear poroelasticity, which is first order in strain. The difference is insensitive to $T$ for $T\gtrsim {}1$, but scales as $T^{-1}$ for faster periods. Decreasing $T$ leads to increasingly strong localisation at the left boundary, and hence increasingly large deformations, even for small $A$, as expected from the early time and very fast analyses above (§§ 3.2–3.4). We explore this localisation in more detail in § 5.

As expected, figure (4) shows good agreement between the numerical and nonlinear quasi-static solutions for large periods, worsening as $T$ decreases. The difference scales with $T^{-1}$ for $T\gtrsim {}1$ and with $T^{-1/2}$ for shorter periods (left panel); the difference also scales with $A$ (right panel), consistent with the scaling of the non-quasi-static terms in (3.10). Based on these results, we distinguish between ‘slow loading’ (SL; $T\gtrsim {}{\rm \pi}$), where spatial variations in porosity are relatively small, and ‘fast loading’ (FL; $T\lesssim {}0.1{\rm \pi}$), where spatial variations in porosity are relatively large. For very slow loading ($T\gg {}1$), the porosity is uniform and the response is quasi-static (see § 3.5). For very fast loading ($T\ll {}1$), the oscillations are localised near the left boundary and the right portion of the material is in a state of static compression (see § 3.4).

5.1. Impact of loading period

To visualise the distinct poromechanical responses for FL and SL, and the transition between them, we plot $u_s$, $\Delta {\phi _f}$ and $q_f$ (all normalised) over one cycle in the periodic regime for $\phi _{f,0}=0.75$, $A=0.1$, and four different values of $T$: two for FL (left two columns) and two for SL (right two columns) (figure 5).

Figure 5. Evolution of normalised $u_s$ (a–d) and $\Delta {\phi _f}$ (e–h) at times $t=nT$ to $(n+1)T$ in increments of $0.1T$ (dark to light), where $n$ is an integer, during one cycle in the periodic regime for $A=0.1$, $\phi _{f,0}=0.75$, and $T=0.03{\rm \pi}$ (a,e,i,m), $0.1{\rm \pi}$ (b, f, j,n), ${\rm \pi}$ (c,g,k,o) and $10{\rm \pi}$ (d,h,l,p). As in figure 2, we plot all fields in (a–h) against the Lagrangian spatial coordinate $X=x-u_s(x,t)$ for clarity. In (i–l), we plot normalised $q_f$ against $t$ for 10 different values of $X$ from $0$ to $1$ (dark to light). We distinguish between the loading half of the cycle ($\dot {a}>0$; solid curves) and the unloading half of the cycle ($\dot {a}<0$; dashed curves). In (m–p), we plot phase portraits of $\sigma ^\prime (a,t)$ against $a(t)$ from $t=0$ to $t=20T$, emphasising the last cycle (dotted black).

Figure 5(a–d) show that, for FL (i.e. $T=0.03{\rm \pi}$ and $0.1{\rm \pi}$), the displacement is highly nonlinear in $X$ (and in $x$), with substantial differences between the loading and unloading phases. As $T$ increases, transitioning into SL (i.e. $T=1{\rm \pi}$ and $10{\rm \pi}$), the displacement is increasingly linear in $X$ and converges toward the quasi-static solution, which is fully determined by the instantaneous value of $a$ and is thus symmetric between loading and unloading. For all values of $T$, the displacement is of characteristic size $A$ at the left and vanishes at the right.

For SL, $\Delta {\phi _f}$ is uniform in space and varies only in time, per the quasi-static solution. The material is in a uniform state of compression, fully determined by the instantaneous value of $a$ and thus symmetric between loading and unloading and fully relaxed at the beginning/end of each cycle. As $T$ decreases, transitioning into FL, $\Delta {\phi _f}$ becomes increasingly localised near the left boundary and also increasingly asymmetric between loading and unloading (figure 5e). For FL, the material experiences a substantial amount of tension in the left portion of the domain during unloading, despite the overall mean compression. Tension emerges for FL because this regime is, by definition, one where the rate of loading is much faster than the poroelastic relaxation time, so the left boundary must pull the material to the left during unloading. The right portion of the domain experiences an overall more limited range of porosities and remains compressed throughout the cycle, never reaching a state of tension or even full unloading. The latter feature is also visible in the corresponding displacement fields.

The fluid flux $q_f=\phi _fv_f$ is particularly relevant to the transport of solutes since it drives advection. Fluid leaves the domain during the loading phase of the cycle ($q_f(x=a,t)=q_f(X=0,t)<0$ when $\dot {a}>0$) and enters the domain during unloading. For SL, $q_f$ is entirely in phase with the loading, but opposite in sign. For FL, the peak value of $q_f$ in the interior exhibits a lag relative to the peak value of $\dot {a}$, and this lag increases with $x$. Note that the fluid flux is orders of magnitude larger for FL than for SL because the rate of loading is orders of magnitude faster, but this variation is largely scaled out by normalisation with $q_f^\ast (\phi _{f,0},A,T)\propto {}T^{-1}$, per (3.21) (figure 5i–l).

The extreme values of $q_f$ at the left boundary are larger in loading than in unloading, but progressively more symmetric as $T$ increases. This asymmetry is a result of the kinematic and constitutive nonlinearity of large deformations. We show in Appendix E that this asymmetry originates in the nonlinear kinematics of large deformations, meaning that it emerges from the nonlinear model during large deformations even with linear elasticity and constant permeability. This asymmetry is then strongly amplified by deformation-dependent permeability (relative to constant permeability) and slightly suppressed by Hencky elasticity (relative to linear elasticity). The latter occurs because Hencky elasticity stiffens in compression, resulting in a larger poroelastic diffusivity and therefore weaker localisation during loading (see also Appendix A).

The asymmetry between loading and unloading is also highlighted by the evolution of $\sigma ^{\prime }(a,t)$ (figure 5m–p). For FL, $\sigma ^\prime (a,t)$ exhibits hysteresis: for the same value of $a(t)$, the value of $\sigma ^\prime (a,t)$ is considerably higher in magnitude during loading than during unloading (and tensile during much of the unloading phase). This hysteresis decreases as $T$ increases, such that, for SL, $\sigma ^\prime (a,t)$ is modest in magnitude, fully compressive and symmetric between loading and unloading (i.e. non-hysteretic); as a result of the latter, the phase portrait for SL is a single curve (rather than a loop). These phase portraits also illustrate the convergence to the periodic regime: in all cases, the overall response grows less extreme as the transient component decays exponentially over the first $\sim {}T^{-1}$ cycles.

Since a key difference between SL and FL is that the deformation is uniformly distributed in the former and localised toward the left in the latter, we study the emergence of this localisation to further quantify the transition from SL to FL. In figure 6, we plot the normalised maximum and minimum values of $\Delta {\phi _f}$, $q_f$ and $\sigma ^\prime$ at the left boundary ($X=0$) and then either at the right boundary ($X=1$) for $\Delta {\phi _f}$ and $\sigma ^\prime$ or at the material mid-point ($X=1/2$) for $q_f$. The latter is necessary because $q_f$ vanishes at $X=1$. We consider a range of periods $T$ at two fixed values of $A$ for $\phi _{f,0}=0.75$.

Figure 6. We illustrate the transition from FL to SL as $T$ increases by plotting the normalised maximum and minimum values $\Delta {\phi _f}$ (a), $q_f$ (b) and $\sigma ^\prime$ (c) against $T$ for $A=0.1$ (solid lines) and $A=0.001$ (dotted lines) during the periodic regime. We show the maximum values (dark colours) and minimum values (light colours) of all three quantities at the left boundary ($X=0$, blue curves), of $\Delta {\phi _f}$ and $\sigma ^\prime$ at the right boundary ($X=1$, red curves), and of $q_f$ at the material midpoint ($X=1/2$, red curves).

Figure 6 shows that, for SL, the deformation is uniform, varying only in time. The maxima and minima of both $\Delta {\phi _f}$ and $\sigma ^{\prime }$ remain separate, but their left and right values converge. The flux $q_f$ remains linear in space, even for SL (see § 3.5). As $T$ decreases, the deformation is progressively localised at the left, such that the right is increasingly static. At the right, the maximum and minimum values of $\Delta {\phi _f}$ and $q_f$ converge to zero while the maximum and minimum values of $\sigma ^\prime$ converge to weak compression, as also illustrated in figure 5.

We find that a small-amplitude deformation ($A=0.001$) and a large-amplitude deformation ($A=0.1$) exhibit qualitatively similar features. However, large deformations lead to an increasingly strong asymmetry between the maxima and minima of all quantities at the left boundary as $T$ decreases. In particular, the curves are biased downward, such with higher magnitudes reached in loading (minima) than in unloading (maxima). This asymmetry is does not occur for small deformations, for which the maxima and minima of all quantities are symmetric in magnitude. Note also that the maxima of $\Delta {\phi _f}$ and $\sigma ^{\prime }$ increase as $T$ decreases for both values of $A$, whereas the maximum of $q_f$ is independent of $T$ for small deformations but decreases with $T$ for large deformations, consistent with figure 5.

5.2. Impact of loading amplitude

We next examine the impact of loading amplitude $A$ on the poromechanical response. We first assess the interaction between amplitude and period by plotting $\Delta {\phi _f}$ against $X$ at mid-cycle for two different values of $T$ (one for SL and one for FL) and for five different values of $A$, ranging from small to large deformations ($A=0.002$ to $0.2$; figure 7).

Figure 7. Profiles of $\Delta {\phi _f}$ at mid-cycle during the periodic regime for $\phi _{f,0}=0.75$ and for five values of $A$ ranging from small to large deformations (dark to light), each for two different periods: $T=0.1{\rm \pi}$ (FL; solid) and $T=10{\rm \pi}$ (SL; dashed). Panels (a) and (b) are non-normalised and normalised, respectively, showing that this normalisation scales out the leading-order impact of $A$ on $\Delta {\phi _f}$.

Figure 7(a) shows that, for SL, the non-normalised value of $\Delta {\phi _f}$ is uniform in $X$ and increases with $A$, as expected. For FL, the non-normalised value of $\Delta {\phi _f}$ is increasingly non-uniform with $X$ as $A$ increases. Relative to the SL case, the deformation is increasingly amplified at the left boundary and suppressed at the right boundary (recall that the mean value of $\Delta {\phi _f}$ is independent of $T$). The right panel shows the same results, but now normalised by $\Delta \phi _f^M$ (3.18). By definition, all of the SL curves collapse onto $\Delta {\phi _f}/\Delta \phi _f^M=-1$ ($\Delta \phi _f^M$ is the magnitude of the quasi-static change in porosity at mid-cycle). The FL curves also nearly collapse onto a master curve, suggesting that this normalisation scales out the leading-order impact of $A$, even for large deformations for FL. The largest deviations from this collapse are near the left boundary, where nonlinearities are particularly pronounced.

In figure 8, we focus on FL by fixing $T=0.1{\rm \pi}$ and plotting the evolution of $\Delta {\phi _f}$, $q_f$ and $\sigma ^\prime (a,t)$ for three different values of $A$, ranging from small to large deformations.

Figure 8. Evolution of normalised $\Delta {\phi _f}$, $q_f$ and $\sigma ^\prime (a,t)$ for $\phi _{f,0}=0.75$, $T=0.1{\rm \pi}$, and $A=0.01$ (a,d,g), $0.1$ (b,e,h) and $0.2$ (c, f,i). In (a–c), we plot $\Delta {\phi _f}$ against the Lagrangian spatial coordinate $X=x-u_{s}(x,t)$ at times $t=nT$ to $(n+1)T$ in increments of $0.1T$ (dark to light) during one cycle in the periodic regime. In (d–f), we plot $q_f$ against $t$ for ten different values of $X$ from $0$ to $1$ (dark to light); we distinguish between loading ($\dot {a}>0$; solid curves) and unloading ($\dot {a}<0$; dashed curves). In (g–i), we plot phase portraits of $\sigma ^\prime (a,t)$ against $a(t)$ from $t=0$ to $t=20T$, emphasising the last cycle (dotted black lines).

We find that increasing $A$ amplifies the asymmetry in the extreme values of $\Delta {\phi _f}$ and $q_f$ at the left boundary in loading and unloading, as noted in figure 6. The normalisation scales out the leading-order impact of $A$ on all of these quantities, as noted in figure 7; the non-normalised values of all three quantities would vary by two orders of magnitude across this range of $A$. As $A$ increases, $\sigma ^\prime (a,t)$ exhibits more hysteresis and larger normalised magnitudes, in accordance with the amplified asymmetry.

We further explore the impact of $A$ in figure 9 by plotting the normalised maxima and minima of $\Delta {\phi _f}$, $q_f$ and $\sigma ^\prime$ on the left and on the right, as in figure 6, but now against $A$ for two different values of $T$, showing the smooth transition from small deformations ($A \lesssim 0.01$) to large deformations ($A \gtrsim 0.01$). For small deformations, the normalised maxima and minima of all quantities at the left and at the right become independent of $A$, and the extreme values of $\Delta {\phi _f}$ and $q_f$ become symmetric. For SL, $\Delta {\phi _f}$ and $\sigma ^\prime$ become uniform in space, such that their normalised left and right values are equal and depend only weakly on $A$ for the largest deformations shown here (e.g. $A\gtrsim {}0.1$). For FL, the values of all three quantities at the left become increasingly asymmetric and biased downward as $A$ increases.

Figure 9. We illustrate the transition from small to large deformations as $A$ increases by plotting the normalised maximum and minimum values $\Delta {\phi _f}$ (a), $q_f$ (b) and $\sigma ^\prime$ (c) against $A$ for $T=0.1{\rm \pi}$ (solid lines) and $T=10{\rm \pi}$ (dotted lines) during the periodic regime. We show the maximum values (dark colours) and minimum values (light colours) of all three quantities at the left boundary ($X=0$, blue curves), of $\Delta {\phi _f}$ and $\sigma ^\prime$ at the right boundary ($X=1$, red curves) and of $q_f$ at the material midpoint ($X=1/2$, red curves).

5.3. Impact of initial porosity

Finally, we consider the initial porosity $\phi _{f,0}$. In figure 10, we plot the distribution of $\Delta {\phi _f}$ at mid-cycle for three different values of $\phi _{f,0}$ for fixed $A=0.01$ and for two different values of $T$, one for FL ($T=0.1{\rm \pi}$) and SL ($T=10{\rm \pi}$). For both SL and FL, the non-normalised values of $\Delta {\phi _f}$ (left panel) increase in magnitude as $\phi _{f,0}$ decreases because the same change in total volume is a larger portion of the initial fluid volume. The right panel shows that normalisation by $\Delta {\phi _f^M}$ scales out the leading-order impact of varying $\phi _{f,0}$ on $\Delta {\phi _f}$, again with small deviations at the left boundary for FL, when nonlinearities are most pronounced.

Figure 10. Impact of $\phi _{f,0}$ on $\Delta {\phi _f}$ at mid-cycle in the periodic regime for $A=0.1$ and for two different periods: $T=0.1{\rm \pi}$ (FL; solid) and $10{\rm \pi}$ (SL; dashed). We show results for three values of $\phi _{f,0}$ (increasing from dark to light). Panels (a) and (b) are non-normalised and normalised, respectively, showing that this normalisation scales out the leading-order impact of $\phi_{f,0}$ on $\delta{\phi_f}$.

We next plot the evolution of $\Delta {\phi _f}$ and $\sigma ^\prime$ in time for $A=0.01$ and $T=0.1{\rm \pi}$ for the same three values of $\phi _{f,0}$ (figure 11), confirming that normalisation scales out the primary impact of $\phi _{f,0}$ on $\Delta {\phi _f}$ and on $\sigma ^\prime (a,t)$.

Figure 11. Impact of $\phi _{f,0}$ on the evolution of $\Delta {\phi _f}$ and $\sigma ^\prime (a,t)$ for $A=0.01$, $T=0.1{\rm \pi}$, and $\phi _{f,0}=0.25$ (a,d), $0.5$ (b,e) and $0.75$ (c, f). In (a–c), we plot $\Delta {\phi _f}$ against the Lagrangian spatial coordinate $X=x-u_s$ at times $t=nT$ to $(n+1)T$ in increments of $0.1T$ (dark to light) during one cycle in the periodic regime, and we distinguish between loading ($\dot {a}>0$; solid curves) and unloading ($\dot {a}<0$; dashed curves). In (d–f), we plot phase portraits of $\sigma ^\prime (a,t)$ against $a(t)$ from $t=0$ to $t=20T$, emphasising the last cycle (dotted black lines).

5.4. FL in biological examples

We conclude by considering appropriate parameter ranges for several examples of periodic loading in soft biological tissues (table 1). Based on these values, we calculate the poroelastic timescale $T_{{pe}}$ and the relative dimensionless loading period $\tilde {T}=T/T_{{pe}}$ for each example to understand the ranges of relevance for $\tilde A$, $\tilde T$ and $\phi _{f,0}$.

Table 1. Material and loading parameters for some examples of biological materials.

Table 1 indicates that, for a variety of soft tissues, deformations are in the range of the dimensionless amplitudes $\tilde {A}$ considered here, with nearly all being near the upper end. The range of $\phi _{f,0}$ is slightly wider than the range considered here, but we showed in § 5.3 that this value has little impact on the normalised mechanical response of the material. The range of poroelastic timescales $T_{{pe}}$ and respective loading frequencies suggest that dimensionless loading periods $\tilde {T}$ span a wide range, with many corresponding to FL. These values justify our analysis and underscore the importance of characterising the nonlinearity of large poromechanical deformations during FL in particular.