1 Introduction

The dynamic properties of aircraft wakes are critical to understand due to the danger they pose to downstream aircraft [26]. Aircraft wakes may be idealized by a pair of equal-strength counter-rotating vortices [28], which are subject to a range of complex behavior, including merging, rotation, instability, and breakdown. The dynamics and stability of vortex pairs have been reviewed in detail by [20].

Perturbations of tip vortex formation due to, for example, wing pitching and heaving or atmospheric turbulence can significantly alter the trajectory and evolution of the trailing vortex downstream [13, 14]. Simple harmonic motions of a lifting surface introduce curvature into the vortex system in the form of sinusoidal perturbations, akin to those studied by Widnall and colleagues [32, 33]. Bliss [3] used matched asymptotic analysis to show that vortex curvature causes a single vortex subject to a sinusoidal perturbation to rotate as a rigid body under its own self-induced velocity. Multiple curved vortices interacting with each other will further complicate the dynamics of the idealized vortices and their stability. Therefore, a clear understanding of the response of an isolated vortex or vortex pair system to self-induced or external perturbations is critical to the anticipation and prediction of the evolution of the wakes of aircraft and their control surfaces.

Stability analysis of vortex systems is also largely motivated by flow control, where researchers often seek to find and exploit the most unstable mode of the system to promote vortex dissipation and breakdown. The stability of vortex pairs and more general wake configurations has a rich history in the literature, often stemming from the well-known work of [7]. The Crow instability describes the exponential growth of long-wavelength sinusoidal perturbations along an initially straight counter-rotating vortex pair with equal circulation magnitudes, where the principal parameter is the wavelength of the perturbation. The Crow instability may be described as the balance of the self-induced rotation of the individual vortex and the mutual interaction of the vortices such that the induced rotations cancel, allowing the vortices to grow radially in a fixed plane (often 43–48\(^{\circ }\) from the vertical plane). Configurations that grow faster than the Crow instability are highly sought after, such as in [4], where a periodic array for vortex rings along the length of the vortex pair optimally excites the Crow instability. Bristol [5] introduced a new term to Crow’s equations for vortices with unequal strengths that accounted for the rotation of the vortex pair about its vorticity centroid, which allowed traditional modal approaches to be employed; the interested reader is referred to the systematic review of modal analysis techniques by [31] as related to the linear stability of vortex filament systems.

Crouch [6] used Floquet analysis to study a more general case of the Crow instability by analyzing the stability characteristics of a periodic 4-vortex system, which represented a tip vortex pair coupled to the vortex pair emanating from inboard flaps on a typical wing. The co-rotating flap vortices in the Crouch system were shown to cause vortex perturbations to grow at a rate that is five times faster than the typical Crow instability for a single counter-rotating vortex pair [6]. Rennich [25] compare a vortex filament model approach against Navier–Stokes simulations of 4-vortex configurations and found qualitative agreement between the two methods. A key finding of their study is how inboard flap vortices amplify perturbations along the outer tip vortices, which has been shown to be an effective strategy to accelerate perturbation growth. Fabre [11] further extended this work to 4-vortex system by implementing optimal perturbation analysis on counter-rotating configurations. A theme among the work of researchers who study 4-vortex systems is the pursuit of instability modes that grow faster than the Crow instability of an ordinary trailing vortex pair. The present work notes that a vortex pair in ground effect can be effectively modeled as a 4-vortex system with constraints, which invites the question as to whether the ground-image system would also permit instabilities that grow faster than the Crow instability.

Experimental work performed by [1] has focused on the influence of a wall on the evolution of a vortex pair subject to a long-wavelength perturbation. A key finding of their experimental campaign is that the wall appears to inhibit the development of the long-wavelength instability for their ‘vertical rings mode’ interaction. This mode is when the Crow instability has yet to develop into vortex rings prior to vortex–wall interaction, characterized by a generation of strong axial flow within the vortex cores. The topology of the flow field changes dramatically as secondary vorticity is generated at the wall, where an array of vortex rings are generated and rise up from the wall. Asselin [1] notes that this highly three-dimensional flow is a part of a broader class of vortex dynamics problems, where vorticity cancellation due to viscous wall interaction induces a pressure gradient and axial flow within the vortex core. Morris [21] performed experiments of an analogous problem of a straight vortex pair impinging upon a wavy wall. This work found a notable difference between a perturbed vortex impinging upon a flat wall, where the vortices do not exhibit strong axial flows during the wavy wall interaction.

