Transient growth and symmetrizability in rectilinear miscible viscous fingering

Abstract

The influence of dispersion or equivalently of the Péclet number (Pe) on miscible viscous fingering in a homogeneous porous medium is examined. The linear optimal perturbations maximizing finite-time energy gain is demonstrated with the help of the propagator matrix approach based non-modal analysis (NMA). We show that onset of instability is a monotonically decreasing function of Pe and the onset time determined by NMA emulates the non-linear simulations. Our investigations suggest that perturbations will grow algebraically at early times, contrary to the well-known exponential growth determined from the quasi-steady eigenvalues. One of the over-arching objective of the present work is to determine whether there are alternative mechanisms which can describe the mathematical understanding of the spectrum of the time-dependent stability matrix. Good agreement between the NMA and non-linear simulations is observed. It is shown that within the framework of $L^2$-norm, the non-normal stability matrix can be symmetrizable by a similarity transformation and thereby we show that the non-normality of the linearized operator is norm dependent. A framework is thus presented to analyze the exchange of stability which can be determined from the eigenmodes.

Appendices

Appendix A: Linearized equations in (x, y, t) coordinates

In order to investigate the stability of the system, we introduce infinitesimal perturbations, $f', \; f \in \{u, v, p, c, \mu \}$ to the base state, Eqs. (9)–(10) such that $f = f_\mathrm{b}+ f'$. On substituting these in Eqs. (4)–(6), we have

$$\begin{aligned}&\frac{\partial u'}{\partial x} + \frac{\partial v'}{\partial y} = 0, \quad \frac{\partial p'}{\partial x} = -\mu _\mathrm{b} u' - \mu ', \quad \frac{\partial p'}{\partial y} = -\mu _\mathrm{b} v', \end{aligned}$$

(41)

$$\begin{aligned}&\frac{\partial c'}{\partial t} +\frac{\partial c_\mathrm{b}}{\partial x}u' = \frac{1}{\text {Pe}}\bigg (\frac{\partial ^2c'}{\partial x^2} + \delta \frac{\partial ^2c'}{\partial y^2}\bigg ), \end{aligned}$$

(42)

$$\begin{aligned}&\mu ' = \bigg (\frac{\text {d}\mu }{\text {d}c}\bigg |_{c_\mathrm{b}}\bigg )c'. \end{aligned}$$

(43)

Tan and Homsy [12] observed that the validity of QSSA in (x, y, t) is questionable at short times where the base state changes rapidly. A fundamental problem with QSSA approach is that the concentration eigenfunctions are localized around the base state while the eigenfunctions of the operator $\displaystyle {\partial ^2}/{\partial x^2}, \; x \in (-\infty , \infty )$ are global modes. Hence, they do not provide an appropriate basis for streamwise perturbations. Since the present problem can have a small onset time ($t \approx \mathcal {O}(1)$) for instability, the QSSA as such is ill-approximation to the task of resolving the early-time behavior of perturbations. To overcome this difficulty, we have done our analysis in the self-similar coordinate system $(\xi ,y,t)$, where the streamwise operator

$$\begin{aligned} \frac{\partial ^2}{\partial \xi ^2} + \frac{\xi }{2} \frac{\partial }{\partial \xi } \end{aligned}$$

has eigenfunctions as Hermite polynomials which are localized around the base-state [16, 17, 29].

Appendix B: Propagator matrix and measure of energy norm

The non-autonomous system, Eq. (23), can be solved using the propagator matrix approach. For a chosen time interval $[t_\mathrm{p}, t_\mathrm{f}]$, let $\varPhi (t_\mathrm{p};t_\mathrm{f})$ be a formal solution of Eq. (23), where $\phi _\mathrm{c}(t_\mathrm{f}) = \varPhi (t_\mathrm{p};t_\mathrm{f})\phi _\mathrm{c}(t_\mathrm{p})$, $\phi _\mathrm{c}(t_\mathrm{p})$ being an arbitrary initial condition. Substitute this value of $\phi _\mathrm{c}(t_\mathrm{f})$ in Eq. (23) to get a matrix-valued IVP

$$\begin{aligned} \frac{\text {d}}{\text {d}t}\varPhi (t_\mathrm{p};t_\mathrm{f}) = \mathcal {L}(t_\mathrm{f})\varPhi (t_\mathrm{p};t_\mathrm{f}), \end{aligned}$$

(44)

with initial condition $\varPhi (t_\mathrm{p};t_\mathrm{p}) = \mathcal {I}$, where $\mathcal {I}$ is the identity matrix. Assuming that $\phi _\mathrm{c}$ is square integrable on $\mathbb {R}$ we wish to find the maximum perturbation energy gain, that is,

$$\begin{aligned} G(t_\mathrm{f})=G(t_\mathrm{f}, k, \text {Pe}, R, \delta ):=\max _{\phi _\mathrm{c}(t_\mathrm{p})} \frac{E_{\phi _\mathrm{c}}(t_\mathrm{f})}{E_{\phi _\mathrm{c}}(t_\mathrm{p})} = \max _{\phi _\mathrm{c}(t_\mathrm{p})} \Vert \varPhi (t_\mathrm{p},t_\mathrm{f})\phi _\mathrm{c}(t_\mathrm{p})\Vert _2= \Vert \varPhi (t_\mathrm{p}, t_\mathrm{f})\Vert = \displaystyle \sup _{j} s_j(t_\mathrm{f}), \end{aligned}$$

(45)

where $s_j$’s are the singular values of $\varPhi (t_\mathrm{p}; t_\mathrm{f})$, in other words, the eigenvalues of the self-adjoint matrix $\varPhi ^*(t_\mathrm{p}; t_\mathrm{f}) \varPhi (t_\mathrm{p}; t_\mathrm{f})$, and can be found by singular value decomposition (SVD) of $\varPhi (t_\mathrm{p}; t_\mathrm{f})$ given in Eq. (26). Here $E_{g}(t_\mathrm{f}) := \Vert g(t_\mathrm{f})\Vert _2^2=\int _{-\infty }^{\infty } |g(w,t_\mathrm{f})|^2 \text {d} w$. Since there is no time derivative of velocity appeared in our model equations, the velocity perturbations are slaved to the concentration perturbations [45]. Thus, we have chosen the square of the concentrations as a mathematical measure to identify the optimal amplification of perturbations.

However, this choice of flow quantity and norm is somewhat arbitrary, and there are other choices [35, 38, 46] which may be used in the linear framework. For instance, perturbations maximizing velocity perturbation, solute concentration, or total perturbation energy incorporating both velocity and concentration could be computed. Large linear transient growth could thus be associated with different perturbation structures other than those discussed in the present paper.