2.2.1. Elasticity

For a KGD or radial fracture with radius $R$ in an infinite homogeneous linear elastic solid, the relationship between the fracture aperture ${w}$, which is in elastic equilibrium with the imposed pressure ${p}$, can be represented (Hills et al. Reference Hills, Kelly, Dai and Korsunsky1996; Dontsov Reference Dontsov2016, Reference Dontsov2017) by an integral equation of the form

(2.1)\begin{equation} p(s,t)={-}\frac{E^{\prime }}{2{\rm \pi} R}\int_{0}^{1} M(s,s^{\prime })\,\frac{\partial w}{\partial s^{\prime }}\, {\rm d} s^{\prime}, \end{equation}

where for the symmetric KGD fracture, the kernel $M(s,s^{\prime })$ is

(2.2)\begin{equation} M(\rho ,s)=\frac{s}{s^{2}-\rho^{2}}, \end{equation}

while for the radially symmetric fracture, $M(s,s^{\prime })$ is given by

(2.3)\begin{equation} M(\rho ,s)=\left\{ \begin{array}{@{}ll} \displaystyle \dfrac{1}{\rho }\,K \left(\dfrac{s^{2}}{\rho^{2}}\right)+\dfrac{\rho }{s^{2}-\rho^{2}}\,E \left(\dfrac{s^{2}}{\rho^{2}}\right), & \rho >s ,\\ \displaystyle \dfrac{s}{s^{2}-\rho^{2}}\,E \left(\dfrac{s^{2}}{\rho^{2}}\right), & \rho < s , \end{array} \right. \end{equation}

where $K(\cdot )$ and $E(\cdot )$ are complete elliptic integrals of the first and second kind, respectively.

We observe that the expression in (2.2) has been obtained from the following integral equation with a Cauchy kernel by exploiting symmetry:

(2.4)\begin{equation} p(s,t)={-}\frac{E^{\prime }}{4{\rm \pi} R}\int_{{-}1}^{1}\, \frac{\partial w}{\partial s^{\prime }}\frac{1}{s^{\prime }-s}\, {\rm d} s^{\prime }={-}\frac{E^{\prime }}{ 4{\rm \pi} R} \int_{0}^{2}\frac{\partial \hat{w}}{\partial \hat{s}^{\prime }}\, \frac{1}{\hat{s}^{\prime }-\hat{s}}\, {\rm d} \hat{s}^{\prime}.\end{equation}

Shifting to the tip coordinate $\hat {s}$ in (2.4), the first integral reduces to the second integral.

Singularity structure of the kernel functions:

If we consider source points $s=\rho \pm \varepsilon$ a small distance $\varepsilon$ either side of the receiving point $\rho$, and expand both kernel functions $M$ given in (2.2) and (2.3), we then observe that the dominant behaviour of both kernel functions reduce to that of the Cauchy kernel given in (2.4):

(2.5)\begin{equation} M(\rho ,\rho \pm \varepsilon )={\pm} \frac{1}{2\varepsilon }+ \cdots \sim \frac{1}{2(s-\rho )},\end{equation}

where for the KGD case the next term in the expansion is $1/4$, and for the radial case the next term in the expansion is $O(\log (\varepsilon ))$ – both of which, upon integration over a finite interval, yield a result that is finite. Thus the dominant behaviour for the pressure field can be determined by considering the Cauchy kernels given in (2.4) and (2.5).

Given the dominant behaviour of the Cauchy kernel, we will make use of the following integral in the analysis of the pressure associated with an aperture that is a power law:

(2.6)\begin{equation} \int_{0}^{a}\frac{\hat{s}^{\kappa }}{\hat{s}-\hat{\rho}}\, {\rm d} \hat{s}=f(\hat{\rho},\kappa )\,\hat{\rho}^{\kappa }+a^{\kappa }\sum_{m=0}^{u_{\kappa }} \frac{1}{\kappa -m}\left( \frac{\hat{\rho} }{a} \right)^{m},\end{equation}

where $f(\hat {\rho },\kappa )=-{\rm \pi} \cot ({\rm \pi} \kappa )$ and $u_{\kappa }=\infty$ for $-1<\kappa \notin \mathbb {Z}^{+}=\{ 0,1,\ldots \}$, while $f(\hat {\rho },n )=\log | ({\hat {\rho } -a})/{\hat {\rho } }|$ and $u_{n}=n-1$ for $\kappa =n\in \mathbb {Z}^{+}$. Typically, asymptotic analysis for propagating hydraulic fractures considers a semi-infinite hydraulic fracture in a state of plane strain moving at a constant velocity (Garagash et al. Reference Garagash, Detournay and Adachi2011; Dontsov & Peirce Reference Dontsov and Peirce2015). The relevant integral in that case can be obtained by letting $a\rightarrow \infty$ in (2.6), which yields the important result that a power law with $-1<\kappa <0$ is an eigenfunction of the integral operator on the left-hand side of (2.6), since the summation on the right-hand side of (2.6) vanishes. In § 3, we will see that power laws with $-1<\kappa <0$ fail to yield a dominant balance for a receding hydraulic fracture, so the eigenfunction result needs to be extended to include values of $\kappa \ge 0$. For the integral to converge for $\kappa \ge 0$, it is necessary to restrict to $a<\infty$ and to establish the dominant behaviour of the integral operator in (2.6).

2.2.2. Lubrication

By combining Poiseuille's law with the continuity equation, we obtain the lubrication equation (Savitski & Detournay Reference Savitski and Detournay2002; Dontsov Reference Dontsov2016, Reference Dontsov2017) relating $w(s,t)$ and $p(s,t)$, which, expressed in terms of the stretched coordinate $s$, is

(2.7)\begin{equation} \frac{\partial w}{\partial t}-s\,\frac{\dot{R}}{R}\,\frac{\partial w}{\partial s}=\frac{1}{\mu^{\prime}R^{2}s^{\delta -1}}\,\frac{\partial }{\partial s} \left( s^{\delta -1}w^{3}\,\frac{\partial p}{\partial s} \right) -g(s,t), \end{equation}

where $\delta$ is the dimension parameter defined above, which is 1 for the KGD case, and 2 for the radial case. Fluid loss to the permeable rock is assumed to follow the inverse square root leak-off model

(2.8)\begin{equation} {g}(s,t)=\frac{C'}{\sqrt{t-t_{0}({s},t)}},\end{equation}

in which $t_{0}({s},t)$ denotes the time of first exposure of point $s$ to the fracturing fluid. After arrest, $\dot {R}\le 0$, and as a result previously unfractured rock is no longer added to the hydraulic fracture so that the leak-off term ${g}$ is no longer singular. Indeed, beyond the arrest time, ${g}$ becomes progressively more spatially uniform. Thus in order to perform a local analysis to determine the asymptotic behaviour near the tip after arrest, to first order, ${g}$ will be assumed to be spatially homogeneous in the tip region, and denoted by $\hat {g}_0(t)$.

2.2.3. Initial and boundary conditions for the receding fracture

Initial conditions: At the onset of recession $t=t_r$, the initial conditions are

(2.9a,b)\begin{equation} R(t_r)=R_r,\quad w(s,t_r) = w_r(s), \quad \text{where } s=r/R_r. \end{equation}

Boundary conditions: For coalescent fluid and fracture fronts, the boundary conditions at the crack tip $s=1$ are given by zero fracture opening and zero flux conditions (see Detournay & Peirce Reference Detournay and Peirce2014)

(2.10a,b)\begin{equation} w(1,t)=0,\quad w^{3}\,\boldsymbol{\nabla} p|_{s=1}=0. \end{equation}

4.5.1. Closed-form expression for the sunset solution

The similarity solution corresponding to $\gamma = 1/2$, when expressed as a function of $s$, is given by

(4.17)\begin{equation} w(s,t^{\prime})=g_{0}t^{\prime}( 1-s^{2}),\quad s =r/R , \quad R =\varLambda t^{\prime {1/2}}. \end{equation}

The pressure $p$ can be determined explicitly for the KGD case by substituting (4.17) into (2.1) and using the kernel given in (2.2):

(4.18)\begin{equation} p_{KGD}(s,t)=\frac{t^{\prime}E^{\prime }g_{0}}{2{\rm \pi} R} \left[ 2+s\log \left| \frac{ 1-s}{1+s}\right| \right]. \end{equation}

Similarly, the pressure $p$ for the radial case can be determined numerically by substituting (4.17) into (2.1) and using the kernel (2.3). For the radial sunset solution, the pressure at the wellbore $s=0$ can be shown to be

(4.19)\begin{equation} p_{rad}(0,t)=\frac{t^{\prime} E^{\prime }g_{0}}{2R}. \end{equation}

Since $\gamma = 1/2$ for the sunset solution, it follows from (4.18), (4.19) or (4.10) that the time exponent for the pressure is $\beta =1-\gamma =1/2$. In order to close the solution to the forward problem for a receding hydraulic fracture with initial conditions (2.9a,b), the only remaining parameter to specify is the prefactor $\varLambda$ in the power law for $R$ defined in (4.2c).

In figure 2, we provide IMMA solution results (indicated by solid black curves) for a receding radial hydraulic fracture on the interval $t \in [t_r,t_c)$, i.e. from the time of initiation of recession $t_r$ to the time of closure $t_c$. Because of its practical importance, and since the results for the KGD case are essentially similar, we choose to provide these results for only the radial case. The solution depicted in figure 2 corresponds to the dimensionless parameters $(\phi ^V,\omega )=(2,10^{-6})$ defined in § 4.6. In figure 2(a), we plot the scaled fracture radius $R/R(t_r)$ as a function of the scaled reverse time $t^{\prime }/t_r$ on a log–log scale. The dashed red line represents the log–linear regression using data points sampled near $t\sim t_c$, whose gradient is $0.49$, in close agreement with the exponent $\gamma = 1/2$ for the fracture radius of the sunset solution given in (4.17). In figure 2(b), we plot the scaled wellbore aperture $w(0,t)/w(0,t_r)$ (referenced to the left-hand vertical axis) plotted as a function of the scaled reverse time $t^{\prime }/t_r$. The dashed red line represents a linear regression using data points sampled near $t\sim t_c$. The dash-dotted blue curve (referenced to the right-hand vertical axis) represents the decaying leak-off term $g$. The horizontal dashed blue line represents the estimate of $g_0$ obtained from the gradient of the dashed red line.

Figure 2. (a) Plot of the scaled fracture radius $R/R(t_r)$ (solid black) as a function of the scaled reverse time $t^{\prime }/t_r$. (b) The scaled wellbore aperture $w(0,t)/w(0,t_r)$ (solid black) referenced to the left-hand vertical axis plotted as a function of the scaled reverse time $t^{\prime }/t_r$. The decaying leak-off term g (dash-dotted blue) referenced to the right-hand vertical axis is plotted as a function of $t'/t_r$.

4.5.2. Estimating the time to closure and closing the forward problem

For the receding fracture, $\varLambda$ will be determined by the amount of fluid in the fracture at a given radius $R$ on the way to closure. Indeed, using the initial conditions (2.9a,b), assuming that $g$ is constant over the period of recession, and using (4.4) and (4.6a,b), it is possible to obtain the following estimate of the time from the initiation of recession to closure:

(4.20)\begin{equation} t_{c}-t_r\sim \frac{\delta (2+\delta )}{2g_{0}} \int_{0}^{1}w_{r}(s)\,s^{ \delta -1}\, {\rm d} s.\end{equation}

To illustrate this estimate, consider the example shown in figure 2 and use the average value over the interval $[t_r,t_c)$ of $g\sim 50$. We obtain the estimate $(t_c-t_r)\sim 7\times 10^{-4}$, which is close to that obtained from the numerical solution $(t_c-t_r)\sim 6.13 \times 10^{-4}$. Now $\varLambda$ can be determined from the initial radius $R_r$ and the expression for $R$ in (4.2c).

4.5.3. Estimating the leak-off coefficient $C^{\prime }$

We have seen that the solution of the forward problem comprising the coupled equations (2.1) and (2.7) subject to the initial (2.9a,b) and boundary (2.10a,b) conditions reduces, in the asymptotic limit $t\rightarrow t_c$, to the sunset solution (4.17). However, unlike the initial conditions that define a propagating fracture, which are that $w\equiv 0$ and $R\equiv 0$ at $t=0$, the initial conditions (2.9a,b) for a receding fracture are not known a priori. Indeed, the only feasible way to establish reasonable initial conditions for $R_r$ and $w_r$ is via a numerical solution. This involves modelling a hydraulic fracture that has completed the propagation and deflation at rest phases, and is at the point of transition to recession, in order to establish the initial conditions $R_r$ and $w_r$. In addition, the history embedded in the trigger time function $t_0(s,t)$ needs to be determined in order to be able to evaluate the leak-off function $g$ during recession. From a practical point of view, none of this information can be inferred from quantities that can be measured, so it would seem that there is little practical utility in the sunset solution when posed as a forward initial–boundary value problem since it relies for its initial conditions on a numerical algorithm that, in itself, is quite adequate for determining the receding fracture solution. Thus a preliminary assessment would suggest that the primary significance of the sunset solution is that it establishes the theoretical result that this class of receding fractures reduces to a self-similar form as the fracture approaches closure.

However, the similarity solution is, by definition, an asymptotic solution determined in the limit $t^{\prime } \rightarrow 0$, which is devoid of history dependence and is also dependent not on the detailed functional form of the aperture but rather on the radius and initial fracture volume at the start of recession, as can be seen from the estimate (4.20). It is also significant that the sunset solution emerges only because the flux term (4.11) becomes subdominant to the other terms in (4.1) in the limit $t^{\prime } \rightarrow 0$, which essentially reduces (4.1) to a local kinematic condition relating the receding aperture evolution to the leak-off function $g$. It can be seen from (4.11) that this separation of dynamics from kinematics in the decoupling of elasticity from fluid flow removes the dependence of the sunset solution on two of the four fundamental material parameters defined in § 2.1, namely, the dynamic viscosity $\mu ^{\prime }$ and plane strain modulus $E^{\prime }$. We also observe that because the receding fracture does not break new rock, the solution does not depend directly on the fracture toughness $K'$. It is this decoupling of kinematics from dynamics that is responsible for the particularly simple form of the sunset solution (4.17) and its dependence on the single material parameter $C^{\prime }$ through $g_0$.

Potentially, the isolation of the single fundamental parameter $C^{\prime }$ by the novel sunset solution could lead to a new procedure to estimate $C^{\prime }$ from decaying aperture measurements, which has hitherto not been explored in field or laboratory experiments. We reason as follows. First, we note that since the fracture is shrinking to its initial footprint as it closes, it follows that $t_{0}({s},t)\overset {t\rightarrow t_{c}}{\rightarrow }0$ for all points within the receding fracture footprint. Thus from (2.8), we obtain the following approximation for $g_{0}$:

(4.21)\begin{equation} {g}(s,t)\sim g_{0}(t)\overset{t\rightarrow t_{c}}{\sim } \frac{C^{\prime }}{\sqrt{t}}.\end{equation}

We note that $t$ is available readily since it is simply the time elapsed since the initiation of the fracture. As an illustration of this procedure, consider the numerical solution for the decaying fracture aperture represented by the solid black curve in figure 2(b) to be a proxy for the closing fracture aperture measurement from the field or a laboratory experiment. From the sunset solution (4.17), we see that the gradient of the dashed red line in figure 2(b) provides an estimate of the leak-off velocity close to the closure time, $g_0 \overset {t\rightarrow t_{c}}{\approx }34$, whose value is illustrated by the horizontal dashed blue line referenced to the right-hand axis in that figure. The closure time in this case is $\sqrt {t_c} \sim 2.9 \times 10^{-2}$, which along with (4.21) enables us to obtain the estimate

(4.22)\begin{equation} C^{\prime }\approx g_{0}(t_{c})\,\sqrt{t_{c}}=0.986.\end{equation}

Comparing the estimate provided in (4.22) with the value of $C^{\prime }=1$ used to generate the numerical solution, we see that the error in the estimate of $C^{\prime }$ provided by the sunset solution is less than $2\,\%$.

4.5.4. Asymptotic estimate for the radius $R(t')$

If the wellbore fluid pressure is also monitored (and it is assumed that the far-field stress $\sigma _0$ was known from the closure pressure) and one knows $g_0(t_c)$ from the aperture closure measurements, then (4.18) and (4.19) also provide the following a posteriori asymptotic estimate for the decaying fracture length/radius $R$:

(4.23)\begin{equation} R(t^{\prime })\sim \frac{t^{\prime }E^{\prime }\,g_{0}(t_{c})}{2^{\delta -1} {\rm \pi}^{2-\delta}\,p(0,t^{\prime })}. \end{equation}

4.5.5. Order of magnitude estimate for the permeability $k$

As mentioned in § 2.1, the leak-off coefficient $C'$ has a dual interpretation. For a cake-building fracturing fluid, $C'$ is twice the classical Carter's leak-off coefficient $C_{L}$, but for Newtonian fluids (such as water) that do not plug the pores of the rock, $C'$ can be interpreted in terms of the rock permeability. This interpretation requires that the diffusion length – the thickness of the region adjacent to the fracture wall where the pore pressure is perturbed from its initial or far-field value $p_{0}$ – is small compared to the fracture dimension. Indeed, the early-time solution of the diffusion equation indicates that $C'$ can then be expressed as

(4.24)\begin{equation} C'=\frac{2k(\sigma_{0}-p_{0})}{\mu\sqrt{{\rm \pi} c}},\end{equation}

since the net pressure $p_{f}-\sigma _{0}$ is small compared to $\sigma _{0}-p_{0}$, so that $p_{f}-p_{0}\approx \sigma _{0}-p_{0}$. In the above expression, $k$ is the intrinsic permeability of the rock, $\mu$ is the viscosity of the fracturing fluid assumed to be the same as the viscosity of the pore fluid, and $c$ is the diffusivity coefficient. Noting that $c$ can be written as $k/\mu S$, where $S$ is the storage coefficient, which can be approximated as $\phi /K_{f}$ (with $\phi$ denoting the porosity, and $K_{f}$ the bulk modulus of the fluid), the permeability $k$ can be expressed as

(4.25)\begin{equation} k=\frac{{\rm \pi} C'^{2}\mu K_{f}}{4\phi(\sigma_{0}-p_{0})^{2}}.\end{equation}

This expression provides a means to estimate the permeability $k$ from $C'$, if the other parameters in (4.25) can be estimated. It must be remembered that the permeability varies by many orders in rocks, and that a realistic expectation is that this procedure will result in only an order of magnitude estimate for $k$.

4.6.1. Scaling before shut-in

In this subsubsection, we consider the scalings associated with KGD and radial fractures driven to propagate by the injection of a viscous fluid at a constant rate $Q_{0}$ having dimensions $[ {L^{1+\delta }}/{T}]$. Following Lister (Reference Lister1990), we define three pressures fundamental to the hydraulic fracture process, namely, the pressure drop due to viscous flow in the crack $p_{m}$, the pressure required to open the crack in the elastic medium $p_{e}$, and the crack extension pressure $p_{k}$, whose magnitudes scale as follows:

(4.26a–c)\begin{equation} p_{m}\sim \frac{\mu^{\prime }R^{2-\delta }Q_{0}}{w^{3}},\quad p_{e}\sim \frac{E^{\prime }w}{R}\quad \text{and} \quad p_{k}\sim \frac{K^{\prime }}{R^{1/2}}. \end{equation}

There are also three fundamental volumes that need to be accounted for in the injection and leak-off process, namely, the injected volume $V_{i}$, the volume of fluid contained in the fracture $V_{f}$, and the volume of fluid that has leaked off $V_{c}$, whose magnitudes scale as follows:

(4.27a–c)\begin{equation} V_{i}\sim Q_{0}t,\quad V_{f}\sim wR^{\delta }\quad \text{and} \quad V_{c}\sim C^{\prime }\sqrt{t}\,R^{\delta }. \end{equation}

Following Detournay (Reference Detournay2016), the viscous scaling ($m$-scaling) can be identified by requiring that $p_{m}\sim p_{e}$, while the storage scaling can be identified by requiring that $V_{i}\sim V_{f}$, from which it follows that the length/radius $R_{m}$ and aperture $w_{m}$ scales for the viscous storage scaling are given by, respectively,

(4.28a,b)\begin{equation} R_{m}\sim \left( \frac{E^{\prime }Q_{0}^{3} t^{4}}{\mu^{\prime }}\right)^{{1}/({3\delta +3})}, \quad w_{m}\sim \left( \frac{\mu^{\prime \delta }Q_{0}^{3} t^{3-\delta }}{E^{\prime \delta }}\right)^{{1}/({3\delta +3})},\end{equation}

while the dimensionless toughness $\mathcal {K}_{m}$ and leak-off coefficient $\mathcal {C}_{m}$ become

(4.29a,b)\begin{align} \mathcal{K}_{m}:=\frac{p_{k}}{p_{e}}\sim \left( \frac{K^{\prime 14\delta -10}t^{2(\delta -1)}}{E^{\prime 10\delta -7}\mu^{\prime 4\delta -3} Q_{0}^{2\delta -1}}\right)^{{1}/({14\delta -10})},\quad \mathcal{C}_{m}:=\frac{V_{c}}{V_{i}}\sim \left( \frac{C^{\prime 12\delta-6}E^{\prime 3\delta -2}t^{6\delta -5}}{\mu^{\prime 3\delta -2} Q_{0}^{3\delta }}\right)^{{1}/({12\delta -6})}.\end{align}

The viscous leak-off scaling ($\tilde {m}$-scaling) can be obtained by requiring $V_{i}\sim V_{c}$ instead of $V_{i}\sim V_{f}$, from which it follows that the length/radius $R_{\tilde {m}}$ and aperture $w_{\tilde {m}}$ scales for the viscous leak-off scaling are given by, respectively,

(4.30a,b)\begin{equation} R_{\tilde{m}}\sim \left( \frac{Q_{0}^{2}t}{C^{\prime 2 }} \right)^{{1}/{2\delta }},\quad w_{\tilde{m}}\sim \left( \frac{\mu^{\prime 3\delta -2}Q_{0}^{3 \delta }t}{E^{\prime 3\delta -2}C^{\prime 2}}\right)^{{1}/({ 12\delta -8})}. \end{equation}

We observe from (4.29a,b) that the dimensionless toughness for a fracture driven by a constant injection rate $Q_{0}$ is independent of time for a KGD fracture $\delta =1$, and increases with time for a radial fracture $\delta =2$. The dimensionless leak-off coefficients for both KGD and radial fractures driven by a constant injection rate $Q_{0}$ increase with time. We also observe that the transition time $t_{mk}$ from viscosity- to toughness-dominated propagation (which is valid only in the radial case $\delta =2$, since $\mathcal {K}_{m}$ is time-independent in the KGD case $\delta =1$) and the transition time $t_{m\tilde {m}}$ from viscosity-storage-dominated propagation to leak-off-dominated propagation are given by

(4.31a,b)\begin{equation} t_{mk}=\left( \frac{E^{\prime 10\delta -7}\mu^{\prime 4\delta -3}Q_{0}^{2\delta -1}}{K^{\prime 14\delta -10}}\right)^{{1}/{2(\delta -1)}},\quad t_{m\tilde{m}}=\left( \frac{\mu^{\prime 3\delta -2}Q_{0}^{3\delta }}{ E^{\prime 3\delta -2}C^{\prime 12\delta -6}}\right)^{{1}/({6\delta -5})},\end{equation}

while the propagation regime parameter $\phi$ (valid only for the radial case) is defined to be

(4.32)\begin{equation} \phi = \frac{t_{mk}}{t_{m\tilde{m}}} =\left(\frac{E^{\prime 11}\mu^{\prime 3}C^{\prime 4}Q_{0}}{K^{\prime 14}}\right)^{{9}/{14}}. \end{equation}

4.6.2. Scaling for a radial fracture after shut-in

In practice, the evolution of a hydraulic fracture involves propagation due to the injection of a fluid at flux $Q_{0}$ (considered here to be constant) followed by shut-in at a certain time $t_{s}$, after which the fracture may continue to propagate (depending on the regime of propagation at shut-in) until ultimately it comes to rest either due to excessive leak-off or because the stress intensity factor has dropped below the critical fracture toughness. In the latter case, there is an arrest period during which the stress intensity factor $K$ decreases as the fracture continues to lose volume until $K=0$, at which point transition to recession is initiated and the fracture starts to recede. The appropriate scaling for the dynamics of a hydraulic fracture with a fixed injected volume $V_{0}$ can be obtained directly from those of a fracture driven by a constant flux $Q_{0}$ given in (4.28a,b) and (4.29a,b) by making the simple substitution $Q_{0}=V_{0}/t$. In this case, the length/radius $R_{m}^{V}(t)$ and aperture $w_{m}^{V}(t)$ scaling factors are given by

(4.33a,b)\begin{equation} R_{m}^{V}(t)\sim \left( \frac{E^{\prime }V_{0}^{3}t}{\mu^{\prime }}\right) ^{{1}/({3\delta +3})}\quad \text{and} \quad w_{m}^{V}(t)\sim \left( \frac{\mu^{\prime \delta }V_{0}^{3}}{E^{\prime \delta }t^{\delta }}\right)^{{1}/({ 3\delta +3})}. \end{equation}

Here, we have followed Mori & Lecampion (Reference Mori and Lecampion2021) and used the superscript $V$ to denote the scaling for a fracture with a fixed injected volume $V_{0}$ at time $t$. The corresponding leak-off ($\tilde {m}$-scaling) length/radius $R_{\tilde {m} }^{V}(t)$ and aperture $w_{\tilde {m}}^{V}(t)$ scalings for a hydraulic fracture with a fixed injected volume $V_{0}$ are

(4.34a,b)\begin{equation} R_{\tilde{m}}^{V}(t)\sim \left( \frac{V_{0}^{2}}{C^{\prime 2 }t}\right)^{{1}/{2\delta }}\quad \text{and} \quad w_{\tilde{m}}^{V}(t)\sim \left( \frac{\mu^{\prime 3\delta -2}V_{0}^{3\delta }}{E^{\prime 3\delta -2}C^{\prime 2}t^{3\delta -1}} \right)^{{1}/({12\delta -8})}. \end{equation}

The dimensionless toughness and leak-off parameters for a hydraulic fracture with a fixed injected volume $V_{0}$ are $\mathcal {K} _{m}^{V}(t)=( {t}/{t_{mk}^{V}})^{({4\delta -3})/({14\delta -10})}$ and $\mathcal {C}_{m}^{V}(t)=( {t}/{t_{m\tilde {m}}^{V}} )^{({9\delta -5})/({12\delta -6})}$, where

(4.35a,b)\begin{align} t_{mk}^{V}=\left( \frac{\mu^{\prime 4\delta -3}E^{\prime 10\delta -7}V_{0}^{2\delta -1}}{K^{\prime 14\delta -10}}\right)^{{1}/({4\delta -3})}\quad \text{and} \quad t_{m\tilde{m}}^{V}=\left( \frac{\mu^{\prime 3\delta -2}V_{0}^{3\delta }}{E^{\prime 3\delta -2}C^{\prime 12\delta -6}}\right)^{{1}/({9\delta -5})}. \end{align}

We define the following arrest regime parameter $\phi ^{V}$ for KGD and radial hydraulic fractures with a fixed injected volume to characterize the modes of arrest:

(4.36)\begin{equation} \phi^{V}=\frac{t_{mk}^{V}}{t_{m\tilde{m}}^{V}}= \left( \frac{\mu^{\prime 2\delta +1}E^{\prime 8\delta +5}C^{\prime 4\delta +2}V_{0}}{K^{\prime 10\delta +6}}\right)^{({8\delta -7})/({61\delta -57})}.\end{equation}

We observe that the parameter $\phi ^{V}$ defined in (4.36) has no meaning in the zero toughness case since $t_{mk}^{V}=\infty$.

4.6.3. Characteristic power law for the duration of recession

Since recession does not depend directly on the rock toughness $K^{\prime }$, we expect the recession phase to be largely independent of $\phi ^{V}$. Because the recession process that we consider is driven by leak-off, which persists through propagation, arrest and recession, we choose to define the dimensionless shut-in time $\omega$ as the ratio of the shut-in time $t_{s}$ to the storage–leak-off transition time $t_{m\tilde {m}}$, i.e.

(4.37)\begin{equation} \omega =\frac{t_{s}}{t_{m\tilde{m}}} .\end{equation}

Now making use of (4.31a,b) and (4.35a,b), the following relationship can be established between the fixed injected volume transition time $t_{m\tilde {m}}^{V}$ and the constant injection rate storage–leak-off transition time $t_{m\tilde {m}}$ in terms of the dimensionless shut-in parameter $\omega$:

(4.38)\begin{equation} t_{m\tilde{m}}^{V}=t_{m\tilde{m}}\omega^{{3\delta }/({9\delta -5})} =t_{s}\omega^{({5-6\delta })/({9\delta -5})}, \end{equation}

where the second relationship in (4.38) comes directly from the definition of $\omega$.

Modes of arrest and recession: Using (4.38) and (4.36), we now characterize the way in which the arrest time $t_a$ and recession time $t_r$ are impacted by the relative magnitudes of $\phi ^{V}$ and $\omega$.

1. (i) Consider $\phi ^{V}\ll 1$.

   1. (a) If $\omega \ll 1$, then it follows from (4.38) and (4.36) that

      (4.39)\begin{equation} t_{s}\ \&\ t_{mk}^{V}\ll t_{m\tilde{m}}^{V}\ll t_{m\tilde{m}}, \end{equation}
      so after shut-in, arrest will be determined by $t_{a}\sim t_{mk}^{V}$, and recession will be determined by $t_{r}\sim t_{m\tilde {m}}^{V}$.

   2. (b) If $\omega \gg 1$, then it follows from (4.38) and (4.36) that

      (4.40)\begin{equation} t_{m\tilde{m}}\ \&\ t_{mk}^{V}\ll t_{m\tilde{m}}^{V}\ll t_{s}, \end{equation}
      so after shut-in, arrest and recession will be determined by $t_{a}\sim t_{r}\sim t_{s}$.

2. (ii) Consider $\phi ^{V}\gg 1$.

   • (a) If $\omega \ll 1$, then it follows from (4.38) and (4.36) that

     (4.41)\begin{equation} t_{s}\ll t_{m\tilde{m}}^{V}\ll t_{m\tilde{m}}\ \&\ t_{mk}^{V}, \end{equation}
     so after shut-in, arrest and recession will be determined by $t_{m \tilde {m}}^{V}$, i.e. $t_{a}\sim t_{r}\sim t_{m\tilde {m}}^{V}$.

   • (b) If $\omega \gg 1$, then it follows from (4.38) and (4.36) that

     (4.42)\begin{equation} t_{m\tilde{m}}\ll t_{m\tilde{m}}^{V}\ll t_{s}\ \&\ t_{mk}^{V}, \end{equation}
     so after shut-in, arrest and recession will be determined by $t_{a}\sim t_{r}\sim t_{s}$.

Typical field values of $\phi ^V$ and $\omega$:

We have seen that the two dimensionless parameters $\phi ^V$ and $\omega$ characterize the fracture arrest and the nature of transition to recession, respectively. To establish typical field values for these parameters, assume the following ranges of material parameters that are encountered typically in the field: $E^{\prime }\sim 1\unicode{x2013}30\ \mathrm {GPa}$, $\mu ^{\prime }\sim 10^{-2}\unicode{x2013}10 \mathrm {Pa}\ \mathrm {s}$, $C^{\prime }\sim 10^{-5}\unicode{x2013}10^{-8} \mathrm {m\ s}^{-1/2}$, $K^{\prime }\sim 0.3\unicode{x2013}3 \mathrm {MPa}\, \mathrm {m}^{1/2}$, $t_{s}\sim 3600\ \mathrm {s}$ and $Q_{0}\sim 10^{-5}\unicode{x2013}10^{-3}\ \mathrm {m}^{2}\ \mathrm {s}^{-1}$ for KGD hydraulic fractures, and $Q_{0}\sim 10^{-3}\unicode{x2013}10^{-1}\ \mathrm {m}^{3}\ \mathrm {s}^{-1}$ for radial hydraulic fractures. For these ranges of material and injection parameters, the ranges of the dimensionless parameters for KGD hydraulic fractures are $10^{-23}\lesssim \omega _{{KGD}}\lesssim 10$ and $10^{-11}\lesssim \phi _{{KGD}}^{V}\lesssim 10^{4}$, while for radial hydraulic fractures they are $10^{-10}\lesssim \omega _{{RAD}}\lesssim 3$ and $10^{-10}\lesssim \phi _{{RAD}}^{V}\lesssim 10^{4}$.

4.6.4. Scaling law for the duration of recession $t_{c}-t_{r}$

In this subsubsection, we use the sunset solution to provide the estimate

(4.43)\begin{equation} t_{c}-t_{r}\sim \frac{w_{r}\sqrt{t_{r}}}{C^{\prime }},\end{equation}

where $t_{c}$ is the closure time, $t_{r}$ is the time of initiation of recession, and $w_{r}$ is an estimate of the fracture aperture at the start of recession.

Small shut-in time $\omega \ll 1$:

Since in this limit $t_{r}\sim t_{m\tilde {m}}^{V}$, we use (4.33a,b) to provide the estimate $w_{r}=w_{m}^{V}(t_{m\tilde {m}}^{V})\sim$ $w_{s} \omega ^{{\delta ( 6\delta -5) }/{3( \delta +1) ( 9\delta -5) }}$, where $w_{s}$ is the aperture at shut-in, which from ( 4.28a,b) has the estimate $w_{s}\sim C^{\prime }t_{m\tilde {m} }^{1/2}\omega ^{{1}/({6\delta -3})}$. Now substituting these estimates into (4.43) and using (4.38), we obtain the following scaling for the duration of recession:

(4.44)\begin{equation} t_{c}-t_{r}\overset{\omega \ll 1}{\sim } t_{m\tilde{m}} \omega^{{3\delta}/({9\delta -5})}.\end{equation}

Comparing (4.44) and (4.38), we see that $t_{c}\sim 2t_{r}$, which implies that at the time recession starts, only some fraction (roughly half) of the fluid has leaked off, so that the time $t_{c}-t_{r}$ that the fracture takes to recede is of the same order as the time it took to get to the point of recession $t_{r}$. In order to demonstrate this further, we consider the asymptotics of the efficiency $\eta (t)$ defined to be the ratio of the volume of fluid enclosed in the fracture $V_f(t)$ at any time $t$ to the volume of fluid that has been injected into the fracture $V_i(t)$ until time $t$, i.e. $\eta (t):={V_f(t)}/{V_{i}(t)}$, where $V_{i}=Q_{0}[ t( 1-H(t-t_{s})) +t_{s}\,H(t-t_{s})]$ and $H(t)$ is the Heaviside step function. We now use (4.34a,b) to obtain the following estimate for $\eta _{r}$:

(4.45)\begin{equation} \eta _{r}\sim \frac{w_{\tilde{m}}^{V}(t_{\tilde{m}}^{V})R_{\tilde{m}}^{V}(t_{ \tilde{m}}^{V})^{\delta }}{V_{0}}=1. \end{equation}

Now, we know that $\eta (t)<1$ for $t>0$, so the estimate $\eta _{r}= O(1)$ established in (4.45) should be interpreted as demonstrating that the efficiency at the start of recession, $\eta _r<1$, does not depend on $\omega$ for small values of $\omega$.

Large shut-in time $\omega \gg 1$: Since in this limit $t_{r}\sim t_{s}$, we use (4.30a,b) to provide the estimate $w_{r}\sim w_{s}\sim C^{\prime }t_{m\tilde {m} }^{1/2}\omega ^{{1}/({12\delta -8})}$, where $w_{s}$ is the aperture at shut-in. Now substituting these estimates into (4.43), we obtain the following scaling for the duration of recession:

(4.46)\begin{equation} t_{c}-t_{r}\overset{\omega \gg 1}{\sim }t_{m\tilde{m}} \omega^{({6\delta -3})/({12\delta -8})}.\end{equation}

We use (4.30a,b) to obtain the following estimate for the efficiency $\eta _{r}$:

(4.47)\begin{equation} \eta _{r}\sim \eta (t_{s})\sim \frac{w_{\tilde{m}}(t_{s})\,R_{\tilde{m} }(t_{s})^{\delta }}{V_{0}}=\omega^{({5-6\delta})/({12\delta -8})}. \end{equation}

This establishes the rate at which the efficiency at the start of recession $\eta _r$ decreases to zero for large values of $\omega$. We also observe from (4.38), (4.46) and (4.47) that $\eta _{r}\sim \omega ^{({5-6\delta })/({12\delta -8})}\sim ({t_{c}-t_{r}})/{t_{s}}$.

Comparison of the scaling for the duration of recession with numerical results: The solid black curves in figure 3(a) indicate the duration of recession to storage–leak-off transition time ratios $(t_c-t_r)/t_{m\tilde {m}}$ for KGD fractures $\delta =1$ all plotted as functions of $\omega$ for a range of values of the regime parameter $\phi ^{V}=10^{j}$, $j\in \{-3,-2,-1,-\frac {3}{4},-\frac {1}{2},-\frac {1}{4},0, \frac {1}{4},\frac {1}{2},\frac {3}{4},1,2,3\}$. We observe that the duration of recession in this case, $t_{c}-t_{r}$, is independent of $\phi ^V$. The dashed red line represents the log–linear regression of the case $\phi ^V=1$ using the first few data points, which has gradient $0.746$, in close agreement with the exponent $3/4$ determined by the scaling in the limit $\omega \ll 1$ given in (4.44). The scaling in the limit $\omega \gg 1$ given in (4.46) also yields the exponent $3/4$ when $\delta =1$, which is consistent with the numerical result shown in figure 3(a). The solid black curves in 3(b) indicate the duration of recession to storage–leak-off transition time ratios $(t_c-t_r)/t_{m\tilde {m}}$ for radial fractures $\delta =2$ all plotted as functions of $\omega$ for the range of values of the regime parameter $\phi ^V \in \{0.05,0.1,0.2,0.5,1,1.25, 1.5,2,3\}$. We observe that the duration of recession $t_{c}-t_{r}$ is also independent of $\phi ^V$ for this case. The dashed red line has gradient $0.46$, close to the exponent $6/13$ obtained from the small $\omega$ scaling given in (4.44) when the value $\delta =2$ is used. The dash-dotted blue line has gradient $0.54$, close to the exponent $9/16$ obtained from the large $\omega$ scaling given in (4.44) with $\delta =2$.

Figure 3. The solid black lines indicate the duration of recession to storage–leak-off transition time ratios $(t_c-t_r)/t_{m\tilde {m}}$ plotted as functions of $\omega$ for a range of values of the regime parameter $\phi ^V$. The KGD case $\delta =1$ is plotted in (a), and the dashed red line represents the log–linear regression of the case $\phi ^V=1$ using the first few data points. The radial case $\delta =2$ is plotted in (b), and the dashed red and dash-dotted blue lines represent the log–linear regressions of the case $\phi ^V=1$ using the first/last few data points.

In figure 4(a), we plot the efficiency at the initiation of recession $\eta _r$ as a function of $\omega$ for the KGD case $\delta =1$ for three selected values of the arrest regime parameter $\phi ^{V}=10^{j}$, $j\in \{-2,0, 3\}$. The curve corresponding to $\phi ^{V}=10^{-2}$ is indicated by the ${\bullet }$ symbol at the abscissa $\omega =10^{-12}$, and for the other two curves, $\eta _r$ increases with increasing $\phi ^{V}$. We observe that $\eta _r$ develops a horizontal asymptote for small values of $\omega$, which is consistent with the scaling result (4.45). A log–linear regression, sampled for the few largest values of $\omega$, yields gradient $-0.2$, which is not far from the large $\omega$ scaling estimate (4.46), which yields the exponent $-1/4$ for $\delta =1$. The underestimation of the gradient by the numerical scheme for large $\omega$ is due to the limit on the feasibility of exploring much smaller values of $\eta _r$ associated with large $\omega$ numerically, which would be required for the asymptote to become fully developed. In figure 4(b), we plot the efficiency at the initiation of recession $\eta _r$ as a function of $\omega$ for the radial case $\delta =2$ for three selected values of the arrest regime parameter $\phi ^V \in \{0.05,0.1,3\}$. The curve corresponding to $\phi ^{V}=0.05$ is indicated by the ${\bullet }$ symbol at the abscissa $\omega =10^{-7}$, and for the other two curves, $\eta _r$ increases with increasing $\phi ^{V}$. We observe that $\eta _r$ develops a horizontal asymptote for small values of $\omega$, which is consistent with the scaling result (4.45). A log–linear regression, sampled for the few largest values of $\omega$, yields gradient $-0.38$, which is not far from the large $\omega$ scaling estimate (4.46), which yields the exponent $-7/16$ for $\delta =2$. As in the KGD case, the numerical estimate for this exponent underestimates the exponent established by scaling because of the challenge in exploring values $\eta _r \ll 1$.

Figure 4. The efficiency $\eta _r$ plotted as functions of $\omega$ for three selected values of the regime parameter $\phi ^V$. The KGD case $\delta =1$ is plotted in (a) for $\phi ^{V}=10^{j},\ j\in \{-2 ( \text {indicated by } {\bullet }),0, 3\}$. The radial case $\delta =2$ is plotted in (b) for $\phi ^V \in \{0.05 ( \text {indicated by } {\bullet }),0.1,3\}$.