We study surface gravity waves in a homogeneous ocean where the depth is very much larger than the wavelength. These waves are affected by the physical conditions at the sea surface, most notably by the presence of thin flexible covers like biogenic films, which may cover large areas of the ocean surface (Gade et al., 2006). In addition to these natural films, we find pollutant slicks from petroleum spills or municipal effluents, which change the surface conditions. In cold regions, the surface cover is usually related to the presence of ice, for example, grease ice (Martin and Kauffman, 1981; Sutherland et al., 2019), or densely packed ice rubble in the marginal ice zone (Squire, 1984). Another example is the seasonal accumulation of Sargassum mats in the Tropical Atlantic (Marsh et al., 2022).

It is practically impossible to formulate a unified mathematical theory for explaining the effect on surface waves of floating material like oil, grease ice, or vegetation. We therefore intend to model this interaction by introducing a freely varying parameter R, which can be adapted to the problem in question. This parameter appears in a very general condition at the ocean surface. The condition is called a Robin condition (Gustafson, 1998; Akin, 2005) and R is the Robin parameter. Traditionally, the Robin condition is a relation between a quantity T (temperature, velocity, etc.) and its derivative ∂ T / ∂ x n , where x n is directed normal to the boundary. The homogeneous Robin condition can be written ∂ T / ∂ x n + R T = 0 at the boundary (sometimes R ^ ∂ T / ∂ x n + T = 0 is used). In this way, the Robin condition becomes a weighted combination of Dirichlet and Neumann boundary conditions. It is common in many branches of physics; see, for example, Hahn and Ozisk (2012) for applications to heat conduction, Tyvand and Nøland (2019) for convection in a porous medium, and Weber and Børve (2021) in the case of continental shelf waves with a permeable coastline. The relevance for oceanographic purposes is the role of the fluctuating tangential stress at the surface, which is typically small compared to the fluctuating normal stress (form stress) for uncontaminated surfaces, but can be large when the surface is covered by sea slicks, oil spills, sea ice, and other materials (e.g., Dorrestein, 1951). The framework presented here is mostly relevant for waves in deep water since the bottom friction and wave radiation stresses will typically dominate in shallow water. We will therefore only discuss waves for which k H ≫ 1 , where k is the wave number and H is the local water depth. We will also ignore rotational effects; hence, we assume that ω ≫ f , where ω   is the angular frequency of the waves and f is the Coriolis parameter.

In Figure 1 , we have depicted a general configuration relevant to the discussion here. The x-axis is horizontal and placed along the undisturbed sea surface, while the z-axis is positive upwards. The x-direction is, without loss of generality, in the direction of wave propagation. In the presence of waves, the sea surface is given by z = η ( x , y , t ) .

We shall here use the Robin formalism to describe the tangential conditions at the surface. Thus, we write

where u and w are the horizontal and vertical velocity components. In Equation 1, R is generally complex, i.e.,

We note from Equation 1 that when R = 0 , the viscous shear stress at the surface vanishes, defining a free surface. At the other extreme, when R → ∞ , the horizontal motion is zero at the surface, which models a flexible, but horizontally inextensible top layer (Lamb, 1932). The application of Equation 1 as a general model for the physical conditions at the surface appears to be novel. We will relate the Robin parameter to the rate of growth and decay in the waves, using known examples for elastic monolayers and thin layers of viscous fluids to aid the physical interpretation of our results.

The rest of this paper is organized as follows: In Section 2, we calculate the wave attenuation/growth by applying the Robin condition, and in Section 3, we discuss the solution when the Robin parameter is purely real. The case of a complex Robin parameter is considered in Section 4. Energy considerations due to dilational waves excited in the surface layer is discussed in Section 5. Finally, Section 6 contains some concluding remarks.

We here consider two-dimensional wave motion in a viscous ocean. For a fluid with viscosity ν, the two-dimensional velocity field can be separated into two parts: one potential part φ and one vorticity part ψ ; see Lamb (1932), such that