In complement to these experiments, [8, 9] performed direct numerical simulations and an optimal perturbation analysis on the development of the long-wavelength instability of a vortex pair descending toward a wall. A principal finding of their analysis is the wavenumber shift to higher values than the unbounded Crow instability as the vortex pair approaches the wall. These simulations also identify various modes of the interaction, which are classified by the amplitude of the instability upon wall interaction. Stuart [30] investigated the transient growth associated with the secondary vortices that form as a result of the vortex pair interaction with a viscous wall and found that the optimal perturbation of 4-vortex system inhibits the rebound of primary the vortices from the wall surface.

The stability of vortices near wall boundaries is a canonical vortex dynamics problem relevant to several engineering contexts, such as vortex generators on wings, ship hull appendages, and tip leakage vortices in turbomachinery [1]. Benton  [2] investigated computationally the connected dynamics of a single vortex with its image in a planar boundary and found that the instability of a single vortex passing along an inviscid slip wall behaves in accordance with the linear theory of [7]. Simulations of a single vortex passing along a viscous no-slip also exhibited the Crow instability yet with a damped growth rate. Experiments by [23] examined the evolution of the Crow instability for a vortex pair and a single vortex above a flat plate, where the Crow instability was not observed for a single vortex above a viscous no-slip wall, as the roll up of boundary-layer vorticity dominated the interaction.

The present research effort explores the temporal evolution and stability characteristics of a counter-rotating trailing vortex pair in the presence of an inviscid planar surface. The mathematical model is formulated in Sect. 2, and details of the non-modal stability theory employed in the present work are provided in Sect. 3. The current analysis admits results for the limiting case of vortex pairs near a boundary where the viscous vortex-boundary-layer interaction is negligible, which is discussed in detail in Sect. 4. The inviscid interaction between the vortices and their image in the ground plane detailed in this section elucidates the near-wall behavior observed by [1] with respect to perturbation growth inhibition. Concluding remarks are then presented in Sect. 5.

Fig. 1
figure 1

Three-dimensional representation of the model for a perturbed vortex pair in ground effect. Two wavelengths are shown of an infinite sinusoidal perturbation that extends to \(\pm \infty \) in the \({\hat{e}}_x\)-direction. The ground effect constraints are held along the entire length of the vortices. The coordinate system local to an individual vortex is shown, where the perturbation directions \(\eta \) and \(\zeta \) are shown with the respective angular orientation \(\theta \) of the vortex

Full size image

2 Model formulation

A counter-rotating vortex pair in ground effect is modeled as a pair of vortex filaments of finite core radius \(a\) with their corresponding mirror image across an inviscid ground plane. The governing equations for a general \(n\)-vortex system are adapted from [6] and [11], who employed a vortex filament method to derive the equations of motion for linear vortex perturbations. Figure 1 illustrates the three-dimensional vortex system, where the centers of the physical vortices are separated by a distance \(b\) and lie at a height \(h\) above the ground plane. A sinusoidal perturbation of wavelength \(kb\) is prescribed along the entire infinite length of each vortex, where \(k\) is the perturbation wavenumber. In mathematical terms, the position vector of the nth vortex element with respect to the global coordinate system is

$$\begin{aligned} {\textbf{r}}_n(x,t) = x\,{\hat{e}}_x + \left( Y_n(t) + \eta _n(t)\text {e}^{\text {i}kx}\right) {\hat{e}}_y + \left( Z_n(t) + \zeta _n(t)\text {e}^{\text {i}kx}\right) {\hat{e}}_z, \end{aligned}$$
(1)

where \(Y_n\) and \(Z_n\) are the coordinates of the unperturbed vortex center position at time \(t\), and \(\eta _n\) and \(\zeta _n\) are the perturbation amplitude components with respect to the local coordinate system attached to vortex \(n\). The unperturbed vortices are oriented along the \({\hat{e}}_x\)-direction such that \(x\) is the variable of integration for the Biot–Savart law,

