Past exam of the mathematics course of the University of Cambridge / 2014 / iii / Paper 78 (Fluid Dynamics of Energy Systems)

Effective permeability of a porous medium. Work to leading order in the slender-layer ratio $H/L$. Take the growing lower wedge to have porous permeability $k_1$, with interface $y=Hx/L$, and the upper wedge to have porous permeability $k_2$. The leading pressure is independent of $y$; vertical flow is smaller than horizontal flow by $H/L$. DefineDarcy's law and the fixed two-dimensional volume flux giveHere $q_i$ is Darcy velocity; the parcel speed is $q_i/\varphi$. Integrating the pressure gradient over the length and defining $Q=H{k}_{\mathrm{eff}}\Delta p/(\mu L)$ yieldsThe equal-porous permeability limit is $k_{\mathrm{eff}}=k_1=k_2$. This logarithmic mean comes from parallel layers at each cross-section followed by series addition of their local hydraulic resistances. It is a slender-layer result; a full two-dimensional transmission problem has small end and interface corrections.

Parcel paths and travel times. Pressure equalization does not mean that parcels stay at fixed $y$: mass conservation requires a small vertical flow. Define a streamfunction by $\psi (x,0)=0$ and $\psi (x,H)=Q$. Its leading expression isLet $r=y_0/H$ label a parcel released at the inlet. There $\psi =Qr$. At the inclined interface $\psi /Q=k_1\xi /K(\xi )$, so the parcel crosses from the upper wedge into the lower one atBefore crossing, $y/H=1-(1-r)K/k_2$; afterwards, $y/H=rK/k_1$. These expressions show why integrating the speed along a horizontal line would give the wrong parcel time.

Put $d=k_1-k_2$ and ${\tau }_V=\varphi HL/Q$, the pore-volume throughput time. Integrating $dt=\varphi dx/q_i$ on the two portions of the path givesIts derivative with respect to ${\xi }_c$ is ${\tau }_{V}dK/(k_1{k}_2)$, so it is monotone. The limiting streamline times at the lower and upper boundaries are ${\tau }_V(k_1+k_2)/(2k_1)$ and ${\tau }_V(k_1+k_2)/(2k_2)$ respectively. HenceThe printed expression assumes $k_1\ge k_2$. The absolute value is needed for a nonnegative maximum difference without that ordering. For equal porous permeabilities every parcel has time ${\tau }_V$.

Oil recovery. In the ideal passive-displacement model with $k_1/k_2=10$, the earliest and latest travel times are $0.55{\tau }_V$ and $5.5{\tau }_V$, a spread of $4.95{\tau }_V$. Preferential paths through the high-porous permeability wedge therefore give early water breakthrough while oil on slower paths remains unswept. Continued injection sends much water through paths already swept; complete displacement requires several pore volumes. The flux-weighted residence time in a porous layer remains ${\tau }_V$. Real waterflooding also depends on phase mobilities, relative permeabilities, capillary pressure and mixing, so these numbers illustrate heterogeneity rather than a quantitative two-phase recovery prediction.

Choose $U$ and $Q$ as the speed and source strength in the same kinematic velocity convention. With $r=(x^2+y^2{)}^{1/2}$ and $\vartheta =\arg (x+iy)$, superposition of uniform potential flow and a two-dimensional point source giveswhere $(u,v)=\nabla \Phi =({\psi }_y,-{\psi }_x)$. The angular branch cut records the volume flux from the source; the logarithm's arbitrary reference length only changes an additive constant in the velocity potential.

The stagnation point lies on the negative axis, where $u=U+Q/(2\pi x)=0$. Thus the limiting upstream reach of source fluid isThe stagnation point is approached asymptotically, rather than reached in finite time by a parcel on the upstream axis. The dividing streamlines through it have $\psi =\pm Q/2$ when the angle is taken in $(-\pi ,\pi )$. Far downstream $\vartheta \to 0$, so these streamlines approach $y=\pm Q/(2U)$. Consequently the source fluid occupies the Rankine half-body between them, with downstream width $2b=Q/U$. In the upper half-plane its boundary can also be parametrized by $y=b(1-\vartheta /\pi )$, $x=y\cot \vartheta$, $0<\vartheta <\pi$.

For arrival at the downstream axis, the direct streamline $y=0$, $x>0$ has speed $U+Q/(2\pi x)$. Its transit time from the ideal point source isThe first source fluid at a downstream cross-section lies on this symmetry axis: for positive $x$, the horizontal source contribution $Qx/[2\pi (x^2+y^2)]$ is maximal at $y=0$, and off-axis paths cannot advance through positive $x$ more quickly. A far-away parcel passing the inlet line at the same release time has transit time $x/U$. The time advance is thereforeThis is an advance, not the source parcel's absolute transit time. If a source parcel is released at time $t_0$ while the distant comparison parcel crosses $x=0$ at time zero, the advance relative to that parcel is $\tau -t_0$; the printed formula uses simultaneous releases. If $U,Q$ instead denote Darcy velocity and bulk-area source flux in a medium of uniform porosity $\varphi$, their ratio still gives the same $b$, but both parcel times acquire a factor $\varphi$. No porosity is specified here, so the displayed time formula fixes the kinematic convention.

