A full-field perturbation approach is modified for an ice-covered ocean and applied to estimating narrowband long-range reverberation caused by roughness of the ice-water interface. First-order approximation of the approach is used which requires the roughness amplitudes be small compared to the acoustic wavelength. To obtain the zeroth-order Green's function and transmission loss field used in the reverberation model, elastic parabolic equation solutions are generated in range-independent environments. Ice is represented by an isospeed layer on top of a linear transition layer. Effects of ice properties are discussed and demonstrated by comparing reverberation calculated for different ice layer thicknesses and wave speeds for typical ice values.

## 1. Introduction

In order to analyze reverberation in complex ice-covered environments, the latest models of scattering are used in conjunction with recently benchmarked full-field parabolic equation solutions for propagation loss. These solutions are well known to be accurate and efficient for several types of environments, in particular, all fluid environments^{1} and environments where the ice layer and ocean bottom can be treated as elastic media.^{2}

A fundamental property of parabolic equation solutions is the assumption that the backscattered acoustic field is small enough to be ignored, which makes analysis of reverberation due to backscattering from rough interfaces difficult. This paper addresses the problem by applying an analytic full-field perturbation approach^{3} with the aim to demonstrate a practical way of calculating the reverberation in Arctic environments caused by a small (compared to acoustic wavelength) roughness on the interface between the ice layer and water column.

Single scattering results based on fluid model parabolic equation (FPE) solutions, where the ice is treated as a fluid layer, will be compared with those based on elastic model parabolic equation (EPE) solutions that include an elastic ice layer. The EPE solutions for range-independent waveguides have been recently benchmarked against wavenumber integration and normal mode methods.^{2} Combining these solutions with a perturbation approach described in Ref. 3 allows one to estimate the impact of varying characteristics of the ice layer, such as thickness and elastic material parameters, on long-range reverberation at low enough frequencies where the approach is applicable. In particular, this means large features such as ice keels are not accounted for in this paper, assuming this is a separate subject for future development. This could combine our approach with that used in previous work where reverberation due to ice keels is a focus, see, e.g., Refs. 4–6 and references therein.

## 2. A first-order reverberation model

We use a first-order full-field perturbation approach to reverberation that considers the total field (acoustic pressure) *p* as a sum of the unperturbed field *p*_{0} and a first order (single scattered) field *p*_{1}. This total pressure field is expressed through corresponding components of the full-field Green's function

with a coefficient *D* introduced so that the Green's function of the unperturbed medium near the source represents a spherical wave of unit magnitude, $\u2009|r\u2212r\u2032|G0(r,r\u2032)\u2192r\u2192r\u20321$. Analytical results and numerical solutions for zeroth-order Green's functions are currently available for complicated background (unperturbed) media, such as range-dependent shallow water or deep-water waveguides, and are based on various propagation models, such as parabolic equation, normal modes, or ray approximations. Then an integral expression for average intensity of the first order (single-scattered) field near the source (i.e., for monostatic configuration) $I1=\u27e8|p1|2\u27e9$, can be derived as follows:

where integration is performed over the unperturbed interface, and $MR$ is a local full-field scattering coefficient attributed to this interface. While full details for scattering from a rough interface between two fluid media appear elsewhere,^{3} we give a specific expression for this coefficient for monostatic geometry which can be presented as follows:

where $\Gamma $ is a contrast factor at the interface, and $mj=\rho j/\rho 0$ and $nj=c0/cj$ are relative density and refraction index of the background medium near the interface between fluid medium 1 (water) and fluid medium 2 (called here “fluid ice” or “slush”). The parameters $\rho 0,\u2009c0,\u2009\u2009and\u2009\u2009k0$ are reference values for density, sound speed, and wavenumber, respectively, and $\Phi (r\u22a5,q)$ is the local power spectrum given by

For monostatic configurations the spectrum is taken at the Bragg backscattering wave vector $q=2k$. The magnitude of the spatially varying roughness height $\zeta (r\u22a5)$ is assumed to be small compared to the wavelength so that perturbation theory can be applied. More details can be found in Ref. 3.

