Wind-Driven Waves on the Air-Water Interface

Abstract

An ocean swell refers to a train of periodic or nearly periodic waves. The wave train can propagate on the free surface of a body of water over very long distances. A great deal of the current study in the dynamics of water waves is focused on ocean swells. These swells are typically created initially in the neighborhood of an ocean storm, and then the swell propagates away from the storm in all directions. We consider a different kind of wave, called seas, which are created by and driven entirely by wind. These waves typically have no periodicity, and can rise and fall with changes in the wind. Specifically, this is a two-fluid problem, with air above a moveable interface, and water below it. We focus on the local dynamics at the air-water interface. Various properties at this locality have implications on the waves as a whole, such as pressure differentials and velocity profiles. The following analysis provides insight into the dynamics of seas, and some of the features of these intriguing waves, including a process known as white-capping.

1. Introduction

In this paper, a “wave” refers to a coherent, localized region in which the fluid flow is irregular. The localized region remains irregular as the fluid flow continues. Please note that this definition of a wave does not suggest any sort of periodicity, as one might see in a swell. Section 2.1 describes the background flow, on the air-water interface before any waves appear. This is the simplest situation, where the interface is a completely flat surface. The wind blows in the air above the interface, and if the water below is moving, it is purely horizontal, due to the flatness of the interface. In addition to the purely horizontal flow, other features in the flow can appear because of small local regions near the interface. In these local regions, several properties can create an irregular flow, and then the irregular flow can create larger, surrounding regions in which the fluid flow is also irregular. Some such local properties include the irregularity in the fluid itself, the irregularity in the pressure of the fluid, and the flow having some vertical component in this small region. The sources of these irregularities are not discussed in depth, but one source can be the diffusing vorticity. The tangential shear stress is continuous across the interface, so the z-derivative of $U\left(\zeta \right(z,t\left)\right)$ must be discontinuous at the interface (see (9)). Thus, vorticity is generated at the interface. This vorticity diffuses into the fluid above the interface and also into the fluid below. This diffusing vorticity can be a source of irregular fluid motion above the interface, or below the interface, or both.

One of the foundational bodies of work in the field of propagating water waves came from J. W. Miles (1957) [1]. Miles viewed the problem as a stability problem, in which waves on the upper surface of the water can grow by drawing energy from a fixed horizontal shear flow in the air above that surface. In contrast, our problem is intrinsically time-dependent: the horizontal shear flows, in the air above the interface and in the water below it, each evolve slowly in time as vorticity (created at the interface due to viscous effects) diffuses away from the interface. An important discovery of Miles was that his problem contains a “critical layer”, where the horizontal speed of the fluid matches the propagation speed of the waves on the interface. There is no fixed “critical layer” in the sense of Miles in our model, because the shear flow in each fluid, as well as the propagating waves, evolve in time. Also, in order to observe the effects of a critical layer in the sense of Miles, there must be some connection between the speed of the waves on the interface and the speed of the air in the critical layer. When working with seas instead of swell, this connection cannot be made.

Another significant advancement in the field came from B.D. Dore (1977) [2]. Dore showed that most of the dissipation of the travelling wave comes from the air. Before, most work assumed a vacuum above the water, since it was believed that the dynamic viscosity of air is so low that its effects are negligible. A comprehensive overview of the current progress in understanding wind-driven waves can be found in Chapter 2 of “Wave modelling—The state of the art” (2007) [3], while most of the other chapters focus on some aspect of the dynamics of swells. Dore’s work introduces the effects of viscosity in the air, which is one step closer to the model we propose, but Dore is considering the dynamics of swell. Our work uses the viscosity of both the water and the air.

Section 2.2 explores the linearized equations for the perturbed motion about the background flow, and the resulting corrections to the quantities introduced in Section 2.1. Although small, these corrections are important, and together they define the set of linearized equations for the fluid flow. A few physical results arise from the linearized equations. One clear example is the “moveable interface”, first mentioned in the first paragraph of Section 2.2. The idea is that in the presence of fluid motion, including the linearized motion, the interface is allowed to drift within a small range of where is it expected to be. Thus, the moveable interface is the exact location of where the air and the water meet at any time and at any x-location. Let $\eta (x,t)$ denote the height of the movable interface at horizontal location (x) and time (t). Given a fixed x-location, $\eta (x,t)$ evolves over time. At any one point in time, the highest point of the wave is known as the wave crest. If the vertical velocity of the fluid is positive and strongest near the crest, then the height of the crest will continue to grow. As it grows, the horizontal wind (which increases in speed as one moves higher) drives the crest to rise even higher. At some point, the speed of the wind is so much higher than that of the wave that the entire crest is torn off. This process is known as white-capping, and has not been well understood in the past: “white-capping, or dissipation in deep water, is the third basic physical process that governs the evolution of wind waves in the open oceans. It is the least understood part of wave evolution” [3]. Section 4 provides some intuition into how the process of white-capping follows from the analysis of the equations given.

A major difference between all of the mentioned models and ours is that other works consider the types of waves that appear in swells while we consider the waves within seas. The lack of periodicity in seas implies a lack of many quantities typically analysed through experimentation, such as wavelength, phase speed, etc. Thus, experimentally measuring the seas is a much more difficult task. In addition, much of the work that has been done to analyze the dynamics of seas numerically or through simulation uses models with turbulent flow, a feature that is absent in our model. See Section 2.1 of [3] for more details.

2. Identify the Equations That Define the Linearized Problem

2.1. The Background Flow

We begin with a shear flow, inspired by “Stokes’ first problem” for viscous fluids. A heavier fluid (such as water) lies beneath a lighter fluid (such as air) under the force of a constant gravitational field. The upper fluid extends vertically up to $+\infty$ , while the lower fluid extends vertically down to $-\infty$ . There are no horizontal boundaries, and the flow is two-dimensional. Both fluids are viscous and incompressible. To make the analogy with Stokes’ first problem, assume that for $t<0$, the (horizontal) interface between the two fluids admits no tangential stress, and that well above the interface, the upper fluid moves horizontally with speed $U_{\infty }>0$, while well below the interface, the lower fluid moves horizontally with speed $U_{-}<U_{\infty }$. At the interface itself, the fluid moves horizontally at speed $U\left(0\right)$, with $U_{-}<U\left(0\right)<U_{\infty }$. At $t=0$, the interface suddenly admits tangential stress, and a boundary layer grows in each fluid, as originally found by Stokes.

In each fluid, there is no vertical velocity in the background flow, and no x-variation in the horizontal velocity, so the Navier-Stokes equations in each fluid reduce to

$$\begin{array}{cc}\hfill {\partial }_{t}U& =\nu {\partial }_z^{2}U\hfill \end{array}$$

(1a)

$$\begin{array}{cc}\hfill W& =0\hfill \end{array}$$

(1b)

$$\begin{array}{cc}\hfill \frac{1}{\rho }{\partial }_{z}P& =-g\hfill \end{array}$$

(1c)

where in each fluid, $\nu$ denotes its constant kinematic viscosity, $\rho$ denotes its constant density, $\left\{U\right(z,t),W\}$ are its two components of velocity, and $P\left(z\right)$ is its pressure. The two fluids can have different values of $\nu$, which we denote ${\nu }_u$ in the upper fluid and ${\nu }_l$ in the lower fluid. In addition, the two fluids can have different densities, ${\rho }_u$ and ${\rho }_l$, with

$$0<{\rho }_u\le {\rho }_l$$

(2)

as required for static stability. Equation (1a) admits self-similar solutions in terms of the variables

$$\zeta =\frac{z}{\sqrt{4{\nu }_{u}t}}\ge 0,{\zeta }_l=\frac{z}{\sqrt{4{\nu }_{l}t}}\le 0,$$

(3)

in the upper fluid $(z\ge 0)$ and lower fluid $(z\le 0)$, respectively. Then (1a) reduces to

$$-2\zeta U_u^{\prime }\left(\zeta \right)=U_u^{\prime \prime }\left(\zeta \right),-2{\zeta }_l{U}_l^{\prime }\left({\zeta }_l\right)=U_l^{\prime \prime }\left({\zeta }_l\right)$$

(4)

in each fluid. These ordinary differential equations can be integrated:

$$\begin{array}{c}{U}_u\left(\zeta \right)=c_1{\int }_0^{\zeta }{e}^{-{\sigma }^2}d\sigma +c_2,\zeta >0,\end{array}$$

(5a)

$$\begin{array}{c}{U}_l\left({\zeta }_l\right)=-c_3{\int }_{{\zeta }_l}^0{e}^{-{\sigma }^2}d\sigma +c_2,{\zeta }_l<0,\end{array}$$

(5b)

where $\{c1,c2,c3\}$ are positive constants, to be determined. The two integrals in (5a) and (5b) are each proportional to the error function. The parameter $c_2=U\left(0\right)$ defines the horizontal velocity at the interface, so $U_u\left(\zeta \right)$, the horizontal velocity in the upper fluid is faster than $U\left(0\right)$, while $U_l\left({\zeta }_l\right)$, the horizontal velocity in the lower fluid, moves slower than $U\left(0\right)$. In fact, $\left(U_l\left({\zeta }_l\right)\right)$ could even flow in the opposite direction from $U_u\left(\zeta \right)$. The formulae in (5a) and (5b) are well known. For example, see Carslaw & Jaeger (1959) [4].

$$\mathrm{Fact}\text{:}{\int }_0^{\infty }{e}^{-{\sigma }^2}d\sigma ={\int }_{-\infty }^0{e}^{-{\sigma }^2}d\sigma =\frac{\sqrt{\pi }}{2}$$

(6)

Let $U_{\infty }$ denote the maximum horizontal velocity, which is attained as $\zeta \to +\infty$ , and let $U_{-}$ denote the minimum horizontal velocity, which is attained as ${\zeta }_l\to -\infty$. (5a) and (5b) guarantee that $U\left(\zeta \right)$ is a monotonically increasing function of $\zeta$, and $U_u\left(0\right)=U_l\left(0\right)=c_2$ is an intermediate speed.

Combining these bounds with the results in (5a) and (5b) leads to

$$U_u\left(\zeta \right)=U_{\infty }-c_1[\frac{\sqrt{\pi }}{2}-{\int }_0^{\zeta }{e}^{-{\sigma }^2}d\sigma ].U_l\left({\zeta }_l\right)=U_{-}+c_3[\frac{\sqrt{\pi }}{2}-{\int }_{{\zeta }_l}^0{e}^{-{\sigma }^2}d\sigma ]$$

(7)

The overall range of horizontal velocities is

$$U_{\infty }-U_{-}=(c_1+c_3)\frac{\sqrt{\pi }}{2}$$

(8)

The two formulae in (5a) and (5b) guarantee that the horizontal velocity is continuous at $\zeta =0$ (equivalently, at $z=0$). In addition, physics requires that the tangential shear stress be continuous across the interface:

$${\rho }_u{\nu }_u·{\partial }_z{U}_u\left(0\right)={\rho }_l{\nu }_l·{\partial }_z{U}_l\left(0\right)$$

(9)

where ${\rho }_{*}$ denotes the density of each fluid, and ${\nu }_{*}$ denotes the kinematic viscosity of the same fluid. Unless ${\rho }_u{\nu }_u={\rho }_l{\nu }_l$, ${\partial }_z{U}_{*}\left(\zeta \right)$ is discontinuous at $\zeta =0$.

$$\mathrm{Define}\text{:}0<R:=\frac{{\rho }_u}{{\rho }_l}\le 1,0<\sqrt{V}:=\sqrt{\frac{{\nu }_u}{{\nu }_l}}$$

(10)

The parameter R defines the ratio of the densities of the two fluids, while V defines the ratio of the two kinematic viscosities. In the most common situation, with air above water, at a temperature in the range of 10–20 ${}^{\circ }$C, $R\approx 0.001$, while $3<\sqrt{V}<4$.

An important practical consequence of this information is obtained by making use of the information in (9), (3) and (5a), (5b). The final result can be written as

$$c_3=R\sqrt{V}·c_1$$

(11)

The parameter $c_1$ in (5a) is an arbitrary scaling parameter. Once $c_1$ is fixed, then $c_3$ in (5b) is fixed by (11).

For the common situation with air above the interface and water below the interface, (11) shows that the range of horizontal velocities in the air is much larger than the range of velocities in the water. A diagram of the horizontal velocities in the upper and lower fluids is shown in Figure 1.

An example: In a wind-tunnel tank, the air is blown through the wind tunnel by a fan at one end, and it exits the tunnel at the other end. Near the interface, the air can drag a thin layer of water along, because the air and the water are both viscous fluids. As the water in this thin layer reaches the other end of the tank, it has nowhere for it to go, so there must be a return flow of the water, but not of the air. The diagram in Figure 1 suggests that almost all of the water column travels at one speed, so the return flow of water in the tank apparently should use most of the water column.

Figure 1 is correct in a coordinate system that is moving to the right with a fixed, horizontal speed of $U\left(0\right)$, the horizontal velocity at the interface. The arrows in the figure show that all of the fluid above the interface is travelling to the right faster than $U\left(0\right)$, while all of the fluid below the interface is moving to the right slower than $U\left(0\right)$. In a coordinate system moving to the right with speed $U\left(0\right)$, the fluid below the interface appears to be moving to the left, but that appearance occurs because the coordinate system itself is moving to the right with speed $U\left(0\right)$.

In Section 2.4 and in subsequent sections, equations such as (29a) and (29b) show that the vertical velocity in the upper fluid is affected by integrals over the fluid above the interface, while the vertical velocity in the lower fluid is affected by integrals over the fluid below the interface. Each of these integrals contains a factor of the form ${\partial }_s{U}_u\left(\zeta (s,t)\right)$. Figure 1 suggests that these integrals are significant near the interface, but they become less significant as one moves away from the interface, because ${\partial }_{z}U\to 0$ as one moves away from the interface, either up or down.

In addition to the horizontal flow, (1c) shows that the background pressure distribution is hydrostatic in each fluid, so

$$\begin{array}{c}{P}_u\left(z\right)=P_0-{\rho }_{u}gz,z>0,\end{array}$$

(12a)

$$\begin{array}{c}{P}_l\left(z\right)=P_0-{\rho }_{l}gz,z<0\end{array}$$

(12b)

Equations (1a, 1b, 1c), (3), (5a), (5b), (7), (10) and (11) completely specify the (background) horizontal velocity field, while (12a) and (12b) describes the background pressure field in each fluid. The analysis from here on focuses on linearized perturbations about this time-dependent background flow.

The unknowns in the full problem, in each of the two fluids, are:

$$\begin{array}{cccc}\hfill \mathit{Pressure}-& \tilde{P}(x,z,t;ϵ)=P_0-\rho gz\hfill & & +ϵ·p(x,z,t)+O\left({ϵ}^2\right);\hfill \\ \hfill \mathit{Horizontal}\mathit{velocity}-& \tilde{u}(x,z,t;ϵ)=U\left(\zeta (z,t)\right)\hfill & & +ϵ·u(x,z,t)+O\left({ϵ}^2\right);\hfill \\ \hfill \mathit{Vertical}\mathit{velocity}-& \tilde{w}(x,z,t;ϵ)=\hfill & & +ϵ·w(x,z,t)+O\left({ϵ}^2\right);\hfill \\ \hfill \mathit{Height}\mathit{of}\mathit{the}\mathit{interface}-& \tilde{\eta }(x,z,t;ϵ)=\hfill & & +ϵ·\eta (x,t)+O\left({ϵ}^2\right);\hfill \end{array}$$

(13)

with subscripts ${(·)}_u$ and ${(·)}_l$ to denote the upper and lower fluids as needed, and with $ϵ$, a small dimensionless parameter that justifies the linearization.

2.2. Guide to the Linearized Problem

The linearized equations for the perturbed motion about the background flow in Section 2.1 are given in (14a), (14b) and (14c). Vishal Vasan (private communication) makes the important observation that $u(x,z,t)$, the perturbation of the horizontal velocity, plays only a secondary role in the linearized equations. Therefore the unknown functions needed to describe the linearized motion in each of the two fluids, above and below the moving interface, are the (perturbation) vertical velocity, $w(x,z,t)$, the (perturbation) pressure, $p(x,z,t)$, and the local height of the moving interface, $\eta (x,t)$. The equations can be somewhat simplified by representing these unknown functions in terms of their Fourier transforms in the x-direction. In what follows, the Fourier transform of the pressure in the upper fluid is ${\widehat{p}}_u(k,z,t)·e^{-\left|k\right|z}$, etc.