Inserting into the governing equations, we then obtain for the linear part of the wave field

Assuming that φ , ψ → 0 ,   z → − ∞ , the normalized linear solutions can be written

where C is a constant. In Equations 8 and 9, the complex wave number and frequency are defined by

Here α ,   β are the real (and small) spatial and temporal damping rates. From Equation 6, we obtain that m 2 = κ 2 − i n / ν   . If we now assume that the wavelength is much larger than the thickness of the oscillatory viscous boundary layer, and that the wave growth/damping is slow compared to the wave period, we have | m | ≫ | κ | and | ω ≫ | β | . Then, we have to leading order that

where γ is defined by

Here, the quantity γ − 1 represents the thin viscous boundary-layer thickness at the surface. To obtain Equation 12, we thus assume that k / γ ≪ 1. From Equation 1, we find that

Finally, assuming no fluctuating normal stress at the water surface, the dynamic boundary condition in the z-direction becomes

It is worth noting at this point that capillary waves can be included in the analysis by adding the normal stress due to the equilibrium value of the surface tension in Equation 15; see, for example, Lucassen (1968). The main conclusions in the present study will not change, however, and we ignore capillary waves in our examples here since their inclusion complicates the mathematical expressions without contributing to the physical interpretation of our results. From Equation 15, we find that

Inserting from Equations 8 and 9 into Equation 16, we obtain the dispersion relation for this problem:

where

Utilizing Equations 10–12, we obtain from the real part of Equation 17 to lowest order the obvious result ω 2 = g k . The imaginary part yields to lowest order that

where now c g = ω / ( 2 k ) is the usual deep-water wave group velocity. Here, ℜ denotes the real part of a complex quantity. We define a nondimensional Robin parameter by

and we scale β in Equation 19 by β i n e x = ω k / ( 4 γ ) , which is the temporal damping coefficient for an inextensible surface film (Lamb, 1932). Furthermore, we define β ^ = β / β i n e x and the nondimensional spatial damping rate by α ^ = 4 γ α / k 2 . Using that k / γ ≪ 1 ,   β / ω ≪ 1 , and α / k ≪ 1 , we obtain from Equation 19 that

where

In Equations 22 and 23, A , B ∈( − ∞ , ∞ ) , and it is easily seen that both P and Q are of order unity. When A = B = 0 , we have P = 2 and Q = 0 . Then Equation 21, in dimensional form, reduces to β + c g α = ω k 2 / γ 2 , as first obtained in Jenkins (1986) for a free surface. For larger values of A, B, we find that ( 2 k / γ ) P becomes negligible in Equation 20. A general form of Equation 19 for damped waves has been derived in Weber (2022), relating the right-hand side to the dissipation in the wave motion.

It is natural to start our discussion for the case when R is purely real, i.e., when B = 0 in Equations 22 and 23. In Figure 2 , we have plotted the attenuation rate Q vs. A = R r / γ in this case. We will assume that the waves are already present; hence, we ignore the generation mechanism and focus on the transient development of the existing wave field. In our analysis here, it is the fluctuating part of the tangential surface stress that is of interest, and hence, we will neglect the slowly varying wind stress that contributes to the upper ocean Ekman response, which is a second-order effect in the theoretical framework we apply here (e.g., Weber, 1983).

We note from Figure 2 that the depicted line tends to unity as A becomes infinitely large. Physically, this means that the horizontal velocity becomes zero at the contact with the surface cover, which is often referred to as the inextensible layer limit. Overall, real R is related to wave damping ( Q > 0 ) . However, it is interesting to note that Q < 0 (wave growth) occurs when − 2 ≤ A ≤ 0 . This may be explained in terms of the wave energy balance. From Equation 1, we may write at the surface that

where τ is the horizontal shear stress at the surface. Hence, by multiplying with the real u and average over the wave cycle (denoted by an over-bar), we find for the work per unit time of the shear stress on the fluid that