$$\begin{aligned} {\textbf{U}}_n = \sum _{m\ne n} \frac{\Gamma _m}{4\pi }\int _{-\infty }^{\infty }\frac{{\textbf{r}}_m-{\textbf{r}}_n}{\Vert {\textbf{r}}_m-{\textbf{r}}_n\Vert ^3}\frac{\partial {\textbf{r}}_m}{\partial x_m}\text {d}x_m, \end{aligned}$$
(2)

where \(n\) is understood to be the total number of vortices. The Biot–Savart law computes the induced velocity, \({\textbf{U}}_n\), by each vortex onto the others. Note that an alternative method is used to compute the effects of self-induction when \(m=n\), which is described ahead in this section.

The governing equations of the vortex system are determined by linearizing the Biot–Savart law to first order in the perturbation amplitudes \(\eta _n\) and \(\zeta _n\), where we assume the perturbations are small relative to the vortex separation distance \(b\). The mean vortex center positions are calculated from the kinematics of the temporally varying position vector and the leading-order expansion of the Biot–Savart law in terms of the perturbation amplitudes, and the perturbation amplitude equations are found at the next order. The characteristic length scale is the initial vortex separation \(b\), the time is non-dimensionalized as

$$\begin{aligned} \tau = \frac{\Gamma _0}{2\pi b^2}t, \end{aligned}$$
(3)

and all other non-dimensional quantities are denoted with a hat. Vortices of equal strength and core sizes are assumed in the present investigation, yet the method can permit more general vortex configurations. The resulting non-dimensional equations for the unperturbed vortex center positions are

figure b

and the non-dimensional equations for the perturbation amplitudes are organized in a similar fashion to [11]:

figure c

The right-hand sides of Eqs. 5a and 5b involve modified mutual induction functions \(\phi \) and \(\psi \), originally determined by [7], where the definitions can be found in Eq. 2.4 in [11]. The self-induced rotation, \({\bar{\omega }}\), is calculated in the same manner as [10], inspired by Saffman’s [27] analysis of Kelvin waves on a vortex filament.

Saffman [27, Chapter 12] analyzed oscillations of a vortex column to calculate the self-induced rotation frequency for small disturbances along a uniform rectilinear vortex filament by obtaining roots to the dispersion relation ([27, Eq. 11.3.4]) for Kelvin waves on a vortex filament. The self-induced rotation frequency, \({\bar{\omega }}\), is calculated by obtaining the lowest retrograde root of this dispersion relation for the bending mode and evaluating the frequency from Eq. 12.1.5 in [27, p. 232]. This method of calculating self-induction has been shown to not produce spurious instabilities at short wavelengths when compared to the original self-induction function of [7], as shown by [10]. Figure 2 illustrates the equivalence of the cutoff method and the results from the dispersion relation for low wavenumbers and small core sizes when substituted into Crow’s original equations [7, Eqs. 10–11]. It is understood that \(\alpha \) in Fig. 2 is the highest positive symmetric root from Crow’s [7] Eq. 11, which represents the exponential growth rate of a symmetric, sinusoidal, vortex perturbation.

The equations for the perturbation amplitudes is a linear, non-autonomous system. To consider cases near a wall boundary, we impose the constraints \({\hat{Y}}_1={\hat{Y}}_3\), \({\hat{Y}}_2={\hat{Y}}_4\), and \({\hat{Z}}_1=-{\hat{Z}}_3\), \({\hat{Z}}_2=-{\hat{Z}}_4\) for the vortex center positions. Likewise, \({\hat{\eta }}_1={\hat{\eta }}_3\), \({\hat{\eta }}_2={\hat{\eta }}_4\), and \({\hat{\zeta }}_1=-{\hat{\zeta }}_3\), \({\hat{\zeta }}_2=-{\hat{\zeta }}_4\) for the perturbation amplitudes. Because these constraints must hold for all time to represent the ground effect problem, the non-autonomous system may be reduced from an \(8 \times 8\) system in [11] to a modified \(4 \times 4\) system.

Fig. 2
figure 2

Spectrum of eigenvalues \(\alpha \) with respect to wavenumber employing Crow’s original theory with both a cutoff method and the dispersion relation for Kelvin waves on a Rankine vortex to calculate self-induction. There is agreement between both methods, as they have equivalent values in the low wavenumber limit

Full size image

The model vortex problem may be compactly expressed as a system of ordinary differential equations, where the matrix A depends on the time-dependent vortex center positions:

$$\begin{aligned} \frac{\text{ d }{{{\textbf {x}}}}}{\text{ d }\tau } = {{\textbf {A}}}(\tau ){{\textbf {x}}}, \qquad \qquad {{\textbf {x}}} = ({\hat{\eta }}_1, \ldots , {\hat{\eta }}_n, {\hat{\zeta }}_1,\ldots ,{\hat{\zeta }}_n)^T. \end{aligned}$$
(6)
Fig. 3
figure 3

Induced trajectory of the vortex pair in free space and approaching an inviscid wall. The wall requires a hyperbolic trajectory of the vortex pair, in accordance with classical theory [7]

Full size image

The governing equations are integrated numerically via a fourth-order Runge–Kutta scheme over a range of dimensionless wavenumbers \(kb\) to observe how the perturbations grow with respect to time. Numerical convergence was achieved for a time step of \(\Delta \tau = 0.005\) or smaller. The key parameters of interest in this problem are the initial dimensionless vortex height \(h/b\) and the fundamental time scale \(\tau \), which is the time it takes for an undisturbed vortex pair to translate a distance \(b\) in the \(-{\hat{e}}_z\) direction. Note that [19] has previously calculated the theoretical trajectory of a pair of point vortices, where the vortices descend on a hyperbolic trajectory away from each other along an inviscid wall. The true trajectory of the vortex pair descending toward a ground plane would rebound and persist above the wall due to viscous fluid effects, as explained by [16]. Figure 3 illustrates representative trajectories of vortex pairs with or without an inviscid ground constraint.

The small finite core size \({\hat{a}} = a/b = 0.098\) is chosen to validate the present method with Crow’s original theory for a free vortex pair. Initial starting heights are selected as \(h/b=5\) and \(h/b=10\) to observe the transient growth of perturbations with respect to distance away from the wall for the constrained cases, with a larger core size of \(a/b = 0.23\), based on the direct numerical simulations of [8].

3 Time-dependent stability and optimal perturbations

The stability of perturbations to time-dependent dynamical systems has been formalized by [12] in their work on generalized stability theory. Generalized stability theory is a non-modal approach that can capture asymptotic and finite-time behaviors of perturbations to time-dependent systems, where the non-autonomous operator may also be non-normal. In fact, it has been shown in [12] that non-normality of a system can lead to transient growth in both autonomous and non-autonomous systems such that the systems are asymptotically unstable. For the present work, the finite-time stability of infinitesimal perturbations is the primary focus, where we are not concerned with the asymptotic stability of the system.

3.1 Non-modal stability theory

Generalized stability theory seeks a fundamental solution matrix \(\Phi _{[\tau ,\tau _0]}\), or the propagator matrix, which maps the state of the system at \({\textbf{x}}(\tau _0)\) to \({\textbf{x}}(\tau )\), where \(\Phi _{[\tau ,\tau _0]}\) satisfies

$$\begin{aligned} \frac{\text {d}{\Phi _{[\tau ,\tau _0]}}}{\text {d}{\tau }} = {\textbf{A}}(\tau )\Phi _{[\tau ,\tau _0]}, \quad \Phi _{[\tau _0,\tau _0]} = {\textbf{I}}. \end{aligned}$$
(7)

This solution is analogous to the familiar matrix exponential solution for autonomous systems. The fundamental solution matrix maps the system at some time \(\tau _0\) to a future time \(\tau \) subject to an initial condition \({\textbf{x}}_0\) and characterizes the system,

$$\begin{aligned} {\textbf{x}} = \Phi _{\left[ \tau ,\tau _0\right] }{\textbf{x}}_0. \end{aligned}$$
(8)

The fundamental solution matrix is analogous to the matrix exponential for autonomous systems, where the solution would be expressed as

$$\begin{aligned} {\textbf{x}} = e^{{\textbf{A}}\tau }{\textbf{x}}_0. \end{aligned}$$
(9)

The asymptotic stability of an autonomous system is characterized by the eigenvalues of the matrix \({\textbf{A}}\), where the amplification of an initial perturbation at some time \(\tau \) is

$$\begin{aligned} \sigma (\tau ) = e^{\alpha \tau }, \end{aligned}$$
(10)