The primary objective in Section 2.3 is to write the pressure in the upper and lower fluid regions explicitly in terms of integrals over the vertical velocity, plus the boundary value of each pressure at the air-water interface. Then the entire linearized problem can be stated in terms of four evolution equations (in time), plus compatibility conditions at the interface. The variables whose evolution in time is given explicitly are as follows.

• ${\widehat{w}}_u(k,z,t)$ is the Fourier transform (in x) of the vertical velocity in the upper fluid.
  Its evolution equation, given in (29a), involves ${\widehat{p}}_u(k,0,t)$, the Fourier transform of the time-dependent pressure in the upper fluid at the interface. A procedure for solving this equation is given in Section 3.1.
• ${\widehat{w}}_l(k,z,t)$ is the Fourier transform of the vertical velocity in the lower fluid.
  Its evolution equation, given in (29b), involves ${\widehat{p}}_l(k,0,t)$, the time-dependent pressure in the lower fluid at the interface. A procedure for solving this equation is given in Section 3.2.
• $\widehat{\eta }(k,t)$ is the Fourier transform of the vertical displacement of the moving interface. Its evolution equation is given in (32).
• $\widehat{\overline{w}}(k,t)$ is the Fourier transform of the vertical velocity at the interface. It is the limiting value of ${\widehat{w}}_u(k,z,t)$ and of ${\widehat{w}}_l(k,z,t)$ as $z\to 0$ in each fluid. $\widehat{\overline{w}}(k,t)$ plays an important role in Section 2.5.

Comment 1.

Our assumption that each fluid is incompressible indicates that ${\widehat{w}}_u(k,z,t)$ and ${\widehat{w}}_l(k,z,t)$ should have the same limit as $z\to 0$. If they have different limiting values as $z\to 0$, then a connected region that includes $z=0$ either develops a vacuum (containing neither fluid) or it contains both fluids simultaneously. Either option violates our assumption of incompressibility. Even so, incompressibility by itself does not guarantee that ${\widehat{w}}_{*}(k,z,t)$ is continuous at the interface. In Section 2.5, Equation (34) provides a necessary condition for ${\widehat{w}}_{*}(k,z,t)$ to be continuous at the interface.

In addition, other variables do not evolve on their own, but evolve in time because they participate in algebraic relations that involve one or more of the four variables listed above. These algebraic relations describe the dynamics at the interface. The other variables are:

• ${\widehat{p}}_u(k,0,t)$ and ${\widehat{p}}_l(k,0,t)$ represent the limiting pressures in the upper and lower fluids, respectively, as $z\to 0$. We do not require that these two pressures coincide.
• Let ${\partial }_z\tilde{w}(k,z,t)$ denote the rate of change of the vertical velocity as a function of height (z) above the fixed interface. If we want to find the rate of change of the vertical velocity at the fixed interface, the symbol changes to ${\partial }_z\tilde{w}{(k,z,t)|}_{z=0}$. As given in (39c), the z-derivative of the vertical velocity at the interface is the common limit of the z-derivative in the upper and in the lower fluids, as $z\to 0$.
• ${\partial }_z^2{\widehat{w}}_{*}(k,z,t)$ must satisfy a jump condition across the interface, which is given in Section 2.5.

2.3. Find the (Linearized) Pressure in Each Fluid

In each fluid, the Navier-Stokes equations in two dimensions, linearized about the background flow in (1a, 1b, 1c), (12a,12b), (13), are the following:

$$\begin{array}{c}{\partial }_{t}u+U\left(\zeta (z,t)\right)·{\partial }_{x}u+w·{\partial }_{z}U\left(\zeta (z,t)\right)+\frac{1}{\rho }{\partial }_{x}p=\nu ·{\nabla }^{2}u,\end{array}$$

(14a)

$$\begin{array}{c}{\partial }_{t}w+U\left(\zeta (z,t)\right)·{\partial }_{x}w+\frac{1}{\rho }{\partial }_{z}p=\nu ·{\nabla }^{2}w,\end{array}$$

(14b)

$$\begin{array}{c}{\partial }_{x}u+{\partial }_{z}w=0.\end{array}$$

(14c)

As noted in Section 2.1, the background flow is purely horizontal, and the background pressure has no connection to the background flow. Linearizing the background flow introduces a perturbation pressure, $p(x,z,t)$, and a perturbation vertical flow, $w(x,z,t)$. This vertical flow is much weaker than the background flow, but it is essential for a well-known physical application, discussed in Section 4, Discussion, at the end of this paper.

The perturbation pressure, $p(x,z,t)$ is defined implicitly by (14a,14b,14c). To find p explicitly, compute {${\partial }_x$ (14a) + ${\partial }_z$ (14b)}, and use (14c) to eliminate unnecessary terms. The result is

$$\begin{array}{c}{\nabla }^{2}p(x,z,t)=-2\rho ·{\partial }_{z}U\left(\zeta (z,t)\right)·{\partial }_{x}w(x,z,t)\end{array}$$

(15)

in each fluid. This is a Poisson equation for $p(x,z,t)$ in a 2-D, semi-infinite domain in each fluid. Vishal Vasan (private communication) points out that, according to (15), $p(x,z,t)$ depends on $w(x,z,t)$, but not on $u(x,z,t)$. In addition, (14b) involves $\{p,w\}$, but not $u(x,z,t)$. Therefore (14b) and (15) define a closed system of equations for two unknowns, $\{p,w\}$. After $w(x,z,t)$ has been obtained by solving those two equations, one constructs $u(x,z,t)$ from (14c).

Comment 2.

This suggests that $u(x,z,t)$ has no effect on the dynamics of the linearized flow. This is nearly correct, but (14c) asserts that at any order in powers of $ϵ$, there always is a $u(x,z,t)$ that accompanies $w(x,z,t)$. In Section 2.5, below, this accompanying $u(x,z,t)$ plays a significant role in the boundary conditions that require the continuity of the normal and tangential stresses across the interface.

It is necessary to impose some restrictions on $\left\{p\right(x,z,t),w(x,z,t\left)\right\}$ in order to solve (15) for $p(x,z,t)$ in terms of $\left\{w\right(x,z,t),U(\zeta ),\rho \}$. In the upper fluid ($z>0$), we require that there are positive constants, $K_1$ and $K_2$, such that:

$$\begin{array}{c}{\int }_{-\infty }^{\infty }\left|p_u(x,z,t)\right|dx\le K_1<\infty \mathrm{for}\mathrm{all}z>0,t>0;\end{array}$$

(16a)

$$\begin{array}{c}{\int }_{-\infty }^{\infty }\left|w_u(x,z,t)\right|dx\le K_2<\infty \mathrm{for}\mathrm{all}z>0,t>0;\end{array}$$

(16b)

Similar restrictions apply in the lower fluid for $z<0$ and for $z\to -\infty$, and these restrictions remain valid as $z\to 0$. As $z\to +\infty$, $p_u(x,z,t)$ and $w_u(x,z,t)$ are each uniformly bounded in x and t, and

$$p_u\to 0,w_u\to 0,{\partial }_x{w}_u\to 0.$$

(17)

Next, we solve (15) by the method of Variation of Parameters in the upper fluid and also in the lower fluid.

Comment 3.

The method of Variation of Parameters appears repeatedly in this analysis.

The first step in this method is to ignore the right-hand side of (15), and to solve the homogeneous problem (i.e., the Laplace equation). It follows from (16a) that the solution of the homogeneous problem, $p_h(x,z,t)$, has a bounded Fourier transform, and that the solution of the homogeneous version of (15) can be written as

$$\begin{array}{c}{p}_{h,u}(x,z,t)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{\widehat{p}}_{h,u}(k,t)e^{ikx}{e}^{-\left|k\right|z}dk\end{array}$$

(18a)

for $z>0$, and as

$$\begin{array}{c}\hfill p_{h,l}(x,z,t)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{\widehat{p}}_{h,l}(k,t)e^{ikx}{e}^{\left|k\right|z}dk\end{array}$$

(18b)

for $z<0$.

Comment 4.

A consequence of (18a, 18b) is that in the absence of fluid motion (so ${\partial }_{x}w(x,z,t)=0$ for all finite z in (15)), for every fixed, real-valued $k\ne 0$, the equilibrium pressure field decays exponentially to zero away from $z=0$, like $e^{-\left|k\right|z}$ for $z>0$ and like $e^{\left|k\right|z}$ for $z<0$. This fact appears repeatedly in Section 3.1 and Section 3.2 of these notes, where solutions of equations derived in this Part I are found. In addition, ${\widehat{p}}_{h,*}(k,t)$ is determined uniquely by $p_{h,*}(x,0,t)$, and we have assumed that $p_{h,*}(x,0,t)$ also satisfies (16a).

Comment 5.

The facts in [Comment 4] are irrelevant for our current purpose, because (18a, 18b) is simply a step along the way to the solution of the nonhomogeneous problem.

The second step in Variation of Parameters is to allow the free parameters in the problem to vary in z for $t>0$. Thus, for $z>0$, we replace ${\widehat{p}}_{h,u}(k,t)$ with ${\widehat{p}}_u(k,z,t)$, whose z-dependence is to be determined from (15). Then the functions on the right side of (18a), (18b) are generalized to

$$\begin{array}{c}{p}_u(x,z,t)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{\widehat{p}}_u(k,z,t)·e^{ikx}·e^{-\left|k\right|z}dk\end{array}$$

(19a)

for $z>0$, and to

$$\begin{array}{c}\hfill p_l(x,z,t)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{\widehat{p}}_l(k,z,t)·e^{ikx}·e^{\left|k\right|z}dk\end{array}$$

(19b)

for $z<0$.

Comment 6.

Based on an astute observation from Diane Henderson (private communication), note that the Fourier transform of $p_u$ is not ${\widehat{p}}_u(k,z,t)$, but instead is ${\widehat{p}}_u(k,z,t)·e^{-\left|k\right|z}$, and that the Fourier transform of $p_l$ is ${\widehat{p}}_l(k,z,t)·e^{\left|k\right|z}$.

It follows from (19a), (19b) that for $z>0$,

$$\begin{array}{c}{\nabla }^2{p}_u(x,z,t)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }[{\partial }_z^2{\widehat{p}}_u-2\left|k\right|{\partial }_z{\widehat{p}}_u]e^{ikx}{e}^{-\left|k\right|z}dk;\end{array}$$

(20a)

and for $z<0$,

$$\begin{array}{c}\hfill {\nabla }^2{p}_l(x,z,t)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }[{\partial }_z^2{\widehat{p}}_l+2\left|k\right|{\partial }_z{\widehat{p}}_u]e^{ikx}{e}^{\left|k\right|z}dk.\end{array}$$

(20b)

It is also necessary to write $w(x,z,t)$ in terms of its Fourier transform (in x), in each fluid. It follows form (16b) that w has a uniformly bounded Fourier transform in each fluid, so we write

$$w_{*}(x,z,t)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }[{\widehat{w}}_{*}(k,z,t)·e^{ikx}]dk,$$

(21)

from which it follows that

$${\partial }_x{w}_{*}(x,z,t)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }[ik{\widehat{w}}_{*}(k,z,t)·e^{ikx}]dk.$$

(22)

Substituting (20a), (20b) and (22) into (15) yields

$$\begin{array}{c}[{\partial }_z^2{\widehat{p}}_u-2\left|k\right|{\partial }_z{\widehat{p}}_u]e^{-\left|k\right|z}=-2{\rho }_u[{\partial }_z{U}_u\left(\zeta (z,t)\right)·ik{\widehat{w}}_u]\end{array}$$

(23a)

for $z>0$, and

$$\begin{array}{c}\hfill [{\partial }_z^2{\widehat{p}}_l+2\left|k\right|{\partial }_z{\widehat{p}}_l]e^{\left|k\right|z}=-2{\rho }_l[{\partial }_z{U}_l\left({\zeta }_l(z,t)\right)·ik{\widehat{w}}_l]\end{array}$$

(23b)

for $z<0$. Equivalently,

$$\begin{array}{c}\hfill {\partial }_z[{\partial }_z{\widehat{p}}_u·e^{-2\left|k\right|z}]e^{\left|k\right|z}=-2ik{\rho }_u[{\partial }_z{U}_u\left(\zeta (z,t)\right)]·[{\widehat{w}}_u]\end{array}$$

(23c)

for $z>0$, and

$$\begin{array}{c}\hfill {\partial }_z[{\partial }_z{\widehat{p}}_l·e^{2\left|k\right|z}]e^{-\left|k\right|z}=-2ik{\rho }_l[{\partial }_z{U}_l\left({\zeta }_l(z,t)\right)]·[{\widehat{w}}_l]\end{array}$$

(23d)

for $z<0$. Recall that $\zeta (z,t)$ and ${\zeta }_l(z,t)$ are defined in (3).

The solution of each of the second-order, ordinary differential equations (in z) in (23c), (23d) has two constants of integration. (15) is an elliptic P.D.E, so in each fluid, one of the constants is determined as $\left|z\right|\to \infty$ and the other is determined as $z\to 0$. It follows that

$$\begin{array}{c}{\partial }_z{\widehat{p}}_u(k,z,t)·e^{-2\left|k\right|z}=2ik{\rho }_u{\int }_z^{\infty }[{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}]ds\end{array}$$

(24a)

for $z>0$, and

$$\begin{array}{c}\hfill {\partial }_z{\widehat{p}}_l(k,z,t)·e^{2\left|k\right|z}=-2ik{\rho }_l{\int }_{-\infty }^z[{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}]ds\end{array}$$

(24b)

for $z<0$.

Please note that ${\widehat{w}}_u(k,s,t)$ is bounded as $s\to \infty$ by (16b) and that if $4{\nu }_{u}t>1$, then

$$\begin{array}{c}\hfill |{\partial }_s{U}_u\left(\zeta (s,t)\right)·e^{-\left|k\right|s}|<|K_3|·e^{\frac{-s^2}{4{\nu }_{u}t}}·\frac{1}{\sqrt{4{\nu }_{u}t}}<\left|K_3\right|·e^{\frac{-s^2}{4{\nu }_{u}t}}\end{array}$$

(24c)

for some constant, finite $|K_3|$, so the integral in (24a) decays faster than exponentially as $z\to \infty$, which guarantees that ${\partial }_z{\widehat{p}}_u(k,z,t)\to 0$ as $z\to \infty$. After multiplying (24a) by $e^{2\left|k\right|z}$, one can integrate in z again, to obtain

$$\begin{array}{c}{\widehat{p}}_u(k,z,t)={\widehat{p}}_u(k,0,t)\\ +2ik{\rho }_u{\int }_0^z[e^{2\left|k\right|q}\left({\int }_q^{\infty }\{{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}\}ds\right)]dq.\end{array}$$

(25a)

Similarly,

$$\begin{array}{c}{\widehat{p}}_l(k,z,t)={\widehat{p}}_l(k,0,t)\\ +2ik{\rho }_l{\int }_z^0[e^{-2\left|k\right|q}\left({\int }_{-\infty }^q\{{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}\}ds\right)]dq.\end{array}$$

(25b)

Multiply (25a) by $e^{-\left|k\right|z}$, and multiply (25b) by $e^{\left|k\right|z}$ to obtain

$$\begin{array}{c}{\widehat{p}}_u(k,z,t)·e^{-\left|k\right|z}={\widehat{p}}_u(k,0,t)·e^{-\left|k\right|z}\\ +2ik{\rho }_u·e^{-\left|k\right|z}{\int }_0^z[e^{2\left|k\right|q}\left({\int }_q^{\infty }\{{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}\}ds\right)]dq.\end{array}$$

(26a)

for $z>0$, and