The imposed pressure gradient and Darcy's law give a linear velocity profile. Use $U$ as the mean pore velocity in the contaminant advection-diffusion equation; thenThe last conditions express zero transverse solute flux under the standard closed-layer interpretation. Transverse exchange through the layer boundaries would change this cell problem and requires additional boundary data. If the prescribed mean speed is a Darcy velocity, its corresponding pore velocity is that speed divided by porosity.

The given parameter is the ratio of transverse mixing time to travel time:Thus transverse dispersion mixes the contaminant across the layer many times during its passage along the long layer. Different fluid speeds still matter: repeated exchange between slow and fast portions produces Taylor dispersion, rather than transport solely at each parcel's original speed. The averaged description applies after $t\gg H^2/D_T$ and on longitudinal scales large compared with the distance travelled during that mixing time.

Write $C=\overline{c}$ and, to leading correction order, $c=C+\chi (y)C_x+\cdots$, with $\overline{\chi }=0$. Taking $C_t\simeq -U{C}_x$ in the local equation gives the cell problemWith $\eta =y/H$, its solution isAveraging the local advection-diffusion equation now yields $C_t+U{C}_x=[D_L-\overline{(u-U)\chi }]C_{xx}$. An integration by parts in the cell problem gives $-\overline{(u-U)\chi }=D_T\overline{({\chi }^{\prime }{)}^2}$. Since ${\int }_0^1({\eta }^2-\eta {)}^{2}d\eta =1/30$, the effective equation isThis is Taylor dispersion in a linear porous-layer velocity profile. For a localized pulse away from the inlet and outlet, its mean position advances at $U$ and its longitudinal variance grows as $2D_{\mathrm{eff}}t$.

For the inlet step, use the semi-infinite inlet approximation $x>0$, $C(0,t)=1$, $C(x,0)=0$ and $C\to 0$ as $x\to \infty$ at fixed $t$. The resulting constant-concentration inlet solution isOne derivation is to Laplace transform in time: $\tilde{C}(x,p)=p^{-1}\exp [(U-\sqrt{U^2+4D_{\mathrm{eff}}p})x/(2D_{\mathrm{eff}})]$. Inverting gives the two complementary error functions. Their sum is one at $x=0$; for fixed $x>0$ both terms vanish as $t\downarrow 0$, and substitution verifies the averaged equation. Retaining only the first term gives the familiar moving error-function front far from the inlet when longitudinal advection dominates diffusion, but does not satisfy the inlet condition exactly.

For a literal finite layer $0<x<L$, a downstream boundary condition is also needed once the outlet influences the solution. The PDF's step-inlet request specifies only the initially clean region $x>0$, so the formula above is the semi-infinite/long-layer solution, not a claim of a unique finite-interval solution with unspecified outlet data.

Saturation equation. Introduce intrinsic permeability of a porous medium $K$ and porosity $\varphi$. The relative porous permeabilities depend on the wetting fluid saturation. Define the phase mobilities and fractional flowUse the conventional capillary pressure $p_c=p_{nw}-p_w$, generally decreasing as wetting saturation increases. Neglect gravity in this horizontal model. The two Darcy fluxes are $q_w=-K{\lambda }_w{p}_{w,x}$ and $q_n=-K{\lambda }_n(p_{w,x}+p_c^{\prime }(s)s_x)$. Eliminate the common pressure gradient using $q_w+q_n=Q$ to obtainThe wetting-phase mass conservation law is consequentlyFor $p_c^{\prime }\le 0$, $D_c\ge 0$ and the capillary pressure provides nonlinear diffusion. The initially uniform saturation is $s_0$. Pure wetting-fluid injection prescribes $q_w(0,t)=Q$ and $q_n(0,t)=0$; in the zero-capillarity, zero-residual-nonwetting idealization this corresponds to inlet saturation one. With residual nonwetting fluid, replace one by the appropriate maximum accessible saturation.

Shock formation and its speed. Neglect capillarity first. Smooth saturation values travel along characteristic curves at speed $(Q/\varphi )F^{\prime }(s)$. When the trailing values have larger characteristic speeds than the values ahead, the curves intersect, and the Buckley-Leverett equation requires a shock selected as an entropy solution. This is typical for the increasing-convex part of a physical fractional-flow curve; shocks are not inevitable for every possible $F$.