where \(\alpha \) is the largest real eigenvalue of \({\textbf{A}}\). For a non-autonomous system, stability is not quantified by the eigenvalues of the time-dependent operator \({\textbf{A}}\), as the eigenvalues may change at each instant in time. Instead, the optimal perturbation may be deduced by taking the singular value decomposition (SVD) of the fundamental solution matrix at some time instant of interest \(\tau \), where the amplification \(\sigma \) is the highest singular value. The right singular vector of the SVD corresponds to the initial condition that would lead to the optimal amplification. Further details of the non-modal stability theory are also described by [24].

Fig. 4
figure 4

Confirmation of non-modal (SVD) stability analysis with modal Crow method for a counter-rotating free vortex pair. The methods agree on the peak eigenvalue or its equivalent at \(\alpha _{{\max }} = 0.83\) at \(kb = 0.73\), where \(\alpha \) is the real component of the complex eigenvalue defined by Crow. \(\ln \left| \sigma \right| /\tau \) approximates \(\alpha \) to make a direct comparison between the Crow method and the present non-modal approach. The vortex core size is \({\hat{a}} = a/b = 0.098\)

Full size image

4 Results and discussion

The case of a free vortex pair studied by Crow [7] reduces to a linear autonomous system where the asymptotic behavior of initial perturbations are measured by the largest real eigenvalue \(\alpha \) of this system. To relate the optimal perturbation of a free vortex pair to the results from Crow, the maximum eigenvalue is approximated as

$$\begin{aligned} \alpha = \lim _{\tau \rightarrow \infty }\frac{\ln \left| \sigma (\tau )\right| }{\tau }, \end{aligned}$$
(11)

where \(\sigma \) is the highest singular value (amplification of an initial perturbation) obtained from the SVD of the fundamental solution matrix. Equation 11 is the same definition from Farrell and Ioannou [12] for the first Lyapunov exponent, which quantifies the asymptotic exponential rate of growth or decay of time-dependent systems. The amplification of a perturbation, \(\sigma \), for a non-autonomous system is of primary interest here, where we will use Eq. 11 to facilitate a comparison with Crow’s work, which deals with rates of growth.

To compare the eigenvalue approximation Eq. 11 from the SVD approach to Crow’s original work, we select \(a/b=0.098\) to match his parameters for the trailing vortex core radius of a real aircraft wake vortex. Figure 4 illustrates the agreement between both methods, where there is agreement in the maximally growing mode at \(kb=0.73\). However, it is necessary to integrate forward in time substantially to match the full wavenumber spectrum, roughly 80 non-dimensional time steps. The reason for the difference at initial times would be in part due to the non-normality of the system, where transient growth may not correspond with the asymptotic behavior defined by the eigenvalues. For finite-time stability analysis, this long time requirement is not a relevant concern, as we are interested in the maximally growing mode at a particular instant in time.

Two illustrative examples are provided to compare the evolution of the long-wavelength instability for a vortex pair in ground effect to an unconstrained vortex pair. The vortex core radius is chosen as \(a/b = 0.23\) in line with a subset of the parameters from [9], but it should be noted that the choice of vortex core size is not the primary interest in this study. Modifying the vortex core size (for the equal strength, equal size configurations studied here) will modify at what wavenumber the optimal perturbations occur yet retains the overall behavior of the system. With the goal of comparing the unconstrained free vortex pair to the constrained ground effect case, it is only necessary to be consistent in a choice of core size. We consider an initial height of \(h/b = 5\) and \(h/b=10\) and compare the evolution of the vortex systems with respect to an unconstrained free vortex pair.

Fig. 5
figure 5

Vortex trajectories and optimal wavenumber spectrum for an initial height \(h/b=10\) at \(\tau = 5.5\) (top row), \(\tau = 9.5\) (middle row), and \(\tau = 10.5\) (bottom row). The right column are the amplifications \(\sigma \) of an initial perturbation scaled by the corresponding maximum amplification of the unconstrained free vortex pair case \(\sigma _{\max }\). The vortex core size is \({\hat{a}} = 0.23\)

Full size image
Fig. 6
figure 6