$$\begin{array}{c}1-1{\widehat{p}}_l(k,z,t)·e^{\left|k\right|z}={\widehat{p}}_l(k,0,t)·e^{\left|k\right|z}\\ +2ik{\rho }_l·e^{\left|k\right|z}{\int }_z^0[e^{-2\left|k\right|q}\left({\int }_{-\infty }^q\{{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}\}ds\right)]dq.\end{array}$$

(26b)

for $z<0$.

An important objective of this section is to find explicit formulae for the perturbation pressure in the upper and lower fluids, respectively. The Fourier transform in x of $p_u(x,z,t)$ is ${\widehat{p}}_u(k,z,t)·e^{-\left|k\right|z}$, given by (26a); the Fourier transform of $p_l(x,z,t)$ is ${\widehat{p}}_l(k,z,t)·e^{\left|k\right|z}$, given by (26b).

Given the explicit formulae in (26a), (26b) for the Fourier transforms of the perturbation pressures in the upper and lower fluids, respectively, the inverse Fourier transforms of these two equations give respectively the perturbation pressures in the upper and lower fluids. One shows by direct calculation that the pressures whose Fourier transforms are given in (26a), (26b) satisfy (15) and (16a) and (16b).

What are needed in (14b) are $\frac{1}{{\rho }_u}{\partial }_z({\widehat{p}}_u(k,z,t)·e^{-\left|k\right|z})$ in the upper fluid (for $z\ge 0$), and $\frac{1}{{\rho }_l}{\partial }_z({\widehat{p}}_l(k,z,t)·e^{\left|k\right|z})$ in the lower fluid (for $z\le 0$). These can be computed directly from (26a), (26b). For $z\ge 0$,

$$\begin{array}{c}\frac{1}{{\rho }_u}{\partial }_z({\widehat{p}}_u(k,z,t)·e^{-\left|k\right|z})=-\frac{\left|k\right|}{{\rho }_u}{\widehat{p}}_u(k,0,t)·e^{-\left|k\right|z}\\ -2ik\left|k\right|·e^{-\left|k\right|z}{\int }_0^z[e^{2\left|k\right|q}{\int }_q^{\infty }\{{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}\}ds]dq\\ +2ik·e^{-\left|k\right|z}·e^{2\left|k\right|z}{\int }_z^{\infty }[{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}]ds,\end{array}$$

(27a)

and for $z\le 0$,

$$\begin{array}{c}\hfill 1-1\frac{1}{{\rho }_l}{\partial }_z({\widehat{p}}_l(k,z,t)·e^{\left|k\right|z})=\frac{\left|k\right|}{{\rho }_l}{\widehat{p}}_l(k,0,t)·e^{\left|k\right|z}\\ \hfill +2ik\left|k\right|·e^{\left|k\right|z}{\int }_z^0[e^{-2\left|k\right|q}{\int }_{-\infty }^q\{{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}\}ds]dq\\ \hfill -2ik·e^{\left|k\right|z}·e^{-2\left|k\right|z}{\int }_{-\infty }^z[{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}]ds.\end{array}$$

(27b)

2.4. An Equation for $\widehat{w}(k,z,t)$ in Each Fluid

The starting point for this part of the analysis is (14b). Taking the Fourier transform (in x) of (14b) and using (21) yields one evolution equation for ${\widehat{w}}_u(k,z,t)$ in the upper fluid ($z\ge 0$),

$$\begin{array}{c}{\partial }_t{\widehat{w}}_u+ik·U_u\left(\zeta (z,t)\right)·{\widehat{w}}_u+\frac{1}{{\rho }_u}{\partial }_z({\widehat{p}}_u(k,z,t)·e^{-\left|k\right|z})={\nu }_u({\partial }_z^2{\widehat{w}}_u-k^2{\widehat{w}}_u),\end{array}$$

(28a)

and a second evolution equation for ${\widehat{w}}_l(k,z,t)$ in the lower fluid ($z\le 0$),

$$\begin{array}{c}\hfill {\partial }_t{\widehat{w}}_l+ik·U_l\left({\zeta }_l(z,t)\right)·{\widehat{w}}_l+\frac{1}{{\rho }_l}{\partial }_z({\widehat{p}}_l(k,z,t)·e^{\left|k\right|z})={\nu }_l({\partial }_z^2{\widehat{w}}_l-k^2{\widehat{w}}_l).\end{array}$$

(28b)

Use (27a) to replace $\frac{1}{{\rho }_u}·{\partial }_z({\widehat{p}}_u(k,z,t)·e^{-\left|k\right|z})$ with equivalent integrals over ${\widehat{w}}_u$, and use (27b) to replace $\frac{1}{{\rho }_l}·{\partial }_z({\widehat{p}}_l(k,z,t)·e^{\left|k\right|z})$ with equivalent integrals over ${\widehat{w}}_l$, and obtain (for $z\ge 0$)

$$\begin{array}{c}{\partial }_t{\widehat{w}}_u+ik·U_u\left(\zeta (z,t)\right)·{\widehat{w}}_u\\ -2ik\left|k\right|·e^{-\left|k\right|z}{\int }_0^z[e^{2\left|k\right|q}{\int }_q^{\infty }\{{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}\}ds]dq\\ +2ik·e^{\left|k\right|z}{\int }_z^{\infty }[{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}]ds\\ =\frac{\left|k\right|}{{\rho }_u}·{\widehat{p}}_u(k,0,t)·e^{-\left|k\right|z}+{\nu }_u({\partial }_z^2-k^2){\widehat{w}}_u,\end{array}$$

(29a)

and obtain (for $z\le 0$)

$$\begin{array}{c}{\partial }_t{\widehat{w}}_l+ik·U_l\left({\zeta }_l(z,t)\right)·{\widehat{w}}_l\\ +2ik\left|k\right|·e^{\left|k\right|z}{\int }_z^0[e^{-2\left|k\right|q}{\int }_{-\infty }^q\{{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}\}ds]dq\\ -2ik·e^{-\left|k\right|z}{\int }_{-\infty }^z[{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}]ds\\ =-\frac{\left|k\right|}{{\rho }_l}·{\widehat{p}}_l(k,0,t)·e^{\left|k\right|z}+{\nu }_l({\partial }_z^2-k^2){\widehat{w}}_l.\end{array}$$

(29b)

(29a) describes the dynamics of the upper fluid, which is the region beginning from $z=0$ and extending as $z\to \infty$. On the right-hand side of (29a), the term ${\widehat{p}}_u(k,0,t)$ provides external pressure, defined at the fixed interface, and it can change over time. As one moves away from the interface, but still in the upper fluid, this external pressure field weakens according to ${\widehat{p}}_u(k,z=0,t)·e^{-\left|k\right|z}$, for $z>0$. This pressure field operates as an external forcing function in (29a). Another external variable is $U_u\left(\zeta (s,t)\right)$, defined for $0<s<\infty$. This variable determines the strength of the external forcing due to the integral terms. Please note that (29b) has a similar structure to (29a): the external pressure field is defined by ${\widehat{p}}_l(k,z=0,t)$, and another other external variable is $U_l\left({\zeta }_l(s,t)\right)$. In the lower fluid ($z<0$), the pressure field weakens according to ${\widehat{p}}_l(k,z=0,t)·e^{\left|k\right|z}$. As shown in Section 3.1 and Section 3.2, exploiting this structure provides an avenue to solve these equations exactly. At this point in the analysis, the velocity field in the upper fluid evolves according to (29a), while the velocity field in the lower fluid evolves according to (29b).

In either fluid, each of these terms provides information about a region away from the interface. ${\widehat{p}}_{*}(k,z=0,t)$ provides information that appears at the interface, while $U_{*}\left(\zeta (s,t)\right)$ provides information from regions above or below $z=0$. Although these terms both provide external forcing, they have no relation to each other. Each forcing term decays exponentially fast, or faster, as one moves away from the interface, in either fluid. Hence, once these two forcing functions are known, then the velocity fields in the upper fluid and in the lower fluid are determined completely by these forcing functions and by the initial data in each fluid. In addition, the pressure fields in the upper and lower fluids are defined in (26a), (26b), each in terms of the velocity field in its own fluid, plus a similar forcing function (at $z=0$) for that fluid. Finally, the perturbation horizontal velocity in each fluid is determined by the vertical velocity in its own fluid, according to (14c).

Comment 7.

In addition to the external pressure field in each fluid, a second external forcing field in the upper fluid, of the form ${\partial }_s{U}_u\left(\zeta (s,t)\right)$, appears in the integral terms in (29a) while a similar external forcing term, of the form ${\partial }_s{U}_l\left({\zeta }_l(s,t)\right)$, appears in the integral terms in (29b).

Comment 7a.

(29a) is an evolution equation (in time) for ${\widehat{w}}_u(k,z,t)$ with $z>0$. ${\widehat{w}}_u(k,z,t)$ is complex valued, so within the term is an amplitude and a phase.

• Keeping only the first line of (29a) provides an approximate equation for the evolution of the phase, with the solution:

  $${\widehat{w}}_u(k,z,t)={\widehat{w}}_u(k,z,0)·e^{-ik{\int }_0^t{U}_u\left(\zeta (z,\tau )\right)d\tau }$$

  (30)
• Keeping the first three lines of (29a) changes the phase of ${\widehat{w}}_u(k,z,t)$ in a complicated way.
• Keeping the t-derivative from the first line of (29a), plus the boundary-pressure term from the 4th line, increases the amplitude of ${\widehat{w}}_u(k,z,t)$, provided ${\widehat{p}}_u(k,0,t)>0$.
• Keeping the t-derivative of ${\widehat{w}}_u(k,z,t)$ from the first line plus the last set of terms in the 4th line creates the heat equation in two spatial dimensions for ${\widehat{w}}_u(k,z,t)$, so this combination introduces the diffusion of ${\widehat{w}}_u(k,z,t)$

Comment 7b.

The mathematical structure for (29b) is very similar to that of (29a), except for some signs.

Comment 7c.

The objectives of Section 3.1 and Section 3.2 are to describe the dynamics of the fluid in the upper fluid (in Section 3.1), and the dynamics of the lower fluid (in Section 3.2). The governing equations for the upper fluid are (44), along with (45a), (45b), (45c). These equations are equivalent to (29a), with somewhat different notation. The governing equations for the lower fluid are (73), along with (74a), (74b), (74c). These equations are equivalent to (29b), with somewhat different notation. The most important results in Section 3.1 and Section 3.2 are families of exact solutions of (44) and (73), respectively. The families of exact solutions for the upper fluid, in Section 3.1, are constructed by solving a sequence of increasingly general truncations of (44). The corresponding family of exact solutions for the lower fluid, in Section 3.2, are also obtained by solving increasingly general truncations of (73). This procedure, based on solving sequential truncations of (44) and/or (73), applies equally well to constructing solutions of (29a) and (29b).

2.5. Linearized Equations at the Interface

The next steps enforce the compatibility of the two vertical velocity fields and of the two pressure fields where they meet, at the interface between the two fluids. The compatibility conditions at the interface are a set of physically plausible rules that constrain the flow. For example, (29a), (29b) are integro-differential equations for the evolution of ${\widehat{w}}_{*}(k,z,t)$ in each fluid, with k as a fixed parameter. These equations need boundary conditions in z before we can solve the initial-value problem in t. The boundary conditions as $z\to +\infty$ in the upper fluid, or as $z\to -\infty$ in the lower fluid, are clear: ${\widehat{w}}_{*}(k,z,t)\to 0$ in each limit. The conditions at the interface are more complicated. The location of the (movable) interface is defined by

$$\begin{array}{c}z=ϵ\eta (x,t)+\mathcal{O}\left({ϵ}^2\right)\end{array}$$

(31a)

where $\eta (x,t)$ is differentiable with respect to both ($x,t$) and we require that

$$\begin{array}{c}\hfill {\int }_{-\infty }^{\infty }\left|\eta (x,t)\right|dx<\infty \end{array}$$

(31b)

For $t\ge 0$, $\eta (x,t)$ evolves according to

$$\begin{array}{c}\hfill {\partial }_t\eta +U\left(0\right)·{\partial }_x\eta =w(x,0,t)\end{array}$$

(31c)

where $U\left(0\right)$ is the background, horizontal velocity, and $w(x,0,t)$ is the vertical velocity, both at the interface (where $z=0$). As discussed in [Comment 1] in Section 2.2, $\widehat{\overline{w}}(k,t)$ is the limit as $z\to 0$ of both ${\widehat{w}}_u(k,z,t)$ and ${\widehat{w}}_l(k,z,t)$, the solutions of (29a) in the upper fluid and of (29b) in the lower fluid. The evolution of $\widehat{\eta }(k,t)$, the Fourier transform of $\eta (x,t)$, is obtained by taking the Fourier transform of (31c):

$${\partial }_t\widehat{\eta }(k,t)+ik·U\left(0\right)·\widehat{\eta }(k,t)=\widehat{\overline{w}}(k,t).$$

(32)

$\widehat{\eta }(k,t)$ is complex-valued, so it has a phase and an amplitude. Setting the left hand side of (32) equal to 0 provides an evolution equation for the phase, with the solution,

$$\widehat{\eta }(k,t)=\widehat{\eta }(k,0)·e^{-ik{\int }_0^t{U}_u\left(\zeta (z,\tau )\right)d\tau }.$$

Please note that this equation is the same as in (30), i.e., the vertical velocity, $\widehat{w}(k,z,t)$, and the height of the movable interface, $\widehat{\eta }(k,t)$, propagate at the same speed and in the same direction. Keeping only the first and last terms in (32) yields an evolution equation for the amplitude of $\widehat{\eta }(k,t)$:

$$\widehat{\eta }(k,t)=\widehat{\eta }(k,0)+{\int }_0^t\widehat{\overline{w}}(k,\tau )d\tau$$

The evolution equation for $\widehat{\overline{w}}(k,t)$ follows from evaluating each of (29a) and (29b) at $z=0$, and then adding the two equations. The result is that the first integral in (29a) and the first integral in (29b) are both eliminated:

$$\begin{array}{c}2\{{\partial }_t\widehat{\overline{w}}(k,t)+ik·U\left(0\right)·\widehat{\overline{w}}(k,t)\}\\ -\{{\nu }_u{\partial }_z^2{\widehat{w}}_u(k,0,t)+{\nu }_l{\partial }_z^2{\widehat{w}}_l(k,0,t)\}+\{{\nu }_u+{\nu }_l\}{k}^2\widehat{\overline{w}}\\ +2ik·{\int }_0^{\infty }[{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}]ds\\ -2ik·{\int }_{-\infty }^0[{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}]ds\\ =|k|·(\frac{{\widehat{p}}_u(k,0,t)}{{\rho }_u}-\frac{{\widehat{p}}_l(k,0,t)}{{\rho }_l}).\end{array}$$

(33)

(33) is obtained by adding (29b) and (29b). Subtracting (29b) from (29a) provides a necessary condition for ${\widehat{w}}_{*}(k,z,t)$ to be continuous across the interface (at $z=0$) for all $t\ge 0$ is

$$\begin{array}{c}{\nu }_l{\partial }_z^2{\widehat{w}}_l(k,0,t)-{\nu }_u{\partial }_z^2{\widehat{w}}_u(k,0,t)+\{{\nu }_u-{\nu }_l\}{k}^2\widehat{\overline{w}}(k,t)\\ +2ik·{\int }_0^{\infty }[{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}]ds\\ +2ik·{\int }_{-\infty }^0[{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}]ds\\ =|k|·(\frac{{\widehat{p}}_u(k,0,t)}{{\rho }_u}+\frac{{\widehat{p}}_l(k,0,t)}{{\rho }_l}).\end{array}$$

(34)

When considering continuity and discontinuity across the interface, there are two issues: What happens at $\mathcal{O}\left(1\right)$? What happens at $\mathcal{O}\left(ϵ\right)$?