Integrating mass conservation across a shock from upstream $s_s$ to downstream $s_0$ gives the Rankine-Hugoniot conditionIf an upstream rarefaction wave joins the shock smoothly in similarity coordinates, its terminal characteristic speed equals $V_s$. Therefore the saturation at that join obeys the fractional-flow tangent constructionThis is the printed equality. It is a tangent-selection condition, not the dimensional speed by itself; the speed includes $Q/\varphi$. An arbitrary shock between prescribed constant states satisfies the jump condition but need not satisfy this additional tangency relation.

Positive capillary diffusion spreads a jump into a transition layer. For a travelling profile $s(x-V_{s}t)$, integration givesThe endpoint states still satisfy the same Rankine-Hugoniot condition. Thus small capillarity primarily gives the front a finite thickness rather than changing the selected limiting speed; vanishing phase mobility at an endpoint can make the regularization degenerate.

The specified reciprocal fractional flow. On its stated branch, $F^{\prime }=1/s^2$ and $F^{\prime \prime }=-2/s^3<0$. If $F$ extends continuously to the initial saturation and the inlet is at one, the decreasing saturation from behind to ahead makes characteristic speeds increase forwards: there is a rarefaction wave, not a compressive shock. Write $a=Q/\varphi$ and $\xi =x/t$. The self-similar solution isThe fan connects continuously to the initial state and broadens linearly in time. This interpretation requires $F(s_0)=2-1/s_0$. A nonnegative physical fractional flow on the entire interval $[s_0,1]$ also requires $s_0\ge 1/2$.

The PDF specifies the formula only for $s>s_0$ and earlier calls $s_0$ residual. If that means an immobile wetting phase with $F(s_0)=0$, the endpoint is an extra constraint and the continuous-branch answer applies without conflict only at $s_0=1/2$. For $0<s_0<1/2$, the printed branch would give negative fractional flow for $s_0<s<1/2$, so it cannot be the complete physical constitutive law. With an additional admissible extension, for example $F=0$ up to $1/2$ and the given positive branch above it, a compound fan–shock is possible. The tangent condition then givesThe fan $s=\sqrt{a/\xi }$ ends at $\xi =V_s$, followed by a jump to $s_0$. This is a conditional physical extension, not data supplied by the question. For $s_0>1/2$ an immobile residual endpoint is likewise incompatible with a continuous version of the stated branch. These distinctions identify what is, and is not, determined by the printed constitutive assumption.

Let $K$ be intrinsic permeability of a porous medium, $\varphi$ the mobile pore fraction, and let $h$ measure thickness normal to the sloping cap. Define cold and reservoir-temperature properties byThe imposed inequality makes ${\delta }_c>0$, so the cold fluid is still buoyant. Assume $K>0$, $\sin \theta >0$, positive densities and viscosities, and the usual positive thermal-expansion and viscosity coefficients. Away from localized fronts, neglect the along-current thickness gradient. Hydrostatic water pressure and Darcy's law then give the upslope Darcy velocitiesTypically $u_h>u_c$: warming increases buoyancy and decreases dynamic viscosity.

Current depths and the moving moving thermal front in a porous current. For this slope-driven porous gravity current, the injected mass flux per unit well length fixes the cold plateau thickness:The moving thermal front in a porous current moves at $V_T=\Gamma u_c$. Since $u_c$ is explicitly a Darcy velocity, the fluid mass stored per unit horizontal area is $\varphi \rho h$. Apply mass conservation across the moving moving thermal front in a porous current:Writing $\gamma =\varphi \Gamma$ and $r=u_c/u_h$, this gives the distal warm depthFor a retarded moving thermal front in a porous current take $0\le \gamma <1$; otherwise this two-plateau geometry needs reconsideration. The factor $(1-\gamma )/(1-\gamma r)$ is important: the distal mass flux is not generally $Q$, because mass is being stored as the cold region replaces the warm region. Setting the two fluxes equal would silently assume a stationary moving thermal front in a porous current.

If the model absorbs porosity into storage and uses pore velocity as its transport speed, the same formulas use $\gamma =\Gamma$. The explicit Darcy convention in the PDF instead gives $\gamma =\varphi \Gamma$; porosity must be specified or absorbed consistently. In the stationary-front limit $\gamma =0$, the two steady depths reduce to $Q/({\rho }_c{u}_c)$ and $Q/({\rho }_h{u}_h)$. In the Boussinesq approximation, use a common reference density in the mass factors while retaining ${\delta }_c,{\delta }_h$ in the driving force.