Vortex trajectories and optimal wavenumber spectrum for an initial height \(h/b=5\) at \(\tau = 2.7\) (top row), \(\tau = 4.1\) (middle row), and \(\tau = 5.5\) (bottom row). The right column are the amplifications \(\sigma \) of an initial perturbations scaled by the corresponding maximum amplification of the unconstrained free vortex pair case \(\sigma _{\max }\). The vortex core size is \({\hat{a}} = 0.23\)

Full size image

Figure 5 shows the growth of an initial perturbation across all wavenumbers considered for an initial height of \(h/b=10\) above the wall compared to the free vortex pair case. Prior to any significant wall interaction, the spectrum of perturbation growth for the ground-constrained case coincides closely with the free vortex pair case. Normalization of the amplification \(\sigma \) with respect to the maximum amplification of the free vortex pair highlights the difference between the two configurations. The normalization with respect to the maximum amplification factor of the free vortex pair occurs at each instant considered, which accentuates the relative difference between the maximum amplification of the free vortex pair and the ground-constrained pair. As the trajectories diverge, the amplification of the ground-constrained vortex pair perturbations is less than half of the free vortex pair. In a real viscous fluid, the Crow instability would dominate the evolution of a vortex pair at an initial height of \(h/b=10\) and convert into vortex rings before the wall interaction [1]. This phenomenon does not occur in our model, as vortex ring creation is not an inviscid process. The dynamics of the inviscid model system agrees with experiments in the context of perturbation growth inhibition, despite the inability to model the creation of vortex rings for such large initial heights. For long time intervals, the perturbation amplitudes grow exponentially as if they were a free vortex pair if they were to follow along their inviscid idealized trajectory. Kliment [18] investigated experimentally vortex pairs in ground effect and noted that viscous effects may be negligible for heights above the wall 0.5–1.5 times the vortex span \(b\), and that the wall effect may be modeled with image vortices. This spacing constraint does not imply that viscous effects are negligible for all times that the vortices are within this spatial region. As the vortices persist above the ground in a real viscous fluid, the boundary layer will grow until distinct secondary vortices are generated and induce a ‘rebound’ in the primary vortex pair causing the vortices to deviate from their idealized hyperbolic path [16]. The present analysis is constrained to trajectories before vortex rebound would typically occur, which occurs approximately 0.5–2.0b spatial units laterally from their original locations depending on the circulation-based Reynolds number of the vortices [22].

Fig. 7
figure 7

Three-dimensional representation of the optimal perturbations for free vortex pair (a) and a vortex pair subject to a ground constraint (b). The perturbations at \(\tau = 2.7\), \(4.1\), and \(5.5\) correspond to instability of interest in Fig. 6. All spatial positions, wavelengths, amplitudes, and angular orientations of the vortices are obtained by integrating Eqs. 5a and 5b, with the optimal perturbation prescribed as the initial condition up until the corresponding instant of interest

Full size image

Figure 6 shows the same qualitative behavior as Fig. 5 for a vortex pair that starts closer to the wall at an initial height \(h/b=5\). The amplification curves are in alignment at early times before the trajectories diverge followed by clear perturbation growth inhibition. It should be noted that the time and spatial dependence of the perturbation spectrum for the free vortex pair case is due to transient growth at early times, where the system approaches its steady-state solution that would be characterized by its eigenvalues. Beyond the clear reduction of perturbation growth in the ground-constrained case, it is evident that the wavenumber of the optimal perturbation shifts to higher values as the vortex pair approaches the wall. This observation indicates that wavenumbers higher than the Crow instability (optimal wavenumber for the unconstrained case) should be selected to maximize the growth of the vortex pair at the wall.

Figure 7 further elucidates the structure of the optimal perturbations for a free vortex pair versus a vortex pair subject to a ground constraint. The rotation of the plane of the most unstable initial perturbation at a given instant in time to a nearly horizontal position in Fig. 7 is a striking difference between the final state of the free vortex pair, which maintains its angular position of approximately \(47^{\circ }\). The final amplitudes of the free and ground-constrained configurations are also notably different, with the latter vortex pair possessing an amplitude that is roughly \(40\%\) of the free vortex pair. In a real viscous fluid, the free vortex pair would begin to undergo a reconnection process to form vortex rings and secondary vorticity generation at the ground plane would initiate vortex rebound [1].

The optimal initial condition and wavenumber given by the SVD approach can be directly applied to the system and integrated in time to form a time history of the optimal perturbations. The instantaneous energy gain of the system \(G\) is defined by