• Section 2.1 describes what happens at $\mathcal{O}\left(1\right)$. Equations (5a) and (5b) show that the background, horizontal velocity is continuous across the interface at $\mathcal{O}\left(1\right)$. There is no vertical velocity at $\mathcal{O}\left(1\right)$. (38b), below, shows that the tangential (or shear) stress must be continuous across the interface. Then it follows that ${\partial }_{z}U\left(\zeta \left(z\right)\right)$ must be discontinuous at $z=0$, unless ${\rho }_u{\nu }_u={\mu }_u={\rho }_l{\nu }_l={\mu }_l$. (12a) and 12b) shows that the pressure is still continuous at $\mathcal{O}\left(1\right)$ at the fixed interface (where $z=0$), and (25a), (25b) show that the pressure can be discontinuous at $\mathcal{O}\left(ϵ\right)$ at the fixed interface.
• As noted in [Comment 1] in Section 2.2, assuming that two fluids are incompressible suggests that the vertical velocity is continuous at $z=0$. Even so, incompressibility is not enough to guarantee that ${\widehat{w}}_{*}(k,z,t)$ is continuous across the interface. Equation (34) is also needed for continuity of ${\widehat{w}}_{*}(k,z,t)$ at the interface (at $z=0$).
• We show next that ${\partial }_z{\widehat{w}}_{*}(k,z,t)$ is also continuous at $z=0$, so

  $$\begin{array}{c}{\partial }_z{\widehat{w}}_l(k,z,t)-{\partial }_z{\widehat{w}}_u(k,z,t)\to 0\mathrm{as}z\to 0.\end{array}$$

  (35a)

Proof. 

–

The Fourier transform in x of (14c) on each side of the interface is

$$\begin{array}{c}{\partial }_z{\widehat{w}}_u(k,z,t)+ik{\widehat{u}}_u(k,z,t)=0,\mathrm{for}z>0,\\ {\partial }_z{\widehat{w}}_l(k,z,t)+ik{\widehat{u}}_l(k,z,t)=0,\mathrm{for}z<0.\end{array}$$

(36)

–

In order to have a finite shear stress at the interface, the horizontal velocity must be continuous at the interface, not only at $\mathcal{O}\left(1\right)$, but also at $\mathcal{O}\left(ϵ\right)$. Therefore

$$ik[{\widehat{u}}_u(k,z,t)-{\widehat{u}}_l(k,z,t)]\to 0\mathrm{as}z\to 0$$

(37)

–

Combining (36) and (37) guarantees that ${\partial }_z{\widehat{w}}_u(k,z,t)$ and ${\partial }_z{\widehat{w}}_l(k,z,t)$ approach a common limit as $z\to 0$, as asserted in (35a). □

In addition to continuity of each component of velocity, it is also necessary that the normal and tangential stresses be continuous across the interface. The full, nonlinear formulae for these two conditions can be found on p. 150 of Batchelor’s Introduction to Fluid Mechanics [5]. The full (i.e., not linearized) condition for continuity of tangential stress across the interface (at $z=ϵ\eta$) is

$$\begin{array}{c}{\rho }_u{\nu }_u[4\left({\partial }_x\tilde{\eta }\right)\left({\partial }_z{\tilde{w}}_u\right)+(1-{\left({\partial }_x\tilde{\eta }\right)}^2)({\partial }_z{\tilde{u}}_u+{\partial }_x{\tilde{w}}_u)]\\ ={\rho }_l{\nu }_l[4\left({\partial }_x\tilde{\eta }\right)\left({\partial }_z{\tilde{w}}_l\right)+(1-{\left({\partial }_x\tilde{\eta }\right)}^2)({\partial }_z{\tilde{u}}_l+{\partial }_x{\tilde{w}}_l)]\end{array}$$

(38a)

where again the variables in (38a) are the full nonlinear variables, not merely the perturbation variables. When one expands the variables in (38a) in powers of $ϵ$$(0<ϵ<<1)$, the $\mathcal{O}\left(1\right)$ terms assert that

$$\begin{array}{c}\hfill {\rho }_u{\nu }_u{[{\partial }_z{U}_u\left(\zeta \left(z\right)\right)]}_{z=0}={\rho }_l{\nu }_l{[{\partial }_z{U}_l\left({\zeta }_l\left(z\right)\right)]}_{z=0},\end{array}$$

(38b)

which is a restatement of (9). At $\mathcal{O}\left(ϵ\right)$, (38a) gives

$$\begin{array}{c}\hfill {\rho }_u{\nu }_u{[({\partial }_z{u}_u+{\partial }_x{w}_u)]}_{z=0}={\rho }_l{\nu }_l{[({\partial }_z{u}_l+{\partial }_x{w}_l)]}_{z=0}.\end{array}$$

(38c)

Comment 8.

Keep in mind that the velocity components in (38c) are the components of the linearized velocity field, while those in (38a) represent the full nonlinear flow.

Comment 9.

(37) guarantees that the tangential stress at the interface is finite. (38c) guarantees that the tangential stress is continuous at the interface.

As we show next, a consequence of (37) is that ${\partial }_z^2\widehat{w}(k,z,t)$ is usually not continuous across the interface. Instead, continuity of the tangential stress across the interface leads to a well-defined jump condition for ${\partial }_z^2\widehat{w}(k,z,t)$. A feature of the linearized problem is that $u(x,z,t)$, the horizontal component of the perturbation velocity field, plays no significant role in the linearized dynamics—after finding $w(x,z,t)$ in each fluid, one finds $u(x,z,t)$ a posteriori by using (14c). Even so, there is useful information in (37). From (36):

$$\begin{array}{c}ik\widehat{u}(k,z,t)+{\partial }_z\widehat{w}(k,z,t)=0\Rightarrow \widehat{u}(k,z,t)=\frac{i}{k}{\partial }_z\widehat{w}(k,z,t)\end{array}$$

(39a)

Inserting this into the Fourier transform of (38c) leads to

$$\begin{array}{c}\hfill {\rho }_u{\nu }_u{[{\partial }_z^2{\widehat{w}}_u+k^2{\widehat{w}}_u]}_{z=0}={\rho }_l{\nu }_l{[{\partial }_z^2{\widehat{w}}_l+k^2{\widehat{w}}_l]}_{z=0}.\end{array}$$

(39b)

From [Comment 1] in Section 2.2, we assume

$$\begin{array}{c}\hfill {\widehat{w}}_u{{|}_{z=0}={\widehat{w}}_l|}_{z=0}=\widehat{\overline{w}}(k,t).\end{array}$$

Inserting this into (39b) yields

$$\begin{array}{c}\hfill {\rho }_u{\nu }_u·{\partial }_z^2{\widehat{w}}_u(k,0,t)-{\rho }_l{\nu }_l·{\partial }_z^2{\widehat{w}}_l(k,0,t)=\{{\rho }_l{\nu }_l-{\rho }_u{\nu }_u\}{k}^2\widehat{\overline{w}}(k,t).\end{array}$$

(39c)

(39c) is a jump condition for ${\partial }_z^2{\widehat{w}}_{*}(k,z,t)$ across the interface (at $z=0$). One can show from (39c) that ${\partial }_z^2{\widehat{w}}_{*}(k,z,t)$ can be continuous across the interface only if

$$\begin{array}{c}\hfill {\rho }_u{\nu }_u={\rho }_l{\nu }_l\end{array}$$

(39d)

Return to Batchelor’s book. For normal stress, the full (i.e., not linearized) condition is that on $z=ϵ\eta (x,t)$,

$$\begin{array}{c}{\tilde{P}}_u-\frac{2{\rho }_u{\nu }_u}{[1+{\left({\partial }_x\tilde{\eta }\right)}^2]}[(1-{\left({\partial }_x\tilde{\eta }\right)}^2)\left({\partial }_z{\tilde{w}}_u\right)-\left({\partial }_x\tilde{\eta }\right)({\partial }_z{\tilde{u}}_u+{\partial }_x{\tilde{w}}_u)]=\\ {\tilde{P}}_l-\frac{2{\rho }_l{\nu }_l}{[1+{\left({\partial }_x\tilde{\eta }\right)}^2]}[(1-{\left({\partial }_x\tilde{\eta }\right)}^2)\left({\partial }_z{\tilde{w}}_l\right)-\left({\partial }_x\tilde{\eta }\right)({\partial }_z{\tilde{u}}_l+{\partial }_x{\tilde{w}}_l)]\end{array}$$

(40)

where we revert to the full, nonlinear variables, given in (12a) and (12b), and $\rho \nu =\mu$ in each fluid. When we expand the variables in (39a) in powers of a small $ϵ$$(0<ϵ<<1)$, the leading-order terms assert that the background pressure is continuous across $z=0$, as guaranteed by (12a) and (12b). At $\mathcal{O}\left(ϵ\right)$, there are three contributions on each side of (40): the background pressure ($\tilde{P}$) is evaluated at $z=ϵ\eta$; the perturbation pressures ($p_u(x,z,t)$, $p_l(x,z,t)$) are each evaluated in their respective limits as $z\to 0$; and $\{{\partial }_z{w}_u(x,z,t),{\partial }_z{w}_l(x,z,t)\}$ are each evaluated at $z=0$. The result at $\mathcal{O}\left(ϵ\right)$ can be written formally as

$$\begin{array}{c}({\rho }_l-{\rho }_u)g·\eta (x,t)+\{p_u(x,0,t)-p_l(x,0,t)\}-2\{{\mu }_u{\partial }_z{w}_u(x,0,t)-{\mu }_l{\partial }_z{w}_l(x,0,t)\}=0.\end{array}$$

(41)

The Fourier transform (in x) of this relation is

$$\begin{array}{c}({\rho }_l-{\rho }_u)g·\widehat{\eta }(k,t)+\{{\widehat{p}}_u(k,0,t)-{\widehat{p}}_l(k,0,t)\}-2\{{\mu }_u{\partial }_z{\widehat{w}}_u(k,0,t)-{\mu }_l{\partial }_z{\widehat{w}}_l(k,0,t)\}=0.\end{array}$$

(42)

Comment 10.

(42) is correct if there is no surface tension at the interface. In the presence of surface tension, (42) generalizes to

$$\begin{array}{c}({\rho }_l-{\rho }_u)(g+\sigma k^2)\widehat{\eta }(k,t)=\\ \{{\widehat{p}}_l(k,0,t)-{\widehat{p}}_u(k,0,t)\}-2({\rho }_l{\nu }_l-{\rho }_u{\nu }_u)·{\partial }_z{\widehat{w}}_{*}(k,0,t),\end{array}$$

(43)

where ${\rho }_{*}{\nu }_{*}={\mu }_{*}$, the dynamic viscosity in each fluid, $\sigma$ represents the surface tension, and we have used (35a).

3. Solve the Equations for Flow in Each Fluid

Objectives of Section 3.1 and Section 3.2: In Section 2.4, Equations (29a) and (29b) define the vertical velocity in the upper and lower fluids, respectively. These equations each contain two known forcing functions, $\{U_u\left(\zeta (z,t)\right),{\widehat{p}}_u(k,z=0,t)\}$ in the upper fluid ($z\ge 0$), and $\{U_l\left({\zeta }_l(z,t)\right),{\widehat{p}}_l(k,z=0,t)\}$ in the lower fluid ($z\le 0$). In Section 3.1 and Section 3.2, Equations (44) and (73) are rewrites of (29a) and (29b), respectively. The objective of Section 3.1 is to find general solutions of (44); these solutions describe the vertical velocity in the upper fluid. The objective of Section 3.2 is to find general solutions of (73); these solutions describe the vertical velocity in the lower fluid.

Comment 11.

The analysis in Section 3.1 shows that (44) admits two distinct families of solutions, so seeking “the general solution of (44)” is too restricted. Instead, in Section 3.1 we seek the most general solution within each of the two families of solutions. In Section 3.1, the method of Variation of Parameters is used repeatedly to find a general solution of (44) within each of the two families. The mathematical structures of (44) and (73) are very similar, so we describe the analysis in Section 3.1 in detail, and carry it to the intermediate conclusions of that part of the analysis. In Section 3.2, it is sufficient simply to identify the important parts of the analysis and follow with the corresponding intermediate conclusions. To simplify the logic of the analysis of the vertical motion in the upper fluid, in Section 3.1, and the corresponding logic of the analysis of the vertical motion in the lower fluid, we combine the two sections into one, Section 3.1 and Section 3.2. Consistently in these sections, we analyze the vertical motion in the upper fluid, followed by nearly identical analysis of the motion in the lower fluid.

Comment 12.

The procedure to generate an ordered sequence of intermediate results in each fluid is straight-forward, but rather involved. The starting point in Section 3.1 is equation (44), which is a rewrite of (29a). For Section 3.2, the starting point is (73), which is a rewrite of (29b). Each of (44) and (73) can be written as the sum of four terms. In either case, step 1 in the analysis is to create and then solve a truncated version of the equation in (44) or (73). In step 1, the truncated equation for each fluid contains only the two leftmost terms in (44) or (73). The solution of this truncated equation can be considered a first approximation of a solution of (44) or (73). Like the solution of (44) or (73), the solution of the truncated equation is complex-valued, so it can be written in terms of its real and imaginary parts. Step 2 involves splitting the solution of the truncated problem into real and imaginary parts. Doing so doubles the numbers of equations to be solved, but it simplifies the structure of each equation. Step 3 generalizes the truncated equation in each fluid by re-inserting the (known) pressure term at the interface in each fluid. These pressure terms are also complex-valued, so they can be written in terms of their real and imaginary parts. One can solve these generalized but approximate equations by using Variation of Parameters again, based on the solutions in steps 1, 2. Step 4 brings the integral term in (44) or in (73) back into the approximate equations obtained in step 3, and thereby reproduces the equations in (44) and (73). In this round-about way, we are able to solve the equations given in (44) and (73) exactly, by starting with the approximate solutions obtained from steps 1, 2, 3, then using Variation of Parameters one more time to obtain the exact solutions for (44) and (73), with no approximations.

3.1. Solve (29a), Also Known as (44)

For the upper fluid with $z\ge 0$, one may write (29a) in the following form:

$${\partial }_t{\widehat{w}}_u+D_u{\widehat{w}}_u+i{L}_u{\widehat{w}}_u=P_u(k,z,t),$$

(44)

where

$$\begin{array}{c}{D}_u{\widehat{w}}_u=-{\nu }_u[{\partial }_z^2{\widehat{w}}_u-k^2{\widehat{w}}_u],\end{array}$$

(45a)

$$\begin{array}{c}{P}_u(k,z,t)=\frac{\left|k\right|}{{\rho }_u}·{\widehat{p}}_u(k,0,t)·e^{-\left|k\right|z},\end{array}$$

(45b)

and

$$\begin{array}{c}\hfill i{L}_u{\widehat{w}}_u(k,z,t)=ik·U_u\left(\zeta (z,t)\right)·{\widehat{w}}_u(k,z,t)\\ \hfill -2ik\left|k\right|·e^{-\left|k\right|z}{\int }_0^z[e^{2\left|k\right|q}{\int }_q^{\infty }\{{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}\}ds]dq\\ \hfill +2ik·e^{\left|k\right|z}{\int }_z^{\infty }[{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{w}}_u(k,s,t)·e^{-\left|k\right|s}]ds.\end{array}$$

(45c)

Both of the differential operators on the left side of (44) are real-valued, while ${\widehat{w}}_u(k,z,t)$ and ${\widehat{p}}_u(k,0,t)$ are typically complex-valued functions. The integral operator, $i{L}_u(·)$, is purely imaginary, so it maps a real-valued function into a purely imaginary function, and vice-versa. The real and imaginary parts of ${\widehat{w}}_u(k,z,t)$ interact with each other only through the integral term in (44), given explicitly in (45c). Step 1 in solving (44) is to create and then solve a truncated version of (44), which contains only the first two terms of (43):

$${\partial }_t{\tilde{w}}_u+D_u{\tilde{w}}_u={\partial }_t{\tilde{w}}_u-{\nu }_u({\partial }_z^2{\tilde{w}}_u-k^2{\tilde{w}}_u)=0$$

(46)

where ${\tilde{w}}_u(k,z,t)$ denotes a (complex-valued) solution of (46), but it need not solve (44). ${\tilde{w}}_u(k,z,t)$ is the simplest non-trivial approximation of ${\widehat{w}}_u(k,z,t)$. Subsequent analysis provides a sequence of increasingly accurate approximations of ${\widehat{w}}_u(k,z,t)$, ending with an explicit formula for ${\widehat{w}}_u(k,z,t)$. ${\tilde{w}}_u(k,z,t)$ is complex-valued, so we can write ${\tilde{w}}_u(k,z,t)$ in terms of its real and imaginary parts:

$${\tilde{w}}_u(k,z,t)={\tilde{\alpha }}_u(k,z,t)+i{\tilde{\beta }}_u(k,z,t)$$

(47)

where ${\tilde{\alpha }}_u(k,z,t)$ and ${\tilde{\beta }}_u(k,z,t)$ are each real-valued functions of real-valued parameters $\{k,z,t\}$. This is Step 2.

Comment 13.

Step 2 is brought into the analysis after step 1, but it can be deferred until after step 3.

The differential operator in (46) is purely real, so ${\tilde{\alpha }}_u(k,z,t)$ and $i{\tilde{\beta }}_u(k,z,t)$ are the real and imaginary parts of ${\tilde{w}}_u(k,z,t)$, each satisfying (46) on their own:

$$\begin{array}{cc}\hfill {\partial }_t{\tilde{\alpha }}_u+D_u{\tilde{\alpha }}_u& ={\partial }_t{\tilde{\alpha }}_u-{\nu }_u({\partial }_z^2{\tilde{\alpha }}_u-k^2{\tilde{\alpha }}_u)=0,\hfill \end{array}$$

(48a)

$$\begin{array}{cc}\hfill {\partial }_t{\tilde{\beta }}_u+D_u{\tilde{\beta }}_u& ={\partial }_t{\tilde{\beta }}_u-{\nu }_u({\partial }_z^2{\tilde{\beta }}_u-k^2{\tilde{\beta }}_u)=0.\hfill \end{array}$$

(48b)

Note from (45a) that $D_u{\tilde{w}}_u=0$ is proportional to the Fourier transform (in x) of the two-dimensional Laplace equation, so (46) can be viewed as the Fourier transform (in x) of the heat equation in two spatial dimensions.

Comment 14.

The Fourier cosine transform and the Fourier sine transform both play roles in the analysis that follows, and we treat them separately.

One can verify by direct substitution that (46) admits solutions in the form of a Fourier cosine transform,

$$\begin{array}{cc}\hfill {\tilde{w}}_{u,c}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\tilde{\omega }}_{u,c}(k,m)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm,\hfill \end{array}$$

(49a)

where ${\tilde{\omega }}_{u,c}(k,m)$ is an arbitrary, complex-valued function of $(k,m)$, and $(k,m)$ are each real-valued. Because ${\tilde{\alpha }}_{u,c}(k,z,t)$ and $i{\tilde{\beta }}_{u,c}(k,z,t)$, the real and imaginary parts of ${\tilde{w}}_{u,c}(k,z,t)$, each satisfy (46), it follows that ${\tilde{\alpha }}_{u,c}(k,z,t)$ and $i{\tilde{\beta }}_{u,c}(k,z,t)$ each have representations in the form of (49a):

$$\begin{array}{cc}\hfill {\tilde{\alpha }}_{u,c}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\tilde{A}}_{u,c}(k,m)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm,\hfill \end{array}$$

(49b)

$$\begin{array}{cc}\hfill i{\tilde{\beta }}_{u,c}(k,z,t)& ={\int }_{-\infty }^{\infty }[i{\tilde{B}}_{u,c}(k,m)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm,\hfill \end{array}$$

(49c)

where ${\tilde{A}}_{u,c}(k,m)$ and ${\tilde{B}}_{u,c}(k,m)$ are both real-valued functions, to be determined. It follows from (47) that

$$\begin{array}{cc}\hfill 1-1{\tilde{\omega }}_{u,c}(k,m)& ={\tilde{A}}_{u,c}(k,m)+i{\tilde{B}}_{u,c}(k,m).\hfill \end{array}$$

(49d)

Evaluating (49b) at $t=0$ shows that ${\tilde{A}}_{u,c}(k,m)$ is proportional to the Fourier cosine transform (in z) of the real part of ${\tilde{w}}_u(k,z,0)$, which relates back to initial data for the 2-D heat equation on a half-space: $\{-\infty <x<\infty ,0\le z<\infty \}$. A similar statement holds for ${\tilde{B}}_{u,c}(k,m)$, in (49c).

Note: (47) also admits solutions of the form

$$\begin{array}{cc}\hfill {\tilde{w}}_{u,s}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\tilde{\omega }}_{u,s}(k,m)·sin\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm.\hfill \end{array}$$

(50a)

Therefore, the real and imaginary parts of ${\tilde{w}}_u(k,z,t)$ each have solutions for (48a), (48b) of the same form:

$$\begin{array}{cc}\hfill {\tilde{\alpha }}_{u,s}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\tilde{A}}_{u,s}(k,m)·sin\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm,\hfill \end{array}$$

(50b)

$$\begin{array}{cc}\hfill i{\tilde{\beta }}_{u,s}(k,z,t)& ={\int }_{-\infty }^{\infty }[i{\tilde{B}}_{u,s}(k,m)·sin\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm.\hfill \end{array}$$

(50c)

As in (49d),

$$\begin{array}{cc}\hfill 1-1{\tilde{\omega }}_{u,s}(k,m)& ={\tilde{A}}_{u,s}(k,m)+i{\tilde{B}}_{u,s}(k,m).\hfill \end{array}$$

(50d)

where ${\tilde{A}}_{u,c}(k,m),$${\tilde{A}}_{u,s}(k,m)$, ${\tilde{B}}_{u,c}(k,m)$, ${\tilde{B}}_{u,s}(k,m)$ are all real-valued functions, to be determined.

We can now see that Equation (29a) determines the time-evolution of ${\widehat{w}}_{u,*}(k,z,t)$, a Fourier transform (in x) of the vertical velocity in the upper fluid, where $z\ge 0$. The vertical velocity field itself is necessarily real-valued, but its Fourier transform (in x) is usually complex-valued; it is purely real-valued only in special cases.

Comment 15.

From here on, ${\tilde{\alpha }}_{u,*}(k,m)$ and other functions with similar subscripts (such as ${\widehat{w}}_{u,*}(k,z,t)$) are used to represent either ${\tilde{\alpha }}_{u,c}(k,m)$ or ${\tilde{\alpha }}_{u,s}(k,m)$, in situations where either variable will serve.

Equations (48a) and (48b) accept solutions of either form, so (48a) and (48b) provide no a priori reason to choose one family of solutions over the other. The two families of solutions differ in terms of their boundary data at $z=0$:

–

A function ${\tilde{\alpha }}_{u,c}(k,z,t)$ with a representation of the form in (49b) can have nonzero values at $z=0$, but its first z-derivative must vanish there;

–

A function ${\tilde{\alpha }}_{u,s}(k,z,t)$ with a representation of the form in (50b) must vanish at $z=0$, but its first z-derivative can have nonzero values there.

–

A function ${\tilde{\alpha }}_{u,*}(k,z,t)$ with nonzero values at $z=0$ and also with nonzero values for its first derivative at $z=0$ necessarily has a representation that is a sum of the form in (49b), plus the form in (50b).

–

Similar statements hold for ${\tilde{\beta }}_{u,*}(k,z,t)$, for $z\ge 0$.

In addition to the two families of solutions of (46), found in (49a,49b,49c,49d) and (50a,50b,50c,50d), there are also limiting cases of (46), which also have physical significance. Recall from Section 2.2 that $\widehat{\overline{w}}(k,t)$ is defined to be the common limit of ${\widehat{w}}_u(k,z,t)$ and ${\widehat{w}}_l(k,z,t)$ as $z\to 0$, and recall from Section 2.5 that both ${\widehat{w}}_{*}(k,z,t)$ and ${\partial }_z{\widehat{w}}_{*}(k,z,t)$ are continuous in z at $z=0$, so apparently the functions in (49a) and (50a) represent two different parts of ${\widehat{w}}_u(k,z,t)$ in the upper fluid, and each must match up with the corresponding function in the lower fluid, in order for each of ${\widehat{w}}_{*}(k,z,t)$ and ${\partial }_z{\widehat{w}}_{*}(k,z,t)$ to be continuous in z at $z=0$.

For either (49a) or (50a), if ${\tilde{\omega }}_{u,c}(k,m)$ and ${\tilde{\omega }}_{u,s}(k,m)$ are both absolutely integrable in m on $-\infty <m<\infty$, then the Riemann-Lebesgue Lemma guarantees that ${\widehat{w}}_u(k,z,t)\to 0$ and ${\widehat{p}}_u(k,z,t)\to 0$ as $z\to \infty$, as required by (16a). Then (48a) and (48b) both show that the initial data in the upper fluid simply diffuse away in 2-D, with one important exception. The exceptional case is

$$\begin{array}{cc}\hfill {\tilde{w}}_u(k,z,t)& =C_u\left(k\right)·e^{-\left|k\right|z}.\hfill \end{array}$$

(51a)

In this case,

$$\begin{array}{cc}\hfill D_u{\tilde{w}}_u(k,z,t)& =0,{\partial }_t{\tilde{w}}_u(k,z,t)=0,\hfill \end{array}$$

(51b)

so the function in (51a) solves (46) trivially because each term in (46) vanishes; there is no diffusion. Please note that $C_u\left(k\right)$ is typically complex-valued.

Comment 16.

The result in (51a) can be interpreted as follows. For every fixed, real-valued $k\ne 0$, the energy in any localized fluid motion imposed at $t=0$ necessarily diffuses away for $t>0$, and the solution of (46) tends to zero as $t\to \infty$. (51a) is the limiting case of a solution of (46) as $t\to \infty$, for any $C_u\left(k\right)$.

As we show next, in Step 3, the imposed pressure term in (44) or equivalently in (52a) and (52b) below, provides a generalization of the function in (51a). The solution of (52a) and (52b) tends to the form given in (51a) as $t\to \infty$ , with the imposed pressure in (52a) and (52b) determining $C_u\left(k\right)$. Denote the generalized forms of $\{{\tilde{w}}_u(k,z,t),{\tilde{\alpha }}_u(k,z,t),{\tilde{\beta }}_u(k,z,t)\}$ respectively by $\{{\tilde{\tilde{w}}}_u(k,z,t),{\tilde{\tilde{\alpha }}}_u(k,z,t),{\tilde{\tilde{\beta }}}_u(k,z,t)\}$.

Step 3a in solving (44) proceeds by including the (known) right-hand side of (44):

$$\begin{array}{cc}\hfill {\partial }_t{\tilde{\tilde{w}}}_u+D_u{\tilde{\tilde{w}}}_u& =P_u(k,z,t),\hfill \end{array}$$

(52a)

or equivalently,

$$\begin{array}{cc}\hfill {\partial }_t{\tilde{\tilde{w}}}_u-{\nu }_u({\partial }_z^2{\tilde{\tilde{w}}}_u-k^2{\tilde{\tilde{w}}}_u)& =\frac{\left|k\right|}{{\rho }_u}·{\widehat{p}}_u(k,0,t)·e^{-\left|k\right|z},\hfill \end{array}$$

(52b)

where ${\tilde{\tilde{w}}}_u(k,z,t)$ denotes a solution to (52b). Please note that (52b) is a linear partial differential equation, with constant coefficients on the left side, and a known forcing function, which depends on $(k,z,t)$, on the right. An equation of this form necessarily has a family of homogeneous solutions, each of which ignores the forcing on the right-hand side, plus a family of particular solutions, each of which depends on the right-hand side.

Comment 17.

Solving (46) and then (52b) is another application of the method of Variation of Parameters, which was used in Section 2.3 to solve (18a, 18b) and then (19a, 19b).

Homogeneous solutions of (52b) are given in (49a) and (50a); particular solutions of (52b) can be found by the method of Variation of Parameters, as follows. In (49a) or (50a), the free “parameters” to be varied are embedded in ${\tilde{\omega }}_{u,c}(k,m)$ or ${\tilde{\omega }}_{u,s}(k,m)$, respectively. In these cases, Variation of Parameters works because additional time dependence is admitted into the problem, so ${\tilde{\omega }}_{u,c}(k,m)\to {\omega }_{u,c}(k,m,t)$ and ${\tilde{\omega }}_{u,s}(k,m)\to {\tilde{\omega }}_{u,s}(k,m,t)$. Therefore, the two families of homogeneous solutions of (52b) acquire slightly different forms:

$$\begin{array}{cc}\hfill {\tilde{\tilde{w}}}_{u,c}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\omega }_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm\hfill \end{array}$$

(53a)

and

$$\begin{array}{cc}\hfill {\tilde{\tilde{w}}}_{u,s}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\omega }_{u,s}(k,m,t)·sin\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm,\hfill \end{array}$$

(53b)

where ${\tilde{\tilde{w}}}_{u,c}(k,z,t)$ and ${\tilde{\tilde{w}}}_{u,s}(k,z,t)$ represent members of two families of solutions of (52b), with ${\omega }_{u,c}(k,m,t)$ and ${\omega }_{u,s}(k,m,t)$ being arbitrary complex-valued functions, to be determined. After substituting either version of (53a) or (53b) into (52b) and then simplifying, one obtains a reduced equation in one of two forms:

$$\begin{array}{c}\hfill {\int }_{-\infty }^{\infty }[{\partial }_t{\omega }_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}dm]=\frac{\left|k\right|}{{\rho }_u}·{\widehat{p}}_u(k,0,t)·e^{-\left|k\right|z}\end{array}$$

(54a)

or

$$\begin{array}{c}\hfill {\int }_{-\infty }^{\infty }[{\partial }_t{\omega }_{u,s}(k,m,t)·sin\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}dm]=\frac{\left|k\right|}{{\rho }_u}·{\widehat{p}}_u(k,0,t)·e^{-\left|k\right|z}.\end{array}$$

(54b)

Either version of (54a), (54b) simplifies if one writes $e^{-\left|k\right|z}$ in terms of the appropriate Fourier transform:

$$\begin{array}{c}\hfill {\int }_0^{\infty }[e^{-\left|k\right|z}·cos\left(mz\right)]dz=\frac{\left|k\right|}{k^2+m^2},\end{array}$$

(55a)

$$\begin{array}{c}\hfill {\int }_0^{\infty }[e^{-\left|k\right|z}·sin\left(mz\right)]dz=\frac{m}{k^2+m^2}.\end{array}$$

(55b)

Substituting these back into (54a), (54b) yields

$$\begin{array}{c}\hfill {\partial }_t{\omega }_{u,c}(k,m,t)=\frac{k^2}{{\rho }_u(k^2+m^2)}·{\widehat{p}}_u(k,0,t)·e^{{\nu }_u(k^2+m^2)t},\end{array}$$

(56a)

$$\begin{array}{c}\hfill {\partial }_t{\omega }_{u,s}(k,m,t)=\frac{\left|k\right|m}{{\rho }_u(k^2+m^2)}·{\widehat{p}}_u(k,0,t)·e^{{\nu }_u(k^2+m^2)t}.\end{array}$$

(56b)

Each equation in (56a) and (56b) can be integrated in t:

$$\begin{array}{c}\hfill {\omega }_{u,c}(k,m,t)-{\omega }_{u,c}(k,m,0)=\frac{k^2}{{\rho }_u(k^2+m^2)}·{\int }_0^t[{\widehat{p}}_u(k,0,\tau )·e^{{\nu }_u(k^2+m^2)\tau }]d\tau ,\end{array}$$

(57a)

$$\begin{array}{c}\hfill {\omega }_{u,s}(k,m,t)-{\omega }_{u,s}(k,m,0)=\frac{\left|k\right|m}{{\rho }_u(k^2+m^2)}·{\int }_0^t[{\widehat{p}}_u(k,0,\tau )·e^{{\nu }_u(k^2+m^2)\tau }]d\tau .\end{array}$$

(57b)

Please note that each equation in (57a) and (57b) also contains a term that represents the starting values at $t=0$: ${\omega }_{u,*}(k,m,0)$. This information is already encoded in (49a) and (50a); it must be added to the respective terms in (57a) and (57b). Substitute all of this into (49a) and (50a) to obtain