Leakage thresholds. A normal thickness $h_i$ generates cap overpressure ${\Pi }_i={\delta }_{i}g{h}_i\cos \theta$. Before leakage modifies the current, the two plateau values areFor the caprock leakage threshold $\Delta p>0$, the candidate critical mass-injection rates are thereforeCold-only leakage occurs for $Q_c^{\ast }\le Q<Q_h^{\ast }$, provided $Q_c^{\ast }<Q_h^{\ast }$. Both plateaus can leak once $Q\ge \max (Q_c^{\ast },Q_h^{\ast })$. In the usual common-density approximation, ${\mu }_c>{\mu }_h$ and $r<1$ ensure $Q_c^{\ast }<Q_h^{\ast }$, giving the expected sequence: no leakage, cold-only leakage, then leakage from both cold and warm regions. Without that approximation the ordering must be checked; the given inequality $\Delta \rho >\alpha \Delta T$ alone does not establish it. These are onset criteria computed on the nonleaking current, not a post-leakage mass budget.

For small slopes $\cos \theta \simeq 1$ and $\tan \theta \simeq \sin \theta \simeq \theta$. If another thickness convention is used, its hydrostatic column and projected flux must be changed consistently; one should not mix a normal thickness with a vertical-pressure formula lacking the cosine.

Cross-current heat conduction. Heat transfer from the warm formation makes temperature vary across the carbon dioxide depth and introduces a warming time controlled by thickness, thermal diffusivity and the surrounding rock's heat capacity. It smooths the sharp thermal adjustment and causes gradual changes in buoyancy, viscosity and the velocity profile. The cold fluid generally warms, becomes more mobile and requires less depth to carry a prescribed mass flux; the enhanced cold-region overpressure and its distinct leakage zone tend to shrink. Warming may occur before a parcel reaches the idealized advective moving thermal front in a porous current, especially for a thin current. A single supplied $\Gamma$ no longer describes all heat transport. Quantitative depths and thresholds then require a coupled temperature equation and thermal boundary data; their exact changes cannot be inferred from the linear property laws alone.

Introduce water dynamic viscosity $\mu$ and the net driving liquid-pressure difference $\Delta P=p_0-p_f>0$. If $p_c>0$ denotes capillary suction relative to an ambient gas pressure $p_a$, then $p_f=p_a-p_c$ and $\Delta P=p_0-p_a+p_c$. If $p_c$ is instead the prescribed liquid front pressure, use $\Delta P=p_0-p_c$. This makes the capillary sign convention explicit.

Hemispherical front and injection flux. A hemispherical surface at radius $r$ has area $2\pi r^2$. Quasi-steady mass conservation makes the total volume flux independent of $r$ in the saturated region, so Darcy's law givesThe changing pore volume is $2\pi \varphi (R^3-R_s^3)/3$, giving $Q=2\pi \varphi R^2\dot{R}$. HenceFor the initially dry idealization $R(0)=R_s$, integration gives the implicit front lawIts right-hand side increases for $R>R_s$, defining one advancing radius and thus the flux implicitly in time. Equivalently, with $r=R/R_s$ and $t_{\ast }=\mu \varphi R_s^2/(k\Delta P)$, $t/t_{\ast }=r^3/3-r^2/2+1/6$. The infinite initial flux is an idealization of zero wetted resistance; inertia, finite inlet geometry and pore-scale effects regularize its earliest stage.

For a thin wetted shell $\ell =R-R_s\ll R_s$, the law reduces toFor $R\gg R_s$, the hemispherical capillary imbibition model instead predictswith $Q/Q_{\infty }=1+R_s/R+O((R_s/R{)}^2)$. Radial spreading makes the hydraulic resistance approach a finite value, while the increasing front area reduces its advance speed.

Comparison with one-dimensional imbibition. Under the same constant-driving-pressure, negligible-gravity assumptions, a planar front at depth $z$ has Darcy flux $q=k\Delta P/(\mu z)$ and storage law $\varphi \dot{z}=q$. ThereforeA fixed-area one-dimensional flux decays as $t^{-1/2}$, whereas the hemispherical total flux tends to a constant and its radius grows as $t^{1/3}$. If purely vertical imbibition includes gravity, upward penetration eventually approaches a capillary-rise height and downward penetration has a different gravity-driven limit. The hemispherical assumption itself neglects that directional gravity effect, so the power-law comparison is within its capillary-dominated regime.

Evaporation-limited steady radius. Take $F_e$ as volumetric water loss per unit front area per unit time. At equilibrium the incoming flux equals $2\pi R^2{F}_e$, givingIf the stated loss is a mass flux, replace $F_e$ by $F_e/{\rho }_w$. Porosity cancels from this equilibrium but affects the time to approach it. The incoming flux decreases with $R$ while the evaporating area increases, so the equilibrium is stable in this model. As $F_e\to 0$, $R_{\ast }\sim [k\Delta P{R}_s/(\mu F_e){]}^{1/2}$; sufficiently large radii can invalidate the negligible-gravity or semi-infinite homogeneous-medium assumptions.