$$\begin{aligned} G = \frac{\left\| {\textbf{x}}\right\| }{\left\| {\textbf{x}}_0\right\| } \end{aligned}$$
(12)

and represents the amplification of an initial perturbation, where \({\textbf{x}}\) are the perturbation amplitude components, and \({\textbf{x}}_0\) is a vector describing an arbitrary initial perturbation of unit magnitude. Note that this definition of \(G\) is different than that of \(\sigma \), as \(\sigma \) is the amplification of a perturbation provided the optimal initial condition is considered for a particular wavenumber at some point in time. It is important to recall that the search for the optimal perturbation seeks the maximally amplified initial condition at each instant in time per wavenumber. This implies that the wavenumber of maximal perturbation growth can be different in time, as shown in Figs. 5 and 6, which can be expected due to the non-autonomous nature of the governing equations.

Integration of the governing equations at particular wavenumbers and initial condition further reveals the time-dependent dynamics of the system. Figure 8a illustrates the growth curves for three different conditions on the optimal perturbation. Specifically, we plot the perturbation energy gain for the optimal perturbation at \(\tau = 2.5\), \(\tau = 5.0\), and \(\tau = 5.5\). The wavenumbers of the associated perturbations are \(kb=0.87\), \(kb=0.91\), and \(kb=0.98\), respectively. The above times are selected to showcase how the optimal perturbation is time-dependent, i.e., optimizing growth for a particular instant in time will not guarantee maximal growth at a later time. Recall that the initial condition is obtained from the right singular vector of the singular value decomposition of the fundamental solution matrix at the time of interest \(\tau \) corresponding to the maximally growing wavenumber.

Fig. 8
figure 8

Energy gain (a) and instability plane angle (b) as a function of time for three wavenumbers and their optimal initial condition for select times of interest. \(kb=0.87\) maximizes the gain at \(\tau = 2.5\), \(kb=0.91\) maximizes energy gain at \(\tau =5.0\), and \(kb=0.98\) maximizes energy gain at \(\tau =5.5\). Vortex core size is \({\hat{a}}=0.23\) for each case

Full size image

For optimal perturbations considered before \(\tau = 5.5\), a decrease in the perturbation amplitudes after the wall interaction is observed, which has also been shown in the experiments of [1]. This effect has been attributed to the rotation of the plane of which the instability exists as the vortices are swept away from each other along the wall. Recall that the Crow instability balances the self-induced rotation of the vortex, the induced strain field from the other vortex, and the mutual interaction between the perturbed vortices, allowing the perturbations to grow radially in a fixed plane [20]. The rotations of the plane of instability are shown in Fig. 8, which all illustrate the rotation induced by the wall. The optimal perturbation at \(\tau = 5.5\) shows a decrease in the rate of amplification, yet is maximally excited compared to the other cases.

Fig. 9
figure 9

Comparison on the effect of initial condition for the wavenumbers \(kb = 0.82\) (a) and \(kb = 0.98\) (b). \(x_{0}^{\text {opt}}\) indicates the optimal initial condition of unit length from the right singular vector of the SVD of propagator matrix, where in this case the chosen times of optimization are \(\tau = 2.5\) (a) and \(\tau = 5.5\) (b). The initial conditions \(x_{0_2}\), \(x_{0_3}\), and \(x_{0_4}\) are randomly generated and scaled to be of unit length. The random initial conditions for the perturbations showcase how for the same wavenumber, optimal energy gain is highly dependent on the initial conditions. The initial height above the wall is \(h/b=5\) with a vortex core size of \({\hat{a}}=0.23\) for each case

Full size image

Beyond choosing the optimal wavenumber at a particular instant in time, it is also necessary to consider the effect of the corresponding initial condition \({\textbf{x}}_0\). The optimal initial condition for a given wavenumber is denoted \({\textbf{x}}_0^{\text {opt}}\), which corresponds to the right singular vector of the singular value decomposition of the fundamental solution matrix at an instant in time \(\tau \). The right singular vector is of unit magnitude and describes the geometric orientation of the initial perturbation. To compare \({\textbf{x}}_0^{\text {opt}}\) to other initial conditions, we randomly generate three other initial conditions of unit magnitude and directly integrate the systems. A similar illustration of this idea is presented by [24] in the context of non-modal growth of perturbations in density-driven convection in porous media. Figure 9 exemplifies that, for a given wavenumber, the geometric initial condition can have significant impact on subsequent perturbation growth. For the configurations studied in this work, all optimal perturbations were geometrically symmetric. This observation is not surprising, as the Crow instability [7] is a symmetric perturbation, where anti-symmetric modes are notably less amplified. For the 4-vortex systems studied by [6] and [11], symmetric configurations were also determined to be maximally growing.