$$\begin{array}{c}\hfill {\tilde{\tilde{w}}}_{u,c}(k,z,t)={\int }_{-\infty }^{\infty }[{\tilde{\tilde{\omega }}}_{u,c}(k,m,0)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm\\ \hfill +{\int }_{-\infty }^{\infty }[\frac{k^2}{{\rho }_u(k^2+m^2)}·cos\left(mz\right)·{\int }_0^t\{{\widehat{p}}_u(k,0,\tau )·e^{-{\nu }_u(k^2+m^2)(t-\tau )}\}d\tau ]dm,\end{array}$$

(58a)

and

$$\begin{array}{c}\hfill {\tilde{\tilde{w}}}_{u,s}(k,z,t)={\int }_{-\infty }^{\infty }[{\tilde{\tilde{\omega }}}_{u,s}(k,m,0)·sin\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm\\ \hfill +{\int }_{-\infty }^{\infty }[\frac{\left|k\right|m}{{\rho }_u(k^2+m^2)}·sin\left(mz\right)·{\int }_0^t\{{\widehat{p}}_u(k,0,\tau )·e^{-{\nu }_u(k^2+m^2)(t-\tau )}\}d\tau ]dm.\end{array}$$

(58b)

One can verify that the functions in (58a) and (58b) each solve (52b), by substituting each into (52b). In addition, one can verify that the first integral in either function also satisfies (46), by substituting that integral into (46). At this point in the analysis, we have found formulae for two complex-valued functions, ${\tilde{\tilde{w}}}_{u,c}(k,z,t)$ and ${\tilde{\tilde{w}}}_{u,s}(k,z,t)$, each of which solves (52b), in part because a third complex-valued function, ${\widehat{p}}_u(k,0,t)$, is embedded in each. Keep in mind that ${\tilde{\tilde{w}}}_{u,c}(k,z,t)$ and ${\tilde{\tilde{w}}}_{u,s}(k,z,t)$ have different boundary data at $z=0$, so both functions are required to represent ${\tilde{\tilde{w}}}_{u,*}(k,z,t)$, the complete solution of (52b). The significance of (58a) and (58b) is that two physical processes are working simultaneously. The first integral in either equation represents diffusion of initial data, while the second integral in either equation represents external forcing due to pressure applied at the fixed interface (at $z=0$).

Step 3b involves one more generalization, of ${\tilde{\tilde{w}}}_{u,*}(k,z,t)$ in order to find exact solutions of (44), instead of only (52b). For this final generalization, we need to separate the real and imaginary parts of ${\tilde{\tilde{w}}}_{u,*}(k,z,t)$, so at this point we split each of (58a) and (58b) into their real and imaginary parts. Some of this splitting is already done, in (47), (49d) and (50d). Separately, the complex-valued pressure in the upper fluid, ${\widehat{p}}_u(k,z,t)$, was first defined in Section 2.3. Here in Section 3.1, only the boundary data, ${\widehat{p}}_u(k,0,t)$, are needed. ${\widehat{p}}_u(k,0,t)$ can be split into its real and imaginary parts: define ${\widehat{\gamma }}_u(k,0,t)$ and ${\widehat{\delta }}_u(k,0,t)$ to be real-valued functions such that

$$\begin{array}{c}\hfill {\widehat{p}}_u(k,0,t)={\widehat{\gamma }}_u(k,0,t)+i{\widehat{\delta }}_u(k,0,t).\end{array}$$

(59a)

Similarly, in (58a), set

$$\begin{array}{c}\hfill {\tilde{\tilde{\omega }}}_{u,c}(k,m,0)={\widehat{\tilde{A}}}_{u,c}(k,m,0)+i{\widehat{\tilde{B}}}_{u,c}(k,m,0),\end{array}$$

(59b)

where ${\widehat{\tilde{A}}}_{u,c}(k,m,0)$ and ${\widehat{\tilde{B}}}_{u,c}(k,m,0)$ are real-valued functions of real-valued parameters $\{k,m\}$, defined by the initial data of the problem. In addition, in (58b), set

$$\begin{array}{c}\hfill {\omega }_{u,s}(k,m,0)={\tilde{A}}_{u,s}(k,m,0)+i{\tilde{B}}_{u,s}(k,m,0).\end{array}$$

(59c)

The method of Variation of Parameters requires generalizations of (59b) and (59c), in order to provide more flexibility in the representation of ${\widehat{\omega }}_{u,c}(k,m,*)$. The desired generalization is

$$\begin{array}{c}\hfill {\widehat{\omega }}_{u,c}(k,m,t)={\widehat{A}}_{u,c}(k,m,t)+i{\widehat{B}}_{u,c}(k,m,t),\end{array}$$

(59d)

where ${\widehat{A}}_{u,c}(k,m,t)$ and ${\widehat{B}}_{u,c}$ are real-valued functions, defined so that the generalization of ${\tilde{\tilde{w}}}_{u,c}$ leads to the exact solutions of (44). Define ${\widehat{w}}_{u,c}(k,z,t)$ so that

$$\begin{array}{c}\hfill {\widehat{w}}_{u,c}(k,z,t)={\int }_{-\infty }^{\infty }[{\widehat{\omega }}_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm\\ \hfill +{\int }_{-\infty }^{\infty }[\frac{k^2}{{\rho }_u(k^2+m^2)}·cos\left(mz\right)·{\int }_0^t\{{\widehat{p}}_u(k,0,\tau )·e^{-{\nu }_u(k^2+m^2)(t-\tau )}\}d\tau ]dm.\end{array}$$

(60a)

Using (59a), the real part of ${\widehat{w}}_{u,c}(k,z,t)={\widehat{\alpha }}_{u,c}(k,z,t)+i{\widehat{\beta }}_{u,c}(k,z,t)$ is

$$\begin{array}{c}\hfill {\widehat{\alpha }}_{u,c}(k,z,t)={\int }_{-\infty }^{\infty }[{\widehat{A}}_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm\\ \hfill +{\int }_{-\infty }^{\infty }[\frac{k^2}{{\rho }_u(k^2+m^2)}·cos\left(mz\right){\int }_0^t\{{\widehat{\gamma }}_u(k,0,\tau )·e^{-{\nu }_u(k^2+m^2)(t-\tau )}\}d\tau ]dm,\end{array}$$

(60b)

while its imaginary part is

$$\begin{array}{c}\hfill i{\widehat{\beta }}_{u,c}(k,z,t)={\int }_{-\infty }^{\infty }[i{\widehat{B}}_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm\\ \hfill +{\int }_{-\infty }^{\infty }[\frac{k^2}{{\rho }_u(k^2+m^2)}·cos\left(mz\right){\int }_0^t\{i{\widehat{\delta }}_u(k,0,\tau )·e^{-{\nu }_u(k^2+m^2)(t-\tau )}\}d\tau ]dm.\end{array}$$

(60c)

Meanwhile, the real part of the generalization of (58b) is

$$\begin{array}{c}\hfill {\widehat{\alpha }}_{u,s}(k,z,t)={\int }_{-\infty }^{\infty }[{\widehat{A}}_{u,s}(k,m,t)·sin\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm\\ \hfill +{\int }_{-\infty }^{\infty }[\frac{\left|k\right|m}{{\rho }_u(k^2+m^2)}·sin\left(mz\right){\int }_0^t\{{\widehat{\gamma }}_u(k,0,\tau )·e^{-{\nu }_u(k^2+m^2)(t-\tau )}\}d\tau ]dm,\end{array}$$

(60d)

while the imaginary part of (58b) is

$$\begin{array}{c}\hfill i{\widehat{\beta }}_{u,s}(k,z,t)={\int }_{-\infty }^{\infty }[i{\widehat{B}}_{u,s}(k,m,t)·sin\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm\\ \hfill +{\int }_{-\infty }^{\infty }[\frac{\left|k\right|m}{{\rho }_u(k^2+m^2)}·sin\left(mz\right){\int }_0^t\{i{\widehat{\delta }}_u(k,0,\tau )·e^{-{\nu }_u(k^2+m^2)(t-\tau )}\}d\tau ]dm.\end{array}$$

(60e)

Step 4 brings the complex operator $i{L}_u(·)$, defined in (45c), back into the problem. This operator distinguishes between the real and imaginary parts of ${\widehat{w}}_{u,*}(k,z,t)$ so now it becomes essential to represent ${\widehat{w}}_{u,*}(k,z,t)$, a solution of (44), in terms of its real and imaginary parts. This was done in (47) for ${\tilde{w}}_u(k,z,t)$, which is the simplest nontrivial approximation of ${\widehat{w}}_{u,*}(k,z,t)$. We assume that a similar formula applies to ${\widehat{w}}_{u,*}(k,z,t)$:

$${\widehat{w}}_{u,*}(k,z,t)={\widehat{\alpha }}_{u,*}(k,z,t)+i{\widehat{\beta }}_{u,*}(k,z,t),$$

(61)

where ${\widehat{\alpha }}_{u,*}(k,z,t)$ and ${\widehat{\beta }}_{u,*}(k,z,t)$ are real-valued functions, chosen so that ${\widehat{w}}_{u,*}(k,z,t)$ in (61) is an exact solution of (44), according to the definitions given in (45a, 45b, 45c). [Note: All of the versions of (60a)–(60e) are consistent with (61).] With those definitions, (45a) guarantees that

$$D_u{\widehat{w}}_{u,*}(k,z,t)=-{\nu }_u({\partial }_z^2{\widehat{\alpha }}_{u,*}-k^2{\widehat{\alpha }}_{u,*})-i{\nu }_u({\partial }_z^2{\widehat{\beta }}_{u,*}-k^2{\widehat{\beta }}_{u,*}).$$

(62)

The complex-valued function ${\widehat{p}}_u(k,z,t)$, which first appears in (45b), is split into real and imaginary parts in (59a). Define the real-valued functions, ${\widehat{\gamma }}_u(k,z,t)$ and ${\widehat{\delta }}_u(k,z,t)$, so that (59a) holds. Then $P_u(k,z,t)$ in (45b) can be written as

$$P_u(k,z,t)=\frac{\left|k\right|}{{\rho }_u}·[{\widehat{\gamma }}_u(k,0,t)+i{\widehat{\delta }}_u(k,0,t)]·e^{-\left|k\right|z}.$$

(63)

Finally, inserting (61) into (45c) splits the integral term in (44) into its real and imaginary parts. The imaginary part of $i{L}_u{\widehat{w}}_u(k,z,t)$ appears when the operator acts on ${\widehat{\alpha }}_{u,*}(k,z,t)$, the real part of ${\widehat{w}}_{u,*}(k,z,t)$, and it generates a purely imaginary result:

$$\begin{array}{c}\mathrm{Im}\left\{i{L}_u{\widehat{w}}_{u,*}(k,z,t)\right\}=ik·U_u\left(\zeta (z,t)\right)·{\widehat{\alpha }}_{u,*}(k,z,t)\\ -2ik\left|k\right|·e^{-\left|k\right|z}{\int }_0^z[e^{2\left|k\right|q}{\int }_q^{\infty }\{{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{\alpha }}_{u,*}(k,s,t)·e^{-\left|k\right|s}\}ds]dq\\ +2ik·e^{\left|k\right|z}{\int }_z^{\infty }[{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{\alpha }}_{u,*}(k,s,t)·e^{-\left|k\right|s}]ds.\end{array}$$

(64a)

Meanwhile, the real part of $i{L}_u{\widehat{w}}_{u,*}(k,z,t)$ appears when the operator acts on $i{\widehat{\beta }}_{u,*}(k,z,t)$, the imaginary part of ${\widehat{w}}_{u,*}(k,z,t)$, and in generates a purely real result:

$$\begin{array}{c}\hfill \mathrm{Re}\left\{i{L}_u{\widehat{w}}_{u,*}(k,z,t)\right\}=-k·U_u\left(\zeta (z,t)\right)·{\widehat{\beta }}_{u,*}(k,z,t)\\ \hfill +2k\left|k\right|·e^{-\left|k\right|z}{\int }_0^z[e^{2\left|k\right|q}{\int }_q^{\infty }\{{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{\beta }}_{u,*}(k,s,t)·e^{-\left|k\right|s}\}ds]dq\\ \hfill -2k·e^{\left|k\right|z}{\int }_z^{\infty }[{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{\beta }}_{u,*}(k,s,t)·e^{-\left|k\right|s}]ds.\end{array}$$

(64b)

An explicit demonstration of the procedure in the upper fluid follows.

Recall from the discussion below (50d) that the vertical velocity in the upper fluid might have a Fourier cosine representation, or a Fourier sine representation, or a linear combination of the two. The boundary conditions at $z=0$ determine which representation is appropriate. For definiteness, we now assume that the boundary conditions require a Fourier cosine representation, but the analysis for either of the other options is quite similar.

Start with the formula in (60b), for the real part of the velocity field. Differentiate both sides of (60b) twice with respect to z, and observe that ${\nu }_u({\partial }_z^2-k^2){\widehat{\alpha }}_{u,c}(k,z,t)$ changes both integral terms by adding a multiplicative term of $-{\nu }_u(k^2+m^2)$ to the integrand of each integral. Then differentiate ((60b) with respect to time, ${\partial }_t{\widehat{\alpha }}_{u,c}(k,z,t)$, and observe that the t-derivative provides the same factor, $-{\nu }_u(k^2+m^2)$, within each integral. It follows that the difference between the two differentiations of (60b) is:

$$\begin{array}{cc}\hfill \{{\partial }_t-{\nu }_u({\partial }_z^2-k^2)\}{\widehat{\alpha }}_{u,c}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\partial }_t{\widehat{A}}_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm\hfill \\ \hfill +{\int }_{-\infty }^{\infty }[\frac{k^2}{{\rho }_u(k^2+m^2)}& ·cos\left(mz\right)·{\widehat{\gamma }}_u(k,0,t)]dm.\hfill \end{array}$$

(65a)

Please note that the left-hand side of (65a) is the differential operator that appears in (44), (46) and (52a).

Separately, note that the second integral on the right-hand side of (65a) can be written as

$$\begin{array}{cc}\hfill \frac{\left|k\right|}{{\rho }_u}·{\widehat{\gamma }}_u(k,0,t)& ·{\int }_{-\infty }^{\infty }[\frac{\left|k\right|}{k^2+m^2}·cos\left(mz\right)]dm.\hfill \end{array}$$

(65b)

Finally, comparing the integral in (65b) with the integral in (55a) and (55b) indicates that the integral in (65b) is simply $\left\{e^{-\left|k\right|z}\right\}$, written in terms of its inverse Fourier cosine transform. Assuming that the assertion below (65b) is correct, then (65a) can be written as

$$\begin{array}{cc}\hfill \{{\partial }_t-{\nu }_u({\partial }_z^2-k^2)\}{\widehat{\alpha }}_{u,c}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\partial }_t{\widehat{A}}_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm\hfill \\ \hfill +& \frac{\left|k\right|}{{\rho }_u}·{\widehat{\gamma }}_u(k,0,t)·e^{-\left|k\right|z}.\hfill \end{array}$$

(66)

Separately, Equation (44) can be split into its real and imaginary parts. The real part of (44) is

$$\{{\partial }_t-{\nu }_u({\partial }_z^2-k^2)\}{\widehat{\alpha }}_{u,*}(k,z,t)+\mathrm{Re}\left\{i{L}_u{\widehat{w}}_{u,*}\right\}=\frac{\left|k\right|}{{\rho }_u}·{\widehat{\gamma }}_u(k,0,t)·e^{-\left|k\right|z}.$$

(67)

This formula is correct for any of the three possible representations of the vertical velocity in the upper fluid, but if the boundary conditions at $z=0$ require a Fourier cosine representation, then (67) becomes

$$\{{\partial }_t-{\nu }_u({\partial }_z^2-k^2)\}{\widehat{\alpha }}_{u,c}(k,z,t)+\mathrm{Re}\left\{i{L}_u{\widehat{w}}_{u,c}(k,z,t)\right\}=\frac{\left|k\right|}{{\rho }_u}·{\widehat{\gamma }}_u(k,0,t)·e^{-\left|k\right|z}.$$

(68)

An important point here is that ${\widehat{\alpha }}_{u,c}(k,z,t)$ is the real part of ${\widehat{w}}_{u,c}(k,z,t)$, so the two functions necessarily have the same representation. Subtracting (68) from the cosine part of (66) leads to

$${\int }_{-\infty }^{\infty }[{\partial }_t{\widehat{A}}_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm=-\mathrm{Re}\left\{i{L}_u{\widehat{w}}_{u,c}(k,z,t)\right\}.$$

(69)

This result becomes more informative if we make use of the cosine part of (64b):

$$\begin{array}{c}{\int }_{-\infty }^{\infty }[{\partial }_t{\widehat{A}}_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm=k·U_u\left(\zeta (z,t)\right)·{\widehat{\beta }}_{u,c}(k,z,t)\\ -2k\left|k\right|·e^{-\left|k\right|z}{\int }_0^z[e^{2\left|k\right|q}{\int }_q^{\infty }\{{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{\beta }}_{u,c}(k,s,t)·e^{-\left|k\right|s}\}ds]dq\\ +2k·e^{\left|k\right|z}{\int }_0^z[{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{\beta }}_{u,c}(k,s,t)·e^{-\left|k\right|z}]ds.\end{array}$$

(70)

Recall from (60b) that ${\widehat{\alpha }}_{u,c}(k,z,t)$, the real part of ${\widehat{w}}_{u,c}(k,z,t)$, consists of two integrals, one of which is ${\int }_{-\infty }^{\infty }[{\widehat{A}}_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm$. It follows from (69) and (60b) that the first integral in (60b) is driven by $\mathrm{Re}\left\{i{L}_u{\widehat{w}}_{u,c}(k,z,t)\right\}$, while the second integral in (60b) is driven by the pressure term, $\frac{\left|k\right|}{{\rho }_u}·{\widehat{\gamma }}_u(k,0,t)·e^{-\left|k\right|z}$. Similar statements hold for $i{\widehat{\beta }}_{u,c}(k,z,t)$ in (60c), ${\widehat{\alpha }}_{u,s}(k,z,t)$ in (60d), and $i{\widehat{\beta }}_{u,s}(k,z,t)$ in (60e).

However, this is only part of the story. The imaginary part of ${\widehat{w}}_{u,c}(k,z,t)$ has its own representation, given in (60c). Following the same logic, but for $i{\widehat{\beta }}_{u,c}(k,z,t)$, one finds the analogue of (69):

$$i{\int }_{-\infty }^{\infty }[{\partial }_t{\widehat{B}}_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm=-\mathrm{Im}\left\{i{L}_u{\widehat{w}}_{u,*}(k,z,t)\right\}.$$

(71)

Using (64a), one finds that the analogue of (70) is

$$\begin{array}{c}{\int }_{-\infty }^{\infty }[{\partial }_t{\widehat{B}}_{u,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_u(k^2+m^2)t}]dm=-k·U_u\left(\zeta (z,t)\right)·{\widehat{\alpha }}_{u,c}(k,z,t)\\ +2k\left|k\right|·e^{-\left|k\right|z}{\int }_0^z[e^{2\left|k\right|q}{\int }_q^{\infty }\{{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{\alpha }}_{u,c}(k,s,t)·e^{-\left|k\right|s}\}ds]dq\\ -2k·e^{\left|k\right|z}{\int }_0^z[{\partial }_s{U}_u\left(\zeta (s,t)\right)·{\widehat{\alpha }}_{u,c}(k,s,t)·e^{-\left|k\right|z}]ds.\end{array}$$

(72)

3.2. Solve (29b), Also Known as (73)

Comment 18.

The steps taken to perform analysis on the lower fluid is nearly identical to those in the upper fluid. Thus, the details in this section are truncated, and refer commonly to setps in Section 3.1.

Step 1 for the lower fluid, with $z\le 0$, begins by writing (29b) in the following form:

$${\partial }_t{\widehat{w}}_l+D_l{\widehat{w}}_l+i{L}_l{\widehat{w}}_l=P_l(k,z,t),$$

(73)

where

$$\begin{array}{c}{D}_l{\widehat{w}}_l=-{\nu }_l({\partial }_z^2-k^2){\widehat{w}}_l\end{array}$$

(74a)

$$\begin{array}{c}{P}_l(k,z,t)=-\frac{\left|k\right|}{{\rho }_l}·{\widehat{p}}_l(k,0,t)·e^{\left|k\right|z},\end{array}$$

(74b)

and

$$\begin{array}{c}i{L}_l{\widehat{w}}_l(k,z,t)=ik·U_l\left({\zeta }_l(z,t)\right)·{\widehat{w}}_l(k,z,t)\\ +2ik\left|k\right|·e^{\left|k\right|z}{\int }_z^0[e^{-2\left|k\right|q}{\int }_{-\infty }^q\{{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}\}ds]dq\\ -2ik·e^{-\left|k\right|z}{\int }_{-\infty }^z[{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{w}}_l(k,s,t)·e^{\left|k\right|s}]ds.\end{array}$$

(74c)

In the same way that (46) was obtained by truncating (44) in the upper fluid, (75) can be obtained by truncating (73) in the lower fluid:

$${\partial }_t{\tilde{w}}_l+D_l{\tilde{w}}_l={\partial }_t{\tilde{w}}_l-{\nu }_l({\partial }_z^2-k^2){\tilde{w}}_l=0,$$

(75)

then by splitting ${\widehat{w}}_l$ into its real and imaginary parts:

$${\tilde{w}}_l(k,z,t)={\tilde{\alpha }}_l(k,z,t)+i{\tilde{\beta }}_l(k,z,t),$$

(76)

and finally introducing the differential operator, $D_l$:

$$\begin{array}{cc}\hfill {\partial }_t{\tilde{\alpha }}_l+D_l{\tilde{\alpha }}_l& ={\partial }_t{\tilde{\alpha }}_l-{\nu }_l({\partial }_z^2{\tilde{\alpha }}_l-k^2{\tilde{\alpha }}_l)=0,\hfill \end{array}$$

(77a)

$$\begin{array}{cc}\hfill {\partial }_t{\tilde{\beta }}_l+D_l{\tilde{\beta }}_l& ={\partial }_t{\tilde{\beta }}_l-{\nu }_l({\partial }_z^2{\tilde{\beta }}_l-k^2{\tilde{\beta }}_l)=0.\hfill \end{array}$$

(77b)

This is Step 2 in the lower fluid.

Next, we verify that (75) admits an exact solution in the form of a Fourier cosine transform:

$$\begin{array}{cc}\hfill {\tilde{w}}_{l,c}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\tilde{\omega }}_{l,c}(k,m)·cos\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm,\hfill \end{array}$$

(78a)

It follows that ${\tilde{\alpha }}_{l,c}(k,z,t)$ and $i{\tilde{\beta }}_{l,c}(k,z,t)$ each have representations in the form of (78a):

$$\begin{array}{cc}\hfill {\tilde{\alpha }}_{l,c}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\tilde{A}}_{l,c}(k,m)·cos\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm,\hfill \end{array}$$

(78b)

$$\begin{array}{cc}\hfill i{\tilde{\beta }}_{l,c}(k,z,t)& ={\int }_{-\infty }^{\infty }[i{\tilde{B}}_{l,c}(k,m)·cos\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm,\hfill \end{array}$$

(78c)

From (76):

$$\begin{array}{cc}\hfill {\tilde{\omega }}_{l,c}(k,m)& ={\tilde{A}}_{l,c}(k,m)+i{\tilde{B}}_{l,c}(k,m)\hfill \end{array}$$

(78d)

Note that (76) admits solutions in the form of a Fourier sine transform:

$$\begin{array}{cc}\hfill {\tilde{w}}_{l,s}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\tilde{\omega }}_{l,s}(k,m)·sin\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm,\hfill \end{array}$$

(79a)

Therefore, the real and imaginary parts of ${\tilde{w}}_{l,s}(k,z,t)$ each have solutions of (77a), (77b) of the same form:

$$\begin{array}{cc}\hfill {\tilde{\alpha }}_{l,s}(k,z,t)& ={\int }_{-\infty }^{\infty }[{\tilde{A}}_{l,s}(k,m)·sin\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm,\hfill \end{array}$$

(79b)

$$\begin{array}{cc}\hfill i{\tilde{\beta }}_{l,s}(k,z,t)& ={\int }_{-\infty }^{\infty }[i{\tilde{B}}_{l,s}(k,m)·sin\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm.\hfill \end{array}$$

(79c)

As in (78d),

$$\begin{array}{cc}\hfill {\tilde{\omega }}_{l,s}(k,m)& ={\tilde{A}}_{l,s}(k,m)+i{\tilde{B}}_{l,s}(k,m)\hfill \end{array}$$

(79d)

where ${\tilde{A}}_{l,c}(k,m),$${\tilde{A}}_{l,s}(k,m)$, ${\tilde{B}}_{l,c}(k,m)$, ${\tilde{B}}_{l,s}(k,m)$ are all real-valued functions, to be determined. Applying the same argument to the case where $z<0$, we see that if ${\tilde{\omega }}_{l,c}(k,m)$ and ${\tilde{\omega }}_{l,s}(k,m)$ are both absolutely integrable in m on $-\infty <m<\infty$, then ${\widehat{w}}_l(k,z,t)\to 0$ and ${\widehat{p}}_l(k,z,t)\to 0$ as $z\to -\infty$. (77a) and (77b) both show that the initial data in the lower fluid simply diffuse away in 2-D, with one exception:

$$\begin{array}{cc}\hfill {\tilde{w}}_l(k,z,t)e& =C_l\left(k\right)·e^{\left|k\right|z}.\hfill \end{array}$$

(80a)

In this case,

$$\begin{array}{cc}\hfill D_l{\tilde{w}}_l(k,z,t)& =0,{\partial }_t{\tilde{w}}_u(k,z,t)=0,\hfill \end{array}$$

(80b)

meaning that there is no diffusion. Please note that $C_l\left(k\right)$ is typically complex-valued.

Step 3a in the lower fluid follows the same logic as in the upper fluid. The imposed pressure term in (73) provides a generalization of the function in (80a), and the solution of (81b) tends to the form given in (80a) as $t\to \infty$ with the imposed pressure in (81b) determining $C_l\left(k\right)$. Denote the generalized forms of $\{{\tilde{w}}_l(k,z,t),{\tilde{\alpha }}_l(k,z,t),{\tilde{\beta }}_l(k,z,t)\}$ respectively by $\{{\tilde{\tilde{w}}}_l(k,z,t),{\tilde{\tilde{\alpha }}}_l(k,z,t),{\tilde{\tilde{\beta }}}_l(k,z,t)\}$.

Begin by including the known right-hand side of (73):

$$\begin{array}{c}\hfill {\partial }_t{\tilde{\tilde{w}}}_l+D_l{\tilde{\tilde{w}}}_l=P_l(k,z,t),\end{array}$$

(81a)

or equivalently,

$$\begin{array}{c}\hfill {\partial }_t{\tilde{\tilde{w}}}_l-{\nu }_l({\partial }_z^2{\tilde{\tilde{w}}}_l-k^2{\tilde{\tilde{w}}}_l)=-\frac{\left|k\right|}{{\rho }_l}{\widehat{p}}_l(k,0,t)·e^{\left|k\right|z},\end{array}$$

(81b)

where ${\tilde{\tilde{w}}}_l(k,z,t)$ denotes a solution to (81a) and (81b). Please note that an equation with the form of (81a) and (81b) has a family of homogeneous solutions, each of which ignores the forcing on the right-hand side, plus a family of particular solutions, each of which depends on the right-hand side.

Homogeneous solutions of (81b) are given in (78a) and (79a); particular solutions of (81b) can be found by the method of Variation of Parameters, as follows.

$$\begin{array}{c}\hfill {\tilde{\tilde{w}}}_{l,c}(k,z,t)={\int }_{-\infty }^{\infty }[{\omega }_{l,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm\end{array}$$

(82a)

and

$$\begin{array}{c}\hfill {\tilde{\tilde{w}}}_{l,s}(k,z,t)={\int }_{-\infty }^{\infty }[{\omega }_{l,s}(k,m,t)·sin\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm,\end{array}$$

(82b)

where ${\tilde{\tilde{w}}}_{l,c}(k,z,t)$ and ${\tilde{\tilde{w}}}_{l,s}(k,z,t)$ represent members of two families of solutions of (81b), with ${\omega }_{l,c}(k,z,t)$ and ${\omega }_{l,s}(k,z,t)$ being arbitrary, complex-valued functions, to be determined. After substituting either version of (82a),(82b) into (81b) and then simplifying, one obtains a reduced equation in one of two forms:

$$\begin{array}{c}\hfill {\int }_{-\infty }^{\infty }[{\partial }_t{\omega }_{l,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm=-\frac{\left|k\right|}{{\rho }_l}·{\widehat{p}}_l(k,0,t)·e^{\left|k\right|z}\end{array}$$

(83a)

or

$$\begin{array}{c}\hfill {\int }_{-\infty }^{\infty }[{\partial }_t{\omega }_{l,s}(k,m,t)·sin\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm=-\frac{\left|k\right|}{{\rho }_l}·{\widehat{p}}_l(k,0,t)·e^{\left|k\right|z}\end{array}$$

(83b)

Either version of (83a), (83b) simplifies if one writes $e^{\left|k\right|z}$ in terms of the appropriate Fourier transform:

$$\begin{array}{c}\hfill {\int }_{-\infty }^0[e^{\left|k\right|z}·cos\left(mz\right)]dz=\frac{\left|k\right|}{k^2+m^2},\end{array}$$

(84a)

$$\begin{array}{c}\hfill {\int }_{-\infty }^0[e^{\left|k\right|z}·sin\left(mz\right)]dz=\frac{m}{k^2+m^2}.\end{array}$$

(84b)

Substituting these back into (83a), (83b) yields

$$\begin{array}{c}\hfill {\partial }_t{\omega }_{l,c}(k,m,t)=-\frac{k^2}{{\rho }_l(k^2+m^2)}·{\widehat{p}}_l(k,0,t)·e^{{\nu }_l(k^2+m^2)t},\end{array}$$

(85a)

$$\begin{array}{c}\hfill {\partial }_t{\omega }_{l,s}(k,m,t)=-\frac{\left|k\right|m}{{\rho }_l(k^2+m^2)}·{\widehat{p}}_l(k,0,t)·e^{{\nu }_l(k^2+m^2)t}.\end{array}$$

(85b)

As with (56a) and (56b), each equation above can be integrated in t:

$$\begin{array}{c}\hfill {\omega }_{l,c}(k,m,t)-{\omega }_{l,c}(k,m,0)=-\frac{k^2}{{\rho }_l(k^2+m^2)}·{\int }_0^t[{\widehat{p}}_l(k,0,\tau )·e^{{\nu }_l(k^2+m^2)\tau }]d\tau ,\end{array}$$

(86a)

$$\begin{array}{c}\hfill {\omega }_{l,s}(k,m,t)-{\omega }_{l,s}(k,m,0)=-\frac{\left|k\right|m}{{\rho }_l(k^2+m^2)}·{\int }_0^t[{\widehat{p}}_l(k,0,\tau )·e^{{\nu }_l(k^2+m^2)\tau }]d\tau .\end{array}$$

(86b)

Substitute the information about the values at $t=0$, encoded in (86a),(86b), into (78a) and (79a) to obtain

$$\begin{array}{c}\hfill {\tilde{\tilde{w}}}_{l,c}(k,z,t)={\int }_{-\infty }^{\infty }[{\tilde{\tilde{\omega }}}_{l,c}(k,m,0)·cos\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm\\ \hfill -{\int }_{-\infty }^{\infty }[\frac{k^2}{{\rho }_l(k^2+m^2)}·cos\left(mz\right)·{\int }_0^t\{{\widehat{p}}_l(k,0,\tau )·e^{-{\nu }_l(k^2+m^2)(t-\tau )}\}d\tau ]dm,\end{array}$$

(87a)

and

$$\begin{array}{c}\hfill {\tilde{\tilde{w}}}_{l,s}(k,z,t)={\int }_{-\infty }^{\infty }[{\tilde{\tilde{\omega }}}_{l,s}(k,m,0)·sin\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm\\ \hfill -{\int }_{-\infty }^{\infty }[\frac{\left|k\right|m}{{\rho }_l(k^2+m^2)}·sin\left(mz\right)·{\int }_0^t\{{\widehat{p}}_l(k,0,\tau )·e^{-{\nu }_l(k^2+m^2)(t-\tau )}\}d\tau ]dm.\end{array}$$

(87b)

Just as with (58a) and (58b), one can verify that the functions above solve (81b), by substituting each into (81b). In addition, one can verify that the first integral in either function also satisfies (75), by substituting that integral into (75).

Step 3b in the lower fluid follows the same steps as in the upper fluid. First, we need to separate the real and imaginary parts of ${\tilde{\tilde{w}}}_{l,*}(k,z,t)$, so at this point we split each of (87a) and (87b) into their real and imaginary parts. Separately, the complex-valued pressure in the lower fluid, ${\widehat{p}}_l(k,z,t)$, was first defined in Section 2.3. Here in Section 3.2, only the boundary data, ${\widehat{p}}_l(k,0,t)$, are needed. ${\widehat{p}}_l(k,0,t)$ can be split into its real and imaginary parts: define ${\widehat{\gamma }}_l(k,0,t)$ and ${\widehat{\delta }}_l(k,0,t)$ to be real-valued functions such that

$$\begin{array}{c}\hfill {\widehat{p}}_l(k,0,t)={\widehat{\gamma }}_l(k,0,t)+i{\widehat{\delta }}_l(k,0,t).\end{array}$$

(88a)

Similarly, in (87b), set

$$\begin{array}{c}\hfill {\tilde{\tilde{\omega }}}_{l,c}(k,m,0)={\widehat{\tilde{A}}}_{l,c}(k,m,0)+i{\widehat{\tilde{B}}}_{l,c}(k,m,0),\end{array}$$

(88b)

where ${\widehat{\tilde{A}}}_{l,c}(k,m,0)$ and ${\widehat{\tilde{B}}}_{l,c}(k,m,0)$ are real-valued functions of real-valued parameters $\{k,m\}$, defined by the initial data of the problem.

The method of Variation of Parameters requires a generalization of (88b), in order to provide more flexibility to the representation of ${\tilde{\tilde{\omega }}}_{l,c}(k,m,*)$. The desired generalization is

$$\begin{array}{c}\hfill {\widehat{\omega }}_{l,c}(k,m,t)={\widehat{A}}_{l,c}(k,m,t)+i{\widehat{B}}_{l,c}(k,m,t),\end{array}$$

(88c)

where ${\widehat{A}}_{l,c}(k,m,t)$ and ${\widehat{B}}_{l,c}(k,m,t)$ are real-valued functions, defined so that the generalization of ${\tilde{\tilde{w}}}_{l,c}$ leads to exact solutions of (73).

Comment 19.

${\widehat{p}}_u(k,0,t)$ in (59a) and ${\widehat{p}}_l(k,0,t)$ in (88a) represent the boundary pressures on either side of the fixed interface, but these two pressures need not agree.

Define ${\widehat{w}}_{l,c}(k,z,t)$ so that

$$\begin{array}{c}\hfill {\widehat{w}}_{l,c}(k,z,t)={\int }_{-\infty }^{\infty }[{\widehat{\omega }}_{l,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm\\ \hfill -{\int }_{-\infty }^{\infty }[\frac{k^2}{{\rho }_l(k^2+m^2)}·cos\left(mz\right)·{\int }_0^t\{{\widehat{p}}_l(k,0,\tau )·e^{-{\nu }_l(k^2+m^2)(t-\tau )}\}d\tau ]dm.\end{array}$$

(89a)

Using (88a), the real part of ${\widehat{w}}_{l,c}(k,z,t)={\widehat{\alpha }}_{l,c}(k,z,t)+i{\widehat{\beta }}_{l,c}(k,z,t)$ is

$$\begin{array}{c}\hfill {\widehat{\alpha }}_{l,c}(k,z,t)={\int }_{-\infty }^{\infty }[{\widehat{A}}_{l,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm\\ \hfill -{\int }_{-\infty }^{\infty }[\frac{k^2}{{\rho }_l(k^2+m^2)}·cos\left(mz\right){\int }_0^t\{{\widehat{\gamma }}_l(k,0,\tau )·e^{-{\nu }_l(k^2+m^2)(t-\tau )}\}d\tau ]dm,\end{array}$$

(89b)

while its imaginary part is

$$\begin{array}{c}\hfill i{\widehat{\beta }}_{l,c}(k,z,t)={\int }_{-\infty }^{\infty }[i{\widehat{B}}_{l,c}(k,m,t)·cos\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm\\ \hfill -{\int }_{-\infty }^{\infty }[\frac{k^2}{{\rho }_l(k^2+m^2)}·cos\left(mz\right){\int }_0^t\{i{\widehat{\delta }}_l(k,0,\tau )·e^{-{\nu }_l(k^2+m^2)(t-\tau )}\}d\tau ]dm.\end{array}$$

(89c)

Meanwhile, the real part of the generalization of (87b) is

$$\begin{array}{c}\hfill {\widehat{\alpha }}_{l,s}(k,z,t)={\int }_{-\infty }^{\infty }[{\widehat{A}}_{l,s}(k,m,t)·sin\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm\\ \hfill -{\int }_{-\infty }^{\infty }[\frac{\left|k\right|m}{{\rho }_l(k^2+m^2)}·sin\left(mz\right){\int }_0^t\{{\widehat{\gamma }}_l(k,0,\tau )·e^{-{\nu }_l(k^2+m^2)(t-\tau )}\}d\tau ]dm,\end{array}$$

(89d)

while the imaginary part of (87b) is

$$\begin{array}{c}\hfill i{\widehat{\beta }}_{l,s}(k,z,t)={\int }_{-\infty }^{\infty }[i{\widehat{B}}_{l,s}(k,m,t)·sin\left(mz\right)·e^{-{\nu }_l(k^2+m^2)t}]dm\\ \hfill -{\int }_{-\infty }^{\infty }[\frac{\left|k\right|m}{{\rho }_l(k^2+m^2)}·sin\left(mz\right){\int }_0^t\{i{\widehat{\delta }}_l(k,0,\tau )·e^{-{\nu }_l(k^2+m^2)(t-\tau )}\}d\tau ]dm.\end{array}$$

(89e)

Just as $i{L}_u(·)$ was re-introduced in the equations for the upper fluid, $i{L}_l(·)$, which is defined in (74c), can be brought back into the problem. This is Step 4. It now becomes essential to represent ${\widehat{w}}_{l,*}(k,z,t)$, a solution of (73), in terms of its real and imaginary parts. This was done in (76) for ${\tilde{w}}_l(k,z,t)$. We assume that a similar formula applies to ${\widehat{w}}_{l,*}(k,z,t)$:

$${\widehat{w}}_{l,*}(k,z,t)={\widehat{\alpha }}_{l,*}(k,z,t)+i{\widehat{\beta }}_{l,*}(k,z,t),$$

(90)

where ${\widehat{\alpha }}_{l,*}(k,z,t)$ and ${\widehat{\beta }}_{l,*}(k,z,t)$ are real-valued functions, chosen so that ${\widehat{w}}_{l,*}(k,z,t)$ in (90) is an exact solution of (73), according to the definitions given in (74a)–(74c). [Note: All of the versions of (89a)–(89e) are consistent with (90).] With those definitions, (74a) guarantees that

$$D_l{\widehat{w}}_{l,*}(k,z,t)=-{\nu }_l({\partial }_z^2{\widehat{\alpha }}_{l,*}-k^2{\widehat{\alpha }}_{l,*})-i{\nu }_l({\partial }_z^2{\widehat{\beta }}_{l,*}-k^2{\widehat{\beta }}_{l,*}).$$

(91)

Define the real-valued functions, ${\widehat{\gamma }}_l(k,z,t)$ and ${\widehat{\delta }}_l(k,z,t)$, so that (88a) holds. Then $P_l(k,z,t)$ in (74b) can be written as

$$P_l(k,z,t)=-\frac{\left|k\right|}{{\rho }_l}·[{\widehat{\gamma }}_l(k,0,t)+i{\widehat{\delta }}_l(k,0,t)]·e^{\left|k\right|z}.$$

(92)

Following the same procedure as in (64a) and (64b), we show the results of applying the operator from (74c):

$$\begin{array}{c}\mathrm{Im}\left\{i{L}_l{\widehat{w}}_{l,*}(k,z,t)\right\}=ik·U_l\left({\zeta }_l(z,t)\right)·{\widehat{\alpha }}_{l,*}(k,z,t)\\ +2ik\left|k\right|·e^{\left|k\right|z}{\int }_z^0[e^{2\left|k\right|q}{\int }_{-\infty }^q\{{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{\alpha }}_{l,*}(k,s,t)·e^{\left|k\right|s}\}ds]dq\\ -2ik·e^{-\left|k\right|z}{\int }_z^{\infty }[{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{\alpha }}_{l,*}(k,s,t)·e^{\left|k\right|s}]ds.\end{array}$$

(93a)

$$\begin{array}{c}\mathrm{Re}\left\{i{L}_l{\widehat{w}}_{l,*}(k,z,t)\right\}=-k·U_l\left({\zeta }_l(z,t)\right)·{\widehat{\beta }}_{l,*}(k,z,t)\\ +2k\left|k\right|·e^{\left|k\right|z}{\int }_z^0[e^{-2\left|k\right|q}\{{\int }_{-\infty }^q[{\partial }_s{U}_l\left({\zeta }_l(s,t)\right)·{\widehat{\beta }}_{l,*}(k,s,t)·e^{\left|k\right|s}\}ds]dq\\ -2k·e^{\left|k\right|z}{\int }_z^{\infty }[{\partial }_s{U}_l\left(\zeta (s,t)\right)·{\widehat{\beta }}_{l,*}(k,s,t)·e^{\left|k\right|s}]ds.\end{array}$$

(93b)

Summary of Section 3.1 and Section 3.2:

• The two comparable analyses of the equations of motion in the upper (Section 3.1) and lower (Section 3.2) fluids simplify the relevant formulae for motion in the upper and lower fluids.
• One important difference between the dynamics in the upper fluid (Section 3.1) vs. that in the lower fluid (Section 3.2) arises if one were to follow the same set of steps outlined above for the equations describing the lower fluid. The second integral on the right side of the lower fluid analogy of (70) and the second integral on the right side of the lower fluid analogy of (72) are both much smaller (numerically) than the corresponding integrals shown for the upper fluid, because the horizontal background flow in the lower fluid is much smaller than that in the upper fluid, as shown in (Figure 1) on page 4. This important difference occurs because the density of air is smaller than that of water by a factor of about 800 (at room temperature). If the densities of the two fluids are more comparable, then the horizontal velocities of the two fluids are also more comparable, and the corresponding integrals in the upper and lower fluids are also more comparable.
• An important part of the structure of these equations is that ${\widehat{w}}_{*}(k,z,t)$ has a real and an imaginary part. Equations (69)–(72) are the only places in this entire analysis where the real and imaginary parts of $w_{u,*}(k,z,t)$ interact with each other. Similarly, the lower fluid analogues of (69)–(72) are the only places where the real and imaginary parts of $w_{l,*}(k,z,t)$ interact with each other.
• Section 3.1 and Section 3.2 describe respectively the dynamics of the fluid above the interface, in Section 3.1, and the dynamics of the fluid beneath the interface, in Section 3.2. In addition to the dynamics in each fluid, the two fluids interact with each other at or near the interface between the two fluids. These interactions are described in Section 2.5, and in more detail in Section 4, which follows.

4. Discussion

4.1. Summary of Results

The analysis carried out in the previous sections provides a mechanism with which to understand wind-driven waves. Beginning in Section 2.4, (29a) and (29b) describe the linearized motion in the upper fluid and lower fluid, respectively. (29a) and (29b) are the fundamental equations that show up in different variations throughout the analysis. Equation (44) is a direct consequence of (29a), with no extra or missing terms. The first two terms in (44) define the diffusion and the term on the right-hand side of (44) is the external pressure. The third term on the left-hand side of (44) describes the interaction between the real and imaginary parts of ${\widehat{w}}_u(k,z,t)$. Please note that while it will be omitted for most of the discussion, analogous processes apply to the lower fluid, defined by (73), a result of (29b).

Turning to the fixed interface at $z=0$, (32) and (33) govern the dynamics near the fixed interface. (33) is obtained by evaluating (29a) and (29b) at $z=0$, and then adding the two, which results in some terms cancelling. In order to maintain continuity in z, the vertical motion in the upper and lower fluids must match at the interface. The term $\widehat{\overline{w}}(k,t)$ is introduced to represent the limit as $z\to 0$ from either above or below. If the top line of (33) is set to 0, the simplified equation can be solved to obtain $\widehat{\overline{w}}(k,t)=\widehat{\overline{w}}(k,0)·e^{-ikU\left(0\right)t}$. The consequence is that the initial shape of the vertical velocity profile (i.e., $\widehat{\overline{w}}(k,0)$) propagates in the x-direction, at a speed approximately equal to $U\left(0\right)$. The shape of the vertical velocity profile does not change as it propagates in the x-direction, which is an important result. Please note that because ${\widehat{w}}_u(k,0,t)={\widehat{w}}_l(k,0,t)$, the vertical velocity must be continuous at the interface. Therefore the vertical velocities in the two fluids at $z=0$ coincide, and together they push the moveable interface either up or down, at the same rate.

The first line of (33), by itself, asserts that $\widehat{\overline{w}}(k,t)$ propagates horizontally, along the fixed interface (where $z=0$) with speed $U\left(0\right)$. These vertical velocities can be viewed as being measured at $z=0$, the location of the fixed interface. Meanwhile, (32) determines how $\widehat{\eta }(k,t)$ evolves in time. The solution of (32) is complex-valued, so it has two parts, an amplitude and a phase. The time-derivative term on the left side of (32) is used twice, once to determine how the phase changes as a function of time (based only on the left side of (32)), and the other use of the time-derivative term is to determine how the amplitude of $\widehat{\eta }(k,t)$ evolves in time, by bringing in the right-hand side of (32). Consider each of those two processes, which operate during the same time-interval. The phase-term, along with the two integral terms in (33), suggest that $\widehat{\eta }(k,t)$ propagates in approximately the same way as $\widehat{\overline{w}}(k,t)$ propagates: horizontally, and approximately along the fixed interface with speed $U\left(0\right)$. Therefore, $\widehat{\overline{w}}(k,t)$, on the right side of (32) controls the amplitude of $\widehat{\eta }(k,t)$, so $\widehat{\overline{w}}(k,t)$ controls the height of the movable interface relative to the fixed interface. In addition, if $\widehat{\overline{w}}(k,t)$ and $\widehat{\eta }(k,t)$ both travel horizontally at the same speed, and if $\widehat{\overline{w}}(k,t)>0$, then $\widehat{\eta }(k,t)$ continues to grow in amplitude as t increases.

4.2. White-Capping of Waves

Finally, this brings us to white-capping of waves in the presence of strong winds. At some initial time ($t=0$), the surface elevation, $\eta (x,t=0)$, has a Fourier transform, $\widehat{\eta }(k,t=0)$, which defines a fixed shape of the surface elevation, as a function of x in some spatial region. At the same time ($t=0$), the vertical motion in the upper and lower fluids, with some shape in each fluid in the x-direction, has its own Fourier transform at $z=0$, $\widehat{\overline{w}}(k,t)$. These two kinds of shapes at $t=0$ might have no relation to each other, so their Fourier transforms might also have no particular connection to each other. If we neglect the sum of the two integral terms in (33), then the vertical velocities near the fixed interface, in each of the two fluids, can push the location of the moveable interface up or down. As long as $\widehat{\overline{w}}(k,t)$ and $\widehat{\eta }(k,t)$ both travel in the same horizontal direction, at approximately the same speed, $U\left(0\right)$, then the height of the moveable interface can increase or decrease significantly over time.

Suppose the crest of a particular wave is pushed up significantly by the vertical velocities near $z=0$, in each fluid. As the crest of the wave is pushed higher in space, it is also pushed into regions of stronger horizontal velocities. (See Figure 1, on p. 4). As the crest of the wave is lifted into stronger winds, other parts of the free surface might not be pushed into such strong winds, so the existing winds can rip the crest of the wave away from the rest of the wave. This is white-capping of waves.