We note that the work is proportional to the square of the horizontal velocity component and independent of the vertical velocity. This is a general result that holds for any value of the Robin parameter. From Equation 25, we realize that to have a positive energy input, we must have R < 0. Furthermore, to have wave growth, this input must be larger than the viscous dissipation D in the fluid. Here,

Due to the strong gradients near the surface, ( ∂ u / ∂ z ) 2 ¯ will dominate in Equation 26. We derive real velocities from Equations 2, 3, 8, 9, and 14, and use that k / γ ≪ 1 . In Figure 3 , we have plotted the dimensionless quantities f 1 = 2 τ u ¯ / ( ρ ν k 2 γ ) and f 2 = 2 D / ( ρ ν k 2 γ ) as functions of A .

We note that in a small window, − 2 ≤ A ≤ 0 , the work by the surface stress is larger than the dissipation, so here we have wave growth. The maximum growth occurs for A = − 1 . At A = − 2 , the two effects balance, so here the growth rate is zero. These results confirm the findings for the attenuation rates depicted in Figure 2 .

Jenkins and Jacobs (1997) studied the wave damping by a thin layer of viscous fluid on top of an infinitely deep fluid of a different viscosity. For a thin, very viscous film of thickness d , they showed that the dynamic boundary condition in the horizontal direction between the two immiscible fluids could be written as

where μ s and μ are the dynamic viscosities of the thin upper surface layer and the infinitely deep bottom layer, respectively. Now, by comparison with Equation 1, we see that in this case

Since R is purely real, the damping in this case will follow the curve for A > 0   in Figure 2 . Here, we ignore any internal motion in the film itself since the analysis is only valid for thin films, with monomolecular films as the limit.

We now study the case R r = γ A ≥ 0 ,   R i = γ B ≥ 0 . In Figure 4 , we have plotted the attenuation rate Q from Equation 23 as a function of B for various values of A .

One classic (and striking) result in Figure 4 is that the maximum nondimensional attenuation rate becomes twice as large as inextensible limit when A = 0 (blue line), i.e., when the Robin parameter is purely imaginary and positive. This is related to the excitation of longitudinal elastic, or dilational, waves at the surface (i.e., Equation 9). Dilational waves in thin surface films efficiently dampen short waves—an intriguing phenomenon that has been the subject of study since antiquity (see, e.g., the historical review by Scott, 1977). From Miles (1967), generalizing an earlier result by Dorrestein (1951), we can write in our notation for an insoluble visco-elastic monolayer at the surface:

where E * is the surface elasticity per unit density, and ζ 1 and ζ 2 are the dilatational and shear viscosities of the film; see also Miles (1991). Comparing with Equation 1, we note that here

We thus see that the Robin boundary condition (Equation 1), with R r , R i > 0 , is related to the viscoelastic properties of the surface cover in that the real part of the Robin parameter represents the viscous effects (as we already have explored), while the imaginary part captures the elastic properties of the cover. This means that R_i , and hence B, must be positive since a negative value of the elasticity is unphysical. Miles (1991) was the first to point out that the coefficient in Equation 24 was a complex quantity.

In the film problem, we take that the surface viscosities are negligibly small (Dorrestein, 1951; Lucassen, 1968). Hence, the Robin parameter is purely imaginary. From Equation 23 we then have that

The maximum value is Q = 2 when B = 2 , as seen in Figure 4 . Introducing the dimensionless parameter δ = γ k 2 E * / ω 2 (Dorrestein, 1951), we have from Equation 32 that B ≡ R i / γ = 2 δ . Hence, we obtain from Equation 21 for pure temporal damping, when 2 k / γ is small, that β ^ m a x = Q m a x = 2 for δ = 1 , which is Dorrestein’s classic result.

Since the wave modes are coupled, Equation 21 yields the damping rate for both the dilational and surface waves.

Jenkins and Jacobs (1997) compare Equations 30 and 27 and conclude that Equation 30 has exactly the same properties as Equation 27 regarding wave damping if we replace ( ζ 1 + ζ 2 ) k 2 / ( ρ ν ) + i E * k 2 / ( ν ω ) in the film case by 4 d k 2 μ s / μ . However, as shown here, this is not true since the elastic film case has a dominant imaginary Robin part yielding increased damping ( Figure 4 , blue line), whereas in the viscous layer case, the Robin parameter is real.

As mentioned, the increased maximum damping for a certain R_i as seen in Figure 4 , is related to the existence of dilational waves (Lucassen, 1968; Weber and Christensen, 2003). Free dilational waves are critically damped through viscous dissipation in the oscillatory boundary layer. For this reason, the linear dispersion relation differs depending on whether the dilational waves are (i) damped in time, (ii) damped in space, or (iii) sustained by an undulating stress at the surface, i.e., forced waves (Weber and Christensen, 2003). Surface waves provide just such a mechanism for sustaining the dilational waves, and the dilational wave frequency relevant here is given by the dispersion relation for forced waves:

In the present problem, where maximum damping occurs for B = R i / γ = 2 , we find from Equation 34), by using Equations 13 and 32, that

Hence, as pointed out in Christensen (2005), maximum damping occurs when the natural frequency of the forced dilational waves exactly coincides with the natural frequency of the surface waves. Physically, the dilational waves, which would, in the absence of surface waves, be nearly critically damped, act as a sink of energy for short gravity waves. This phenomenon is sometimes referred to as “negative resonance”. Similar behaviors appear in more complex setups, for example, the case of finite layer thickness with two elastic interfaces, for which two damping maxima are possible (Ermakov and Khazanov, 2022).

We can gain more insight into the role of dilational waves as an energy sink if we split the horizontal velocity components into two parts u = u ˜ + u ^ , representing the surface wave (Equation 8) and the dilational wave (Equation 9), respectively. The presence of the dilational wave makes the amplitude of the total horizontal velocity u dependent on the elasticity parameter. Maximum damping is associated with a maximum in the horizontal velocity u , and hence a maximum in the work per unit time done by the surface stress (Equation 25). The dissipation is primarily due to the strongly sheared dilational wave motion:

The surface stress is well approximated by τ = ρ υ ∂ u ^ / ∂ z ( z = 0 ) . We now define the work per unit time on the oscillatory boundary layer:

The second term on the right-hand side of Equation 36 is zero because of u ^ and ∂ 2 u ^ / ∂ z 2 always being 90° out of phase, cf. Equations 9 and 13. From Equations 36 and 37, we find that W ^ = D ^ as we would expect for dilational waves (Weber and Christensen, 2003), which means that although the dilational waves are continuously excited by the surface waves, being coupled through the surface condition (Equation 1), they are continuously being dissipated in the oscillatory boundary layer at the same rate.

To model the various effects that may change the physics of the ocean surface, and hence affect growth and decay of surface waves, we have related the shear stress at the water surface and the horizontal velocity using a Robin condition with a complex Robin parameter. Using known results for viscous fluid layers and elastic monolayers, we show how the real part of the Robin parameter is related to viscous or frictional surface effects, while the imaginary part is associated with surface elasticity. In many situations of practical interest, for instance for ice covered waters, large floating mats of seaweeds, or more complex sea slicks with anisotropic horizontal distributions of surfactant, the dynamical conditions at the surface quickly become very complex, and we need to resort to simplified parameterizations of average behavior, particularly on the horizontal scales of numerical wave prediction models. Since it is practically impossible to formulate a unified mathematical theory explaining the effect of surface constituents on surface waves, or even to formulate constitutive relations, we may resort to modeling the interaction between the surface waves and the floating materials by applying the Robin condition and adapting the Robin parameter to the problem in question. Even if a quantitative estimation of the Robin parameter is not possible, a discussion of the problem for limiting values of R may still improve our understanding of the basic physics involved.

The original contributions presented in the study are included in the article/supplementary material. Further inquiries can be directed to the corresponding author.

JW contributed the original idea behind this study. JW and KC contributed equally to the mathematical analysis and writing. Both authors contributed to the article and approved the submitted version.