5 Conclusion

The time-dependent characteristics of a perturbed equal strength and equal size counter-rotating vortex pair interacting with an inviscid wall have been studied as a special case of the equations and methods of [6] and [11]. The non-modal approach of [12] enables the finite-time stability of the temporally varying vortex system to be studied, where we draw connection to the limiting case of the time-independent and non-periodic free vortex pair system by approximating the eigenvalues at long times via Eq. 11. To elucidate the effect of the wall on perturbation growth, the general 4-vortex system is recast as a 2-vortex system with ground constraints. The ground-constrained system is compared against an unconstrained free vortex pair to further highlight the effect of the wall on the development of the cooperative long-wavelength instability of the vortex pair.

There are important differences between the present inviscid analysis and the experimental work of [1] as well as the computational work by [8, 9]. In a real viscous fluid, the vortex pair will undergo viscous dissipation and lose circulation strength during its descent, diminishing the magnitude of the interaction between the vortices. Secondary vorticity generation at the wall boundary will emerge and lead to vortex rebound from the surface. It is understood in the present work that we emphasize the inviscid interaction between the vortex pair and the wall and limit our analysis to where boundary-layer effects can be neglected.

For vortex pairs starting \(h/b=10\) and \(h/b=5\) above the wall, the vortex pair does not significantly deviate from their unconstrained counterparts trajectory until it is approximately a distance of 1.5–2b above the wall. For the initial vortex heights considered, the corresponding amplifications of an initial perturbation across all wavenumbers considered follow the same trend, where we observe a 40–60% reduction in perturbation amplification as the vortices approach the inviscid wall. We also observe that the optimal wavenumber yielding the highest perturbation amplification shifts to a higher wavenumber compared to the free vortex pair case, indicating that the wall influences the conditions for the maximum perturbation growth rate. A similar phenomenon is observed in [9], where they also identify a transient shift in optimal wavenumbers to higher values. The suppression of the long-wavelength instability relative to a free vortex pair is further supported by the experiments of [1], illustrating that perturbation inhibition may at least be partially explained without viscous fluid arguments. As the vortices diverge from their unconstrained trajectories, the magnitude of the mutual induction between vortices is reduced. This reduction can be simply explained by the fact that the mutual interaction functions in the governing equations scale on the inverse of their separation distance.

Direct integration of the model vortex system allows for further investigation on the effect of initial conditions of the present configurations. The present non-modal approach yields both the optimal wavenumber for perturbation gain and the initial condition that leads to maximum growth for any wavenumber at any instant in time. The optimal initial condition for all cases considered corresponded to a symmetric initial condition, where the vortices are originally oriented at 40–50\(^{\circ }\) to the horizontal. The effect of choosing the optimal initial condition is illustrated in Fig. 9a and b, where random initial conditions of unit magnitude were generated for comparison. The initial orientation and symmetry of the vortices influence the resulting magnitude of perturbation growth significantly, yet do not alter the general trend of the energy gain curves. We observe that upon wall interaction, a decrease in perturbation growth rate is observed, with some parameters yielding a temporary decrease in total perturbation amplitude during wall interaction. Choosing the optimal wavenumber at a given instant in time also notably impacts perturbation evolution, where selecting an initial wavenumber higher than the free vortex pair Crow instability results in larger perturbation amplification at the wall. From a control strategy perspective, it would be critical to not only excite a particular wavenumber of interest, but also initialize the perturbation such that it best realizes the corresponding optimal initial condition. It is likely that exciting the optimal perturbation will accelerate vortex breakdown and dissipation, which, for example, is desirable for reducing separation distances at airports during takeoff and landing operations [15, 29]. Despite the limitations of the present inviscid analysis, the above findings can better inform the study of more realistic configurations to study the time-dependent nature of a vortex pair in ground effect.