In addition to the coefficient $MR$, as is seen from Eq. (2), only the magnitude of the Green's function near the interface is needed, which is related to transmission loss (TL) in the unperturbed medium, $TL=20\u2009log10|r0G0|$, where $r0=1$ m is the unit reference distance, so that $TL(r0)=0$ dB. This approach allows for a straightforward calculation of reverberation in those environments where TL estimation methods are available. In this paper, a full-field computational solution obtained from an elastic parabolic equation method is used to obtain TL at the unperturbed interface between the ice and water. Equation (2) is then used to generate the average backscattered intensity caused by perturbations of this interface. When combining this approach with EPE solutions, we must include a thin transition elastic-to-fluid slush layer (detailed in Sec. 3) between elastic ice and the water column since $MR$ was derived for a fluid-fluid interface.

## 3. Parabolic equation solutions

To utilize the reverberation estimates, a full-field numerical solution along the rough interface is more appropriate than methods that only include discrete portions of the acoustic field. The derivation and implementation of parabolic equation solutions are presented in the literature, but this section will introduce the fundamental approaches for FPE and EPE solutions and discuss the ice layer model that allows application of the scattering coefficient in Eqs. (2) and (3).

Underwater acoustic environments with ice cover can be modeled with a fluid ice layer and a fluid sea bottom. The FPE applies the parabolic equation method to obtain an approximation to the full field Helmholtz equation, and is well known to be accurate and efficient for underwater acoustic calculations. In a range independent segment, the azimuthally symmetric wave equation for acoustic pressure $p$ is factored into operators that describe the forward and backward propagating waves. The parabolic approximation ignores the backscattered wave to obtain

where $X$ is a matrix that includes depth-dependent derivatives.^{1} This differential equation is solved by marching from an initial starting field. The $1+X$ operator is approximated using Padé rational functions.^{1}

A key aspect of this approach is the rigorous application of fluid-fluid boundary conditions. For ice covered environments, this means continuity of pressure and vertical particle velocity are enforced at the boundary between the fluid ice layer and the ocean, as well as between the ocean and the ocean bottom. One shortcoming of the FPE in this setting is the application of a pressure release boundary condition at the top of the ice layer, where the physics of elastic media requires a zero traction condition.

For EPE solutions, the acoustic field is determined by the elastic equations of motion in the $(ur,w)$ formulation of elasticity, where $ur=\u2202u/\u2202r$, $u$ is horizontal displacement, and $w$ is vertical displacement. These incorporate the compressional wave speeds $cj$, shear wave speeds $cs,j$, and densities $\rho j$, via the Lamé parameters for each elastic layer *j* = 2, 3 where *j* = 2 is the elastic ice layer and *j* = 3 is the elastic ocean bottom.

The system of elliptic equations that arises is expressed as a vector equation, and factored into a product of parabolic operators that represent outgoing and incoming energy. By applying the assumption that backscattered energy is negligible, the outgoing factor is retained to give

where *L* and *M* are depth-dependent matrix operators^{7} that are approximated by a rational-linear Padé series.

Equation (7) is solved using a marching scheme in range, given an initial condition that can include a source in the water or in one of the elastic layers.^{8} The EPE solutions generated here have been benchmarked against normal mode and wavenumber integration solutions in ice-covered environments. They rigorously apply fluid-elastic boundary conditions at the ice-ocean and ocean-bottom interfaces, as well as the physically accurate zero traction condition at the top of the ice layer.^{2}

Figure 1 shows parabolic equation solutions for a 500 Hz source in a four-layered ice covered environment: (1) The top layer is a 2 m thick homogeneous elastic ice layer with $cp$ = 3200 m/s, $cs=1600$ m/s, $\rho $ = 0.9 g/cm^{3}, compressional attenuation $\alpha p=0.5$ dB/$\lambda ,\u2009$ and shear attenuation $\alpha s=1.5\u2009$ dB/$\lambda $. (2) Underneath the homogeneous elastic ice layer the transition elastic-to-fluid slush layer is comprised of two parts. First is a 0.5 m thick elastic transition layer where $cp$ and $cs\u2009$ change linearly to 1500 and 400 m/s, respectively, and the other parameters are kept the same. Underneath the elastic transition layer we also introduce a thin 0.1 m fluid layer where *c _{p}* = 1500 and $cs=0$ m/s. The small contrast between this thin “fluid slush” layer parameters and “elastic slush” parameters above it allows ignoring the weak reflections from the fluid-elastic interface within the ice-slush layer, while the fluid nature of the interface between this layer and water allows the application of Eq. (2) for estimating reverberation levels.

^{9}The inclusion of these layers is akin to the observed presence of a platelet ice layer in certain environments.

^{10}(3) The third layer is an ocean layer with a canonical Arctic sound speed profile,

^{11}shown in Fig. 1(a), and density

*ρ*= 1.0 g/cm

^{3}. No absorption is included in the water layer. (4) Finally an elastic ocean bottom is considered with $cp$ = 1600 m/s, $cs=500\u2009$ m/s, $\rho $ = 1.0 g/cm

^{3}, $\alpha p$ = 1.0 dB/$\lambda $, and $\alpha s$ = 1.5 dB/$\lambda $. The bottom density and sounds speed were chosen to minimize reflections into the waveguide so they do not overshadow effects of the ice layer.

The source is located 30 m below the water-ice interface. The depth variable begins at $z$ = 0 m at the top of the ice layer so the actual depth of the source is 32 m. Figure 1(a) shows the EPE solution. There is clearly trapped acoustic energy in the upper portion of water column that has multiple interactions with the ice layer caused by the upward refracting profile. The lower panels of Fig. 1 show TL contour plots for the top 10 m of the environment described above. The ice layer profile, the linear transition to water layer parameters, and top 8 m of the water sound speed profile is shown in Fig. 1(c). Figures 1(d) and 1(e) show the FPE and EPE solutions, respectively. In the FPE the shear effects are neglected by assigning $cs=0$, while all other parameters are the same as in the EPE. The FPE solution has an overall higher TL in the ice layer when compared with the EPE solution for the entire propagation range. This is likely due to the transition of acoustic energy into both compressional and shear wave energy as acoustic waves interact with the ice layer. It is also attributed to reflection from the physically accurate zero traction boundary condition at the air-ice interface. The EPE solution also demonstrates interface wave propagation along the ice-ocean interface, consistent with the physics of elastic materials.

## 4. Parameter analysis for reverberation

This section shows results for reverberation levels in terms of normalized (dimensionless) reverberation intensity, $I1r02/\u2009\u2009|D|2$, calculated using Eq. (2) for omnidirectional source and receiver collocated 30 m below the water-ice interface in the 10 km range independent environment used above. In the case of a directional source and/or receiver, these levels would only need an additional scaling factor when applying Eq. (2). The normalized intensity is shown as a function of time by assuming $t=2r/c0$, reasonable for narrowband signals. The range span for integration in Eq. (2) is defined by the pulse duration $\tau $, i.e., $r\u2208(c0t/2,c0t/2+c0\tau /2)$. For narrowband signals with the frequency bandwidth $\Delta f$, it is assumed that $\tau =1/\Delta f$. In calculations below, we use $\tau =0.5$ s, which corresponds to smoothing over a 375 m range span. Note that to keep the focus on the effects of the ice layer thickness and elasticity we consider a fixed input value for the dimensionless full-field scattering coefficient, $MR=0.0035$, the same at the two different frequencies in the following numerical examples. Equation (3) provides the way to consider the possible dependence on frequency of *M _{R}* for any specified roughness spectrum and its parameters, root-mean-square roughness, correlation scale, and others.

Figure 2 compares dimensionless reverberation levels obtained from FPE and EPE solutions to analyze model differences and demonstrate the impact of varying the ice layer thickness on the reverberation at two frequencies, 500 and 2000 Hz. In all four panels the 2000 Hz curves are offset by 10 dB for visual separation. Solutions are generated for the environment described for Fig. 1 above except where noted. Specifically, Fig. 2 shows a comparison between reverberation levels from solutions for 1 m thick ice with a 0.25 m transition layer and those for 2 m thick ice with a 0.5 m transition layer. Material parameters of the ice layer are the same as used in Fig. 1.

For each frequency, as seen in Figs. 2(a) and 2(b), the EPE solution shows a generally insignificant effect of the ice thickness on the overall level, while for the FPE solution the difference in levels at the two frequencies is notable. This is because the level of FPE solution is affected by compressional waves reflected from the ice-air interface whose contribution is significantly reduced at the higher frequency. For EPE solution this reduction effect may be compensated by the contribution of shear waves. Besides, the effect on a “finer” behavior of the time-dependence for the EPE solution is stronger and shows notable differences between the curves near 6 and 10 s, which may correspond to convergence zones as seen in Fig. 1. The 2000 Hz comparisons of the EPE curves have notably different spikes in the reverberation near 1 s, indicating that the ice thickness has an effect on the initial (at short ranges) backscatter from the interface. The EPE curves seem very similar after that, until an apparent phase reversal begins near 10 s, which is not seen for FPE curves. To emphasize the above-mentioned effects, a comparison between EPE and FPE model results on one plot for both frequencies are shown separately when the ice layer thickness is 1 m, Fig. 2(c), and 2 m, Fig. 2(d), respectively. The EPE solutions show clearly less reverberation loss. In addition to the effect discussed above, this can likely be caused by decreased TL levels observed at the ice-ocean interface in Fig. 1 due to the apparent interface wave. The initial spike that occurs in the 2000 Hz curve at 1 s is not reproduced by the FPE. Figure 2(d) shows the same comparison for a 2 m thick ice layer, and leads to consistent observations regarding overall level, wider reverberation, and the reproduction of initial backscatter.

Figure 3 compares fluid and elastic model solutions obtained at 500 and 2000 Hz to show the effect of varying ice layer material parameters for the same environment as described in Fig. 1. The ice layer parameters are controlled in order to establish contributions of compressional or shear speed to reverberation level. Figure 3(a) shows comparisons of FPE solutions for a 2 m ice layer where the compressional speed in the ice layer is either 3200 or 3600 m/s. This typical range of ice layer compressional speed does not result in observable changes in the reverberation. Elastic model comparisons for 2 m ice layers are shown in Fig. 3(b) where the shear speed is held constant at 1600 m/s, but the compressional speed is changed from 3200 to 3600 m/s. This results in a difference in the overall reverberation level, but almost no change in the shape of the curve. This level difference is also observed at several other compressional speed values but not shown here. Similar results are shown in Fig. 3(c), where the compressional speed is held constant at 3200 m/s and the shear wave speed is varied between 1400 and 1600 m/s. The level changes resulting from variations of the elastic wave parameters could result from changes to the amplitude of the interface wave observed in Fig. 1. It is important to note these results are consistent even though one *c _{s}* value is less than the water sound speed at the interface and one value is greater than it. The lower panels of Fig. 3 show the same comparisons as above for a 1 m ice layer, with nearly identical results.

## 5. Concluding remarks

For the canonical environments and frequencies presented here, 500 Hz and 2 kHz, there are notable differences between the acoustic TL for EPE and FPE solutions when the ice layer is present, and these differences contribute to variations in acoustic TL near the ice-water interface. Correspondingly, the differences result in variations in the intensity of reverberation caused by scattering from roughness of this interface. In particular, the EPE solutions show stronger effects on a finer behavior of reverberation and related additional capabilities for discrimination between different ice layer thicknesses, showing contrast in the initial arrivals and curve characteristics at the two frequencies while the FPE solutions only show changes in the overall reverberation level. It thus appears critical to account for excitation of shear waves in the ice layer. It is important, in particular, to accurately represent the zero traction boundary condition applied in the EPE models and excitation of an interface wave at the ice-water interface.

Future work should determine a fluid-elastic scattering coefficient. The reverberation integral is expected to have the same form since it was obtained for an arbitrary rough surface, and only the full field scattering coefficient needs to be modified for a more physically accurate reverberation estimate. Important future considerations would include effects of volume heterogeneity of ice as well, based on previously developed approaches,^{12–14} and also comparisons of the reverberation estimates to data from recent underwater acoustic experiments performed in the Arctic Ocean.

## Acknowledgments

Work supported by grants from the Office of Naval Research to Marist College and University of Washington and by a collaboration grant from the Simons Foundation (Grant No. 279472 to S.D.F.). The authors thank the anonymous reviewers for valuable comments.