Research article 24 Jun 2020
Research article  24 Jun 2020
Spectral attenuation of ocean waves in pack ice and its application in calibrating viscoelastic waveinice models
 ^{1}Nansen Environmental and Remote Sensing Center, Bergen, Norway
 ^{2}Department of Ocean and Resources Engineering, University of Hawaii, Mānoa, HI, USA
 ^{3}Univ. Brest, CNRS, IRD, Ifremer, Laboratoire d'Océanographie Physique et Spatiale (LOPS), IUEM, Brest, France
 ^{4}Department of Civil and Environmental Engineering, Clarkson University, Potsdam, NY, USA
 ^{1}Nansen Environmental and Remote Sensing Center, Bergen, Norway
 ^{2}Department of Ocean and Resources Engineering, University of Hawaii, Mānoa, HI, USA
 ^{3}Univ. Brest, CNRS, IRD, Ifremer, Laboratoire d'Océanographie Physique et Spatiale (LOPS), IUEM, Brest, France
 ^{4}Department of Civil and Environmental Engineering, Clarkson University, Potsdam, NY, USA
Correspondence: Sukun Cheng (sukun.cheng@nersc.no)
Hide author detailsCorrespondence: Sukun Cheng (sukun.cheng@nersc.no)
We investigate a case of ocean waves through a pack ice
cover captured by Sentinel1A synthetic aperture radar (SAR) on 12 October 2015 in the Beaufort Sea. The study domain is 400 km by 300 km, adjacent to a
marginal ice zone (MIZ). The wave spectra in this domain were reported in a
previous study (Stopa et al., 2018b). In that study, the authors divided the
domain into two regions delineated by the first appearance of leads (FAL)
and reported a clear change of wave attenuation of the total energy between
the two regions. In the present study, we use the same dataset to study the
spectral attenuation in the domain. According to the quality of SARretrieved wave spectrum, we focus on a range of wave numbers corresponding to
9–15 s waves from the openwater dispersion relation. We
first determine the apparent attenuation rates of each wave number by pairing
the wave spectra from different locations. These attenuation rates slightly
increase with increasing wave number before the FAL and become lower and more
uniform against wave number in thicker ice after the FAL. The spectral
attenuation due to the ice effect is then extracted from the measured
apparent attenuation and used to calibrate two viscoelastic waveinice
models. For the Wang and Shen (2010b) model, the calibrated equivalent shear
modulus and viscosity of the pack ice are roughly 1 order of magnitude
greater than that in grease and pancake ice reported in Cheng et al. (2017).
These parameters obtained for the extended Fox and Squire model are much
greater, as found in Mosig et al. (2015) using data from the Antarctic MIZ.
This study shows a promising way of using remotesensing data with large
spatial coverage to conduct model calibration for various types of ice
cover.
Highlights. Three key points:

The spatial distribution of wave number and spectral attenuation in pack ice are analyzed from SARretrieved surface wave spectra.

The spectral attenuation rate of 9–15 s waves varies around 10^{−5} m^{2} s^{−1}, with lower values in thicker semicontinuous ice fields with leads.

The calibrated viscoelastic parameters are greater than those found in pancake ice.
 Article
(4253 KB) 
Supplement
(500 KB)  BibTeX
 EndNote
Rapid reduction in Arctic ice in recent decades has become a focal point in climate change discussions (Comiso et al., 2008; Meier, 2017; Rosenblum and Eisenman, 2017; Stroeve and Notz, 2018). The reduction emphasizes the need to better understand the complex interaction between the sea ice, the ocean and the atmosphere. One of these interaction processes is between ocean waves and sea ice. Ocean waves help to shape the formation of new ice covers (Lange et al., 1989; Shen et al., 2001), break existing ice covers (Kohout et al., 2016), modify upperocean mixing (Smith et al., 2018) and potentially compress sea ice through wave radiation stress (Stopa et al., 2018a). In turn, ice covers suppress wave–wind interaction by reducing the fetch. They also alter wave dispersion and attenuation through scattering and dissipation (Squire, 2007, 2018, 2020).
Modeling surface gravity waves in polar oceans requires the knowledge of many source terms. These source terms include wind inputs and dissipation, nonlinear transfer between frequencies, and wave–ice interaction. WAVEWATCH III^{®} (WW3; WAVEWATCH III^{®} Development Group, 2019), one of the most widely used thirdgeneration wind–wave models for global and regional wave forecasts, has implemented several dispersion and dissipation (IC0, IC1, IC2, IC3, IC4 and IC5) as well as scattering (IS1 and IS2) parameterizations, called “switches” in WW3, to estimate the effect of ice on waves. Of the two scattering switches, IS1 redistributes a constant fraction of the incoming wave energy to all directions isotropically. IS2 adopts a linear Boltzmann equation and an estimation of maximum floe diameter due to ice breakup to model wave scattering combined with a creeprelated dissipation. As is discussed in Sect. 3.2, scattering is negligible in the present dataset; hence it will not be considered in the present study. Of the six different dispersion and dissipation parameterizations, IC0 and IC1 are outofdate. IC4 is empirically determined from data obtained in the Southern Ocean. IC2 is based on a theory that states that wave damping is entirely due to the eddy viscosity in the water body beneath the ice cover. The other two parameterization, IC3 and IC5, both assume that wave damping is due to the processes within the ice cover alone. In this study, we will focus on IC3 and IC5 parameterizations. The same method may be applied to calibrate IC2. Both IC3 and IC5 theorize that sea ice can store and dissipate mechanical energy; hence they model the ice cover as a viscoelastic material. The storage property is reflected in the potential and elastic energy, and the dissipative property is in the equivalent viscous damping. The difference between IC3 and IC5 is that IC3 is an extension of the viscous ice layer model with a finite thickness (Keller, 1998) by including elasticity into a complex viscosity (Wang and Shen, 2010b), while IC5 is an extension of the thin elastic plate model (Fox and Squire, 1994) introduced by Mosig et al. (2015) by adding viscosity into a complex shear modulus (equivalent to the complex viscosity via the Voigt model). Below we refer to these two viscoelastic models as WS and FS, respectively. Each of the two models shows a frequencydependent wave propagation determined by the dispersion relation. This relation, which depends on the viscoelastic parameters, specifies how the wave number $k={k}_{\mathrm{r}}+i{k}_{\mathrm{i}}$ is related to the wave frequency f. The complex wave number k contains a real part k_{r}, which determines the wave speed, and an imaginary part k_{i}, which determines the damping rate. To use these models for wave forecasts in icecovered seas, one needs to determine the viscoelastic parameters for all types of ice covers. To derive these parameters using first principles is challenging, as demonstrated by De Carolis et al. (2005), who obtained the viscosity of grease ice using principles in fluid mechanics of a suspension. Alternatively, an inverse method has been adopted to parameterize IC3. Using in situ data from the R/V Sikuliaq field experiment (Thomson et al., 2018), the WS model was calibrated to match the observed wave number and attenuation (Cheng et al., 2017). This calibration was carried out in a marginal ice zone (MIZ) populated predominantly with grease and pancake ice. Although it showed good agreement between the calibrated model and the field data in the frequency band containing most of the wave energy, these calibrated values are limited to the grease and pancake ice types. Further into the ice cover, where more rigid ice with larger floes is present, how the viscoelastic parameters might change is unknown at present. Note that models can only be as robust as the training data. Therefore, using a broader type of sea ice data will result in a more robust model.
Advancements made in remotesensing technology have provided opportunities for such model calibration. Ardhuin et al. (2017) developed a method to conditionally invert from icecovered water wave orbital motion to directional wave spectra from synthetic aperture radar (SAR) images based on the velocitybunching mechanism. The methodology is further improved in Ardhuin et al. (2018) and Stopa et al. (2018b) to study wave state in the icecovered Beaufort Sea using SAR images, which were captured during the R/V Sikuliaq field experiment by the satellite Sentinel1A. Stopa et al. (2018b) retrieved wavenumberdependent spectra under the partially and fully icecovered regions by substantially reducing data contamination by ice features with a similarlength scale as the wavelength. Using the SARretrieved wave spectra, Stopa et al. (2018b) found that the significant wave height attenuated steeply prior to the first appearance of ice leads (denoted as FAL hereafter), with a milder attenuation rate after the FAL. The definition of FAL is described in Appendix A. Furthermore, Monteban et al. (2019) used two overlapping burst images from the Sentinel1 separated by ∼2 s to investigate wave dispersion in the Barents Sea. Based on the method proposed by Johnsen and Collard (2009), this time separation between subsequent images was sufficient to result in a less noisy and higherquality imaginary spectrum, therefore allowing them to obtain spatiotemporal information of the dispersion relation. Their results showed that for long waves (100–350 m) in thin ice (<40 cm), the dispersion relation was the same as in open water.
This study uses the dataset reported in Stopa et al. (2018b). It includes two parts of data analysis: obtaining spectral wave attenuation rates from the retrieved wave data and then using these attenuation rates for waveinice model calibration. This paper is organized as follows: Sect. 2 introduces the wave spectra data from Stopa et al. (2018b) used in this study. From this data, we obtain the dominant wave numbers and the related wave directions over the studied domain. In Sect. 3, we use the directional wave spectra to derive the apparent attenuation rate and then the attenuation rate due to sea ice alone. In Sect. 4, we calibrate two viscoelastic waveinice models using the obtained wavenumberdependent attenuation data. The methodology presented in these two sections is similar to that in Cheng et al. (2017) with modifications to resolve the difference of the spectral data types of the waves. In Cheng et al. (2017), the spectral data were between energy and frequency, while in Stopa et al. (2018b) they were between energy and wave number. Section 5 discusses the characteristics of wavelength and attenuation in pack ice, calibrated viscoelastic parameters and wave attenuation modeling. The final conclusions are given in Sect. 6.
During the R/V Sikuliaq experiment, the Sentinel1A, equipped with a synthetic aperture radar (SAR), acquired six sequential images around 16:50 UTC on 12 October 2015. These SAR images covered a 400 km by 1100 km region including open water, grease and pancake ice and pack ice (refer to Fig. 1 in Stopa et al., 2018b). A large wave event during this time with wave heights exceeding 4 m in the captured region provided quality wave data. From these SAR images, for most of this icecovered region Stopa et al. (2018b) obtained twodimensional wave spectra data E(k_{x}, k_{y}), where E indicates wave energy density and k_{x} and k_{y} are the wave number components in the range and azimuth directions of the satellite track, respectively. Details of this dataset and its retrieval may be found in Stopa et al. (2018b). We convert the twodimensional spectrum at each location into an equivalent wavenumberdirection spectrum E(k_{r},θ), where ${k}_{\mathrm{r}}=\sqrt{{k}_{x}^{\mathrm{2}}+{k}_{y}^{\mathrm{2}}}$ and $\mathit{\theta}=\mathrm{atan}\left({k}_{y}/{k}_{x}\right)$ indicate wave number and direction, respectively. The wave number is discretized from 0.011 to 0.045 m^{−1} with an increment of 0.002 m^{−1}, and the direction is discretized into 360 bins with a 1^{∘} bin width. We then define wavenumberdependent main wave directions for each given spectrum E(k_{r},θ) in the following way: for each wave number k_{r} we fit the corresponding E curve with a Gaussian function (an example is given in Fig. S1 in the Supplement). The mean of this Gaussian function is defined as the main wave direction ${\mathit{\theta}}_{{k}_{\mathrm{r}}}$ for this k_{r}. Furthermore, we define the dominant wave number k_{r,dominant} to be the one corresponding to the maximum directionally integrated wave energy $\int E\left({k}_{\mathrm{r}},\mathit{\theta}\right)\mathrm{d}\mathit{\theta}$.
An overview of the processed wave conditions in terms of the dominant wave number and its main direction is given in Fig. 1a along with the locations of the ice edge and other in situ observations. The associated ice conditions in the region are presented in Fig. 1b, c. Significant spatial variability is observed in both the wave and ice conditions. Figure 1a presents a subregion captured by the SAR images covering a portion of Alaska to the last azimuth position where waves are detectable in the images. Colors indicate the dominant wave number distribution of the retrievals, and arrows indicate their main directions. Cell size is coarsened to 12.5 km×12.5 km to enhance visualization. The ice edge is indicated by the contours of ice concentration (<0.4) from AMSR2 (Advanced Microwave Scanning Radiometer 2; https://seaice.unibremen.de/data/amsr2, last access: 21 June 2020, Spreen et al., 2008). An in situ buoy, AWACI (a subsurface Nortek acoustic wave and current buoy, moored at 150^{∘} W, 75^{∘} N), is marked by a magenta asterisk. Except for the in situ observations from the Sikuliaq ship (green diamond) and several drifting buoys (blue dots) near the ice edge, ice morphology information including ice types and their partial concentrations is absent. Nevertheless, the FAL (red dots) presumably marks the separation between discrete floes and a semicontinuous ice cover with dispersed leads. The ice condition below (before) the FAL was more complex, with thinner ice and lower concentration than that above (after) the FAL.
We select the study domain defined by the azimuth from 450 to 750 km and the range from 0 to 400 km. This domain contains most of the wave data retrieved in the pack ice field. Figure 1b, c show the distributions of ice concentration from AMSR2 and ice thickness from SMOS (Soil Moisture and Ocean Salinity; https://seaice.unibremen.de/data/smos, last access: 21 June 2020, Huntemann et al., 2014), respectively, in the region of interest. As shown in Cheng et al. (2017; Supplement Figs. S6 and S7), these two ice products compared the best with in situ observations in the MIZ. Their accuracies in the pack ice zone are uncertain. The purpose of this work is to retrieve the spectral wave attenuation and use the result to calibrate viscoelastic models in regions dominated by thin pack ice (thickness: <0.3 m). As shown in Fig. 1a, k_{r,dominant} generally declines crossing the FAL towards the north. Before the FAL, k_{r,dominant} increases in the direction of increasing range and decreasing ice concentration but is insensitive to ice thickness variation. After the FAL, k_{r,dominant} decreases in the wavepropagating direction (arrows) associated with the increase in ice thickness, where the ice field is presumably a semicontinuous cover populated with leads.
Figure 2a, b show twodimensional histograms of the main wave direction for each wave number ${\mathit{\theta}}_{{k}_{\mathrm{r}}}$ before and after the FAL, respectively, where wave direction is defined in the meteorological convention (i.e., the direction “from” in degrees clockwise from true north). The gray scale indicates the occurrence frequency of ${\mathit{\theta}}_{{k}_{\mathrm{r}}}$ in each k_{r} bin. We observe a significant change in ${\mathit{\theta}}_{{k}_{\mathrm{r}}}$ crossing the FAL: before the FAL ${\mathit{\theta}}_{{k}_{\mathrm{r}}}$ spreads from 160 to 190^{∘}, while ${\mathit{\theta}}_{{k}_{\mathrm{r}}}$ is more tightly clustered from 180 to 200^{∘} after the FAL. The difference before and after the FAL is most significant for k_{r}<0.035 m^{−1}.
For the subsequent spectral analysis, we further restrict k_{r} to 0.019 m${}^{\mathrm{1}}\le {k}_{\mathrm{r}}\le min\left(\mathrm{2}\mathit{\pi}/{\mathit{\lambda}}_{\mathrm{c}},\phantom{\rule{0.125em}{0ex}}\mathrm{0.045}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{1}}\right)$, where λ_{c} is the azimuth cutoff indicating the minimum resolvable wavelength from the SAR imagery (Stopa et al., 2015; Ardhuin et al., 2017). Below this wavelength, the patterns of icecovered ocean surface roughness from SAR imagery are more related to ice features rather than waves (Stopa et al., 2018b). For k_{r}<0.019 m^{−1}, energy density $E\left({k}_{\mathrm{r}},{\mathit{\theta}}_{{k}_{\mathrm{r}}}\right)$ is small with high spatial variation (Fig. B1 in Appendix B) and hence treated as a noise band and removed from further study.
The corresponding frequency f (wave period T) range is estimated as 0.067 to 0.11 Hz (9 to 15 s) using the openwater dispersion relation (2πf)^{2}=gk_{ow}, where g is the gravitational acceleration, and k_{ow} indicates the wave number in open water. The actual dispersion relation, which is likely dependent on the ice condition, cannot be measured from our instantaneous SAR data. Using accelerometers deployed on ice floes (Fox and Haskell, 2001) and from marine radar or buoy data (Collins et al., 2018), it was found that the wave number in the icecovered region within about 10 km from the ice edge was close to that in open water. Further into the ice cover, Monteban et al. (2019) used two overlapping burst images from the Sentinel1 separated by ∼2 s to investigate wave dispersion in the Barents Sea. The ice concentration and thickness in their study were similar to the present study. The computed wave dispersion within the sea ice for long waves (peak wavelengths: >100 m) was in good agreement with the openwater dispersion relation.
In this section, we first obtain the apparent spectral wave attenuation by pairing the directional spectra at different locations. We then remove the contributions of wave energy between pairs of locations from wind input, wave breaking dissipation and nonlinear transfer between frequency components to extract the spectral attenuation due to ice effects alone. These spectral attenuation rates due to ice are be used in Sect. 4 to calibrate two viscoelastic waveinice models.
3.1 Apparent wave attenuation
We define an apparent wave attenuation for each k_{r} by assuming exponential decay of the wave spectral densities from location A to location B:
where $\overline{\mathit{\theta}}=\frac{{\mathit{\theta}}_{{k}_{\mathrm{r},A}}+{\mathit{\theta}}_{{k}_{\mathrm{r},B}}}{\mathrm{2}}$ is the average of the main wave directions at A and B; D and θ_{AB} are the distance and direction from A to B in the longitude–latitude coordinates, respectively. A selected pair of A and B is named as a pair hereafter. To reduce the uncertainties of naturally present ice and wave variability, a set of quality control criteria are applied to a pair before further analysis:

Ignore substantially oblique waves. The difference between $\overline{\mathit{\theta}}$ and the vector from location A to location B is restricted to $\left\overline{\mathit{\theta}}{\mathit{\theta}}_{AB}\right\le $15^{∘}.

Avoid strong spatial variations of ice condition between A and B. Distance between A and B is restricted to D≤60 km.

As the wave state changes significantly across the FAL as mentioned in Sect. 2, no pair across the FAL is selected. This criterion enables us to detect, if available, the influence of ice morphology on wave attenuation.

Ensure points A and B are both subject to the same wave system. The Pearson correlation coefficient between energy spectra E_{A}(k_{r},θ) and E_{B}(k_{r},θ) is required to be greater than 0.9. The Pearson correlation coefficient is defined as $r=\frac{{\mathrm{\Sigma}}_{i}\left({x}_{i}\overline{x}\right)\left({y}_{i}\overline{y}\right)}{\sqrt{{\mathrm{\Sigma}}_{i}{\left({x}_{i}\overline{x}\right)}^{\mathrm{2}}}\sqrt{{\mathrm{\Sigma}}_{i}{\left({y}_{i}\overline{y}\right)}^{\mathrm{2}}}}$, where x_{i} and y_{i} are the power spectral density (PSD) values at the ith wave number.

Exclude outliers of the spectral attenuation where the wave energy of B is higher or close to that of A. Thus, $\mathit{\alpha}\left({k}_{\mathrm{r}}\right)>{\mathrm{10}}^{\mathrm{6}}$ m^{−1} is required.

A selected pair has to have at least 10 data points of α(k_{r}) to do calibration in Sect. 4.
We obtain 2634 pairs (2194 pairs before the FAL and 440 pairs after the FAL) through the above criteria to calculate the apparent wave attenuation α(k_{r}) by Eq. (1). The results are sorted statistically to show the occurrence frequency of α(k_{r}). The α domain is equally divided into 30 bins, from 10^{−6} to 10^{−4} m^{−1} in log scale, and the k_{r} domain is equally divided from 0.011 to 0.045 m^{−1} with an increment of 0.002 m^{−1} as mentioned earlier. Figure 3a, b show twodimensional histograms of α against k_{r} before and after the FAL, respectively. Gray scales represent the occurrence of α at each combined bin of α and k_{r}. Red curves indicate the most frequent occurrence of α(k_{r}) against k_{r}. It is observed that more data are obtained before the FAL, with a slightly increasing trend of α(k_{r}) versus k_{r} before the FAL, while α(k_{r}) obtained after the FAL are mostly lower and independent of k_{r}.
The range of this apparent spectral attenuation is in agreement with Stopa et al. (2018b), in which the authors selected multiple tracks throughout the ice region and focused on the overall decay of the significant wave height over hundreds of kilometers. The reduction in attenuation crossing the FAL shown in Fig. 3 is also consistent with Stopa et al. (2018b), who reported a drop in attenuation of the significant wave height after the FAL. In Appendix B, we show the results of this longrange attenuation against wave number (Fig. B1). The difference of attenuation obtained by the two methods, one based on short distances (<60 km) to reduce the effect of ice type variability and the other over long distances (∼300 km), is discussed in Sect. 5.
Because of the large study domain and the apparent difference of k_{r,dominant} in the east–west direction before the FAL, as shown in Fig. 1a, it is worthwhile examining the regional variability. Figure 4 displays the results related to α and k_{r} in three subdomains from west to east, bounded by longitudes (150, 151^{∘} W), (148, 149^{∘} W) and (146, 147^{∘} W). The left column shows 467, 233 and 132 pairs selected in the three longitude bins, respectively. Each segment indicates a pair selected in Sect. 3.1, with ends corresponding to the locations A and B and its color indicates the mean ice thickness between the two. In the middle column, the α values corresponding to k_{r}= 0.019, 0.025, 0.031 and 0.039 m^{−1} obtained from these pairs are separated into two groups: before and after the FAL. In each subdomain delineated by the longitude and the FAL, the distribution of α for each k_{r} is presented by a violin plot whose width indicates the probability density distribution of α, and a circle marker inside indicates the median. These violin plots show that α before the FAL is larger than that after the FAL for all k_{r} in all three subdomains. We note that the sample size after the FAL in the easternmost subdomain (146, 147^{∘} W) is the lowest, corresponding to the largest variability of the violin plots. To track wave spectrum evolution from near the ice edge towards the interior ice zone, we collect wave PSDs ($\int E\left({k}_{\mathrm{r}},\mathit{\theta}\right)\mathrm{d}\mathit{\theta})$ at the north end of each selected pairs in the left column. The averaged PSDs of this collection per 0.1^{∘} in latitude are presented by curves in the right column, where colors indicate the latitude. As defined earlier, k_{r} associated with the peak of a PSD curve is k_{r,dominant}. Before the FAL the magnitude of the PSD drops rapidly, while k_{r,dominant} varies slightly as the latitude increases. In contrast, the PSDs change slowly after the FAL, while k_{r,dominant} decreases quickly. Note that at high latitudes, the PSD at k_{r,dominant} (red) is higher than that of the same wave number at low latitudes. We revisit this phenomenon in the discussion section.
3.2 Wave attenuation due to ice effect
Following Eq. (1), the apparent attenuation obtained above is determined by the energy difference between two locations. This apparent attenuation is the result of multiple source terms, including the wind input S_{in}, damping through wave breaking, and swell dissipation S_{ds}, the energy transfer due to nonlinear interactions among spectral components S_{nl} and the dissipation and scattering of wave energy due to ice cover S_{ice}. In this section, we derive the attenuation rate due to S_{ice} from the measured apparent attenuation α.
The radiative transfer equation for surface waves concerning all of the above effects is
where $E=E\left({k}_{\mathrm{r}},\mathit{\theta},x,t\right)$ is the power spectral density depending on wave number, direction and location x; c_{g} is the group velocity; and U is the current velocity. Note that U in this region is below 0.1 m s^{−1} according to OSCAR (Ocean Surface Current Analyses Realtime, ESR, 2009 and Bonjean and Lagerloef, 2002). Because the current speed is at least 1 order of magnitude below the estimated c_{g} using the openwater dispersion relation, we may drop it from Eq. (2). Furthermore, consistent with the fact that the dispersion relation in this study is close to that of the open water, c_{g} is relatively constant. Adopting the exponential wave decay along x, i.e., $E\left({k}_{\mathrm{r}},\mathit{\theta},x\right)=E\left({k}_{\mathrm{r}},\mathit{\theta},x=\mathrm{0}\right){e}^{\mathrm{2}\mathit{\alpha}x}$, we have $\frac{\partial \left({c}_{\mathrm{g}}+U\right)E}{\partial x}\approx {c}_{\mathrm{g}}\frac{\partial E}{\partial x}=\mathrm{2}{c}_{\mathrm{g}}\mathit{\alpha}E$.
Next, we examine the temporal derivative of wave energy $\frac{\partial E}{\partial t}$. It is challenging to calculate $\frac{\partial E}{\partial t}$ from the nearly instantaneous SAR imagery. Instead, we estimate $\frac{\partial E}{\partial t}$ using hourly wave spectra data from two sources around the time stamp of the SAR imagery: the in situ AWACI marked in Fig. 1 and the WW3 simulations of the whole domain (reference model run of Ardhuin et al., 2018). From the AWACI data, we obtain $\frac{\partial E}{\partial t}$ and compare that with ${c}_{\mathrm{g}}\frac{\partial E}{\partial x}$ using SARretrieved wave data at the AWACI site. From the WW3 simulations, we obtain both terms over the whole study domain. Both results consistently show that $\frac{\partial E}{\partial t}$ is at least 2 orders of magnitude below ${c}_{\mathrm{g}}\frac{\partial E}{\partial x}$. Thus, $\frac{\partial E}{\partial t}$ is dropped from Eq. (2).
The other source terms S_{in} and S_{ds} are estimated using formulations from Snyder et al. (1981) and Komen et al. (1984). For S_{nl}, we select the discrete interaction approximation (DIA; Hasselmann and Hasselmann, 1985; Hasselmann et al., 1985) to estimate its value. Note that those parameterizations were formulated from an openwater study. How they might change in the presence of ice covers is an open question. The formulations and associated coefficients used here are described in Cheng et al. (2017). Wherever needed in these formulations, f is approximated by the openwater dispersion relation with the measured k_{r}. Ice concentration is from AMSR2, ice thickness is from SMOS and wind data are from the Climate Forecast System Reanalysis (CFSR; Saha et al., 2010).
Iceinduced wave attenuation is known from the dissipation of wave energy and scattering of waves (e.g., Wadhams et al., 1988; Squire et al., 1995; Montiel et al., 2018). Here we attribute the attenuation entirely to the dissipative process with the following arguments: Ardhuin et al. (2018) reported that in the studied region, ice floe scattering has a weaker effect on wave attenuation compared with other processes, including the boundary layer beneath the ice, inelastic flexing of ice cover and waveinduced ice fracture. The inelastic flexing and ice fracture may be considered as part of the dissipative mechanism within the ice cover already included in the viscous coefficient, but scattering is a redistribution of energy, which must be isolated from the apparent attenuation before using the data to calibrate the viscoelastic models. We estimate the scattering effect based on the study of Bennetts and Squire (2012) as follows: in that study, wave attenuation by floes, cracks and pressure ridges was examined. In the absence of in situ observations, we assume that few and small ridges are present in the studied ice cover (thickness: <0.3 m); hence the effect of ridges is negligible. In our case of long waves propagating through such thin ice cover, Bennetts and Squire (2012) have shown that the floes produce much more attenuation than the cracks. Without in situ observation, WW3 simulations (REF run in Ardhuin et al., 2018) implementing waveinduced fracturing of ice floes gave a range of estimated maximum floe diameter from 70 to 150 m. Using this range of floe diameter, the theoretical scattering results from Bennetts and Squire (2012) indicate that the floe scatteringinduced attenuation is about 10^{−7} m^{−1}, which is negligible compared to the α shown in Fig. 3.
With all the above simplifications and the assumed exponential decay of wave energy, Eq. (2) becomes
which yields
where k_{i} is the attenuation rate due to the ice cover. The occurrence of the k_{i} data is presented by twodimensional histograms of k_{i} against k_{r} in Fig. 5. Not surprisingly, since the effect of wind is low due to the low openwater fraction in the study region and the relatively short distances between the pairs for the nonlinear transfer of energy to accumulate, we observe that k_{i} is very close to the apparent attenuation α. We discuss this iceinduced dissipation further in Sect. 5.
In this section, we present the calibration of the WS model and the FS model. We begin with a brief introduction of the two models then describe the calibration procedure. In both models, k_{r} and k_{i} are solved from the respective dispersion relations for given ice thickness and viscoelastic parameters. The viscoelastic parameters are optimized by minimizing the overall difference of the k_{i}–k_{r} relationship between the models and data obtained in Sect. 3.
4.1 Dispersion relation
The dispersion relations of the two models are written as
For the WS model,
For the FS model,
where H is water depth, σ=2πf is the angular frequency, ρ_{ice} and ρ_{water} are the densities of ice and water, respectively, $k={k}_{\mathrm{r}}+i{k}_{\mathrm{i}}$ is a complex wave number, h is the ice thickness, ${a}^{\mathrm{2}}={k}^{\mathrm{2}}\frac{i\mathit{\sigma}}{{\mathit{\nu}}_{e}}$, S_{k}=sinhkh, S_{a}=sinhah, C_{k}=coshkh, C_{a}=coshah, $N=\mathit{\sigma}+\mathrm{2}i{k}^{\mathrm{2}}{\mathit{\nu}}_{e}$, ${G}_{\mathit{\nu}}=Gi\mathit{\sigma}\mathit{\nu}{\mathit{\rho}}_{\mathrm{ice}}$ and ${\mathit{\nu}}_{e}=\mathit{\nu}+\frac{iG}{{\mathit{\rho}}_{\mathrm{ice}}\mathit{\sigma}}$ and V is Poisson's ratio. Equivalent shear modulus G and kinematic viscosity ν in both models are to be calibrated. In this study, we use H=1000 m for deep water, ρ_{water}=1025 kg m^{−3}, ρ_{ice}=922.5 kg m^{−3} and V=0.3 for ice.
Hereafter, we use superscripts t for the theoretical values and m for the measured data. Specifically, the theoretical wave number and attenuation rate are denoted as ${k}_{\mathrm{r}}^{\mathrm{t}}$ and ${k}_{\mathrm{i}}^{\mathrm{t}}$. They are calculated for each set of G,ν values. Attenuation data due to the ice effect obtained in Sect. 3 are denoted as ${k}_{\mathrm{i}}^{\mathrm{m}}$. The corresponding wave number ${k}_{\mathrm{r}}^{\mathrm{m}}$ is the discretized k_{r} values from 0.011 to 0.045 m^{−1} with an increment of 0.002 m^{−1}.
4.2 Calibration methodology
In this section, we optimize G,ν by fitting the k_{i}–k_{r} relationship from the model via Eq. (5) to the measured values from Eq. (4). Specifically, for given G and ν, we solve arrays of ${k}_{\mathrm{r}}^{\mathrm{t}}$ and ${k}_{\mathrm{i}}^{\mathrm{t}}$ through Eq. (5) for each f densely sampled from 0.0001 to 1 Hz. We then interpolate the results to obtain the theoretical k_{i} at each wave number ${k}_{\mathrm{r}}^{\mathrm{m}}$ from the data in Sect. 3.2. For each of these ${k}_{\mathrm{r}}^{\mathrm{m}}$ we have the measured dissipation rate ${k}_{\mathrm{i}}^{\mathrm{m}}$ shown in Sect. 3.2. We now use an optimization procedure to determine the bestfit parameters G and ν. The objective function for the optimization is defined as the weighted sum of the differences between ${k}_{\mathrm{i}}^{\mathrm{t}}$ and ${k}_{\mathrm{i}}^{\mathrm{m}}$ over ${k}_{\mathrm{r}}^{\mathrm{m}}$, i.e.,
where ∥•∥_{2} is the L2 norm operator and $w\left({k}_{\mathrm{r}}^{\mathrm{m}}\right)$ is a weighting factor to account for the distributions of wave energy and attenuation rate. Cheng et al. (2017) tested two weighting factors, $w=\int E\left(\mathit{\theta},f\right)\mathrm{d}\mathit{\theta}$ and $\int E\left(\mathit{\theta},f\right){f}^{\mathrm{4}}\mathrm{d}\mathit{\theta}$, in calibrating the WS model. In that study, the measured data had a range of f varied from 0.05 to 0.5 Hz, and the attenuation rate varied from 10^{−6} to 10^{−3} m^{2} s^{−1}. The authors found that using $w=\int E\left(\mathit{\theta},f\right)\mathrm{d}\mathit{\theta}$ fitted attenuation best at the most energetic wave band, while $w=\int E\left(\mathit{\theta},f\right){f}^{\mathrm{4}}\mathrm{d}\mathit{\theta}$ performed better to capture significant increasing in k_{i} at high frequencies. No weighting factor could produce a fitting over the entire spectral attenuation curve. For the present study, the variation of spectral attenuation is small, as shown in Fig. 5, with no particular region of emphasis. We thus choose $w=\sqrt{\int E\left(\mathit{\theta},{k}_{\mathrm{r}}\right)\mathrm{d}\mathit{\theta}}$, which has a broad band around the peak energy of the wave field.
We choose a search domain (10^{−7} Pa $\le G\le {\mathrm{10}}^{\mathrm{10}}$ Pa and 10^{−4} m^{2} s${}^{\mathrm{1}}\le \mathit{\nu}\le {\mathrm{10}}^{\mathrm{4}}$ m^{2} s^{−1}) for the WS model. This search domain covers all other reported viscoelastic values for ice covers (e.g., Newyear and Martin, 1999; Doble et al., 2015; Zhao and Shen, 2015; Rabault et al., 2017). For the FS model, it is known that very large G,ν are needed to obtain the level of k_{i} observed (Mosig et al., 2015). We thus choose a large search domain 10 Pa $\le G\le {\mathrm{10}}^{\mathrm{20}}$ Pa and $\mathrm{10}\phantom{\rule{0.125em}{0ex}}\phantom{\rule{0.125em}{0ex}}{\mathrm{m}}^{\mathrm{2}}\phantom{\rule{0.125em}{0ex}}{\mathrm{s}}^{\mathrm{1}}\le \mathit{\nu}\le {\mathrm{10}}^{\mathrm{15}}$ m^{2} s^{−1}. This parameter range is far beyond the measured data from solid ice (Weeks and Assur, 1967). The global optimization procedure is performed by using the genetic algorithm with function ga in MATLAB and the Global Optimization Toolbox R2016a (2016).
4.3 Results
For the WS model, the optimized G,ν are gathered into a small cluster around $G\cong {\mathrm{10}}^{\mathrm{4}}$ Pa and a large cluster around G≅10^{5} Pa, with large variation in ν. The existence of two separate clusters of calibrated results was also found in grease–pancake mixtures near the ice edge (Cheng et al., 2017). Note that multiple solutions (optimized G,ν pairs) could be obtained by optimizing a nonlinear system such as that shown in Eq. (5). Thus, constraints are applied to select solutions that are physically plausible. Note that the solutions around $G={\mathrm{10}}^{\mathrm{4}}$ Pa and solutions with ν near 10^{4} m^{2} s^{−1} lead to the residual from Eq. (6) insensitive to G. It implies that the modeled material is viscousdominant with almost nil elastic effect. Those cases are incompatible with the pack ice condition and are thus removed in the following analysis. Scatter plots of the remaining data points (G,ν) are presented by Fig. S2 in the Supplement. For these remaining data, we apply the bivariate Gaussian distribution to obtain the 90 % probability range of (G,ν) before and after the FAL, respectively. Means and covariances of the related probability density functions are given in Table 1. A slight difference in the distributions of G,ν before and after the FAL is noticed. The mean elasticity is slightly lower before the FAL than after the FAL, and the mean viscosity is slightly higher before the FAL than after the FAL. The results imply that the ice layer behaves more elastic in the inner ice field and more viscous towards the ice edge. It is worth noting that the ranges of both G and ν are about 1 order of magnitude greater than those of the grease and pancake ice (Cheng et al., 2017). As expected, the effective elasticity in the WS model for more solid ice is closer to the elastic modulus of sea ice.
For the FS model, the calibrated (G,ν) are clustered, whereas scatter plots of (G,ν) are given in Fig. S3 in the Supplement. Hence all data points are used in the bivariate Gaussian fitting. The resulting mean and covariance values are also given in Table 1. Notice that the mean values of G, ν from the FS model are extremely larger than the intrinsic values of ice. The distribution of calibrated values is further discussed in the discussion section. Hereafter we only elaborate on the results of the WS model. Figure 6a shows the overall comparison of k_{i}–k_{r} between the measurements (gray) and the WS model (black). Both the median (solid curve) and 90 % boundaries (dashed curves) are in good agreement. In Fig. 6a, we superimpose an empirical model from Meylan et al. (2014), which is also included in WW3 to account for the ice effect as switch IC4. By fitting the wave buoy data from the Antarctic MIZ obtained in 2012, Meylan et al. (2014) proposed a simple perioddependent attenuation rate ${k}_{\mathrm{i}}\left(T\right)=\frac{\mathrm{2.12}\times {\mathrm{10}}^{\mathrm{3}}}{{T}^{\mathrm{2}}}+\frac{\mathrm{4.59}\times {\mathrm{10}}^{\mathrm{2}}}{{T}^{\mathrm{4}}}$, which is converted into a k_{i}–k_{r} relation through the openwater dispersion relation. The empirical formula predicted that attenuation was consistent with the range we obtained here but with a higher sensitivity to k_{r}. Figure 6b shows the normalized wave number $\frac{{k}_{\mathrm{r}}}{{k}_{\mathrm{ow}}}$ against f using the optimized G,ν in the WS model. The deviation of $\frac{{k}_{\mathrm{r}}}{{k}_{\mathrm{ow}}}$ from 1 is less than 5 %, indicating the modeled wave dispersion agrees with the openwater dispersion relation. It is consistent with Monteban et al. (2019), who investigated wave dispersion in ice from the SAR data in the Barents Sea as the studied range of wavelength and ice thickness are similar. The counterpart of the FS model is given in Fig. S4 in the Supplement.
Figure 7a, b show distributions of calibrated G and ν from the WS model in the azimuthrange domain, respectively. Contours represent ice thickness distribution from SMOS. The FAL is marked as red dots. The spatial domain is divided into 12.5 km×12.5 km cells to enhance visualization. Cell color indicates the averaged value of log_{10}G and log_{10}ν of all pairs with midpoints inside the cell. Ignoring the outliers, Fig. 7a shows that G is generally larger after the FAL (thicker ice with leads) than before the FAL, while ν has an opposite trend. Note that the outliers could be from multiple sources, such as noise in the retrieved wave data, spatial variability of ice condition, assumptions made in selecting pairs and calculating attenuation and the adopted genetic algorithm in the model calibration. The counterpart of the FS model is given in Fig. S5 in the Supplement.
In this section, we discuss some key wave characteristics mentioned in the analysis and the behaviors of calibrated model parameters. Some thoughts about waveinice modeling are provided at the end.
5.1 Evolution of wave characteristics
The behaviors of k_{r,dominant} and PSD in Figs. 1 and 4 may come from multiple mechanisms. The decrease in k_{r,dominant} towards the interior ice is in agreement with a similar study in the MIZ (Shen et al., 2018) as well as the lengthening of dominant wave periods reported in many field observations (e.g., Robin, 1963; Wadhams et al., 1988; Marko, 2003; Kohout et al., 2014; Collins et al., 2015). The phenomenon is commonly explained by the lowpass filter mechanism in the literature. That is, higherfrequency waves are preferentially attenuated; thus the power spectrum peak shifts to longer waves. However, though this mechanism is supported by the increasing trend of k_{i} with increasing wave number before the FAL, it is insufficient to explain the observations after the FAL, where the k_{i}–k_{r} curve is practically flat (Fig. 5), and the peak of the PSD shifts toward lower k_{r} with increasing latitude (Fig. 4). A downshift of wave energy towards lower wave number might also be related to nonlinear wave–wave interaction. Through this interaction, spectral wave energy in the main wave frequency could transfer into the lowfrequency components. How to quantify these mechanisms is still an open question, with few data – especially in pack ice – being available and the skills in observations still developing. A piece of critical information is the relationship between the wave period and wave energy, especially in the region after the FAL. Monteban et al. (2019) showed promising applications of overlapping SAR images separated by a sufficient time gap to obtain such information.
The attenuation rate obtained in this study ($\sim {\mathrm{10}}^{\mathrm{5}}$ m^{2} s^{−1}) against wave number shows a slightly increasing trend before the FAL and a nearly flat trend and lower after the FAL. A similar magnitude of attenuation is found in the MIZ in both the Arctic and Antarctic, as reported in previous studies in the same period range, but larger for shorter periods (e.g., Meylan et al., 2014; Doble et al., 2015; Rogers et al., 2016; Cheng et al., 2017). The analysis of attenuation against wave number in Sect. 3 is obtained for pairs of observations over relatively short distances (<60 km). In Appendix B, we present another method to determine the overall attenuation similar to the analysis of the attenuation of the significant wave height shown in Stopa et al. (2018b). By analyzing wave number by wave number over the entire space that the SAR data covers the apparent attenuation coefficient is obtained by fitting hundreds of data points over long distances, ignoring the higher variation of ice thickness (0.01–0.3 m) and the shift of dominant wave number from the PSD curves. The results are shown in Appendix B, Fig. B2. In some cases, we can even have negative values of this attenuation, meaning energy increases as the wave propagates. What we would like to emphasize by showing this result is that over a very long distance, the ice condition can change significantly in addition to nonlinear energy transfer and wind input and dissipation; thus longdistance spectral analysis may become difficult to interpret.
5.2 About model calibration
The covariance values of the calibrated viscoelastic parameters are greater than those found in grease and pancake ice using buoy data (Cheng et al., 2017). In processing the SAR data, there are many challenges to be dealt with. It is extremely difficult to separate (1) sea ice variability, (2) wave height variability and (3) instrument variability (speckle and incoherent SAR noise). All of these influence the variations we see in the wave spectral data. Still, the scatter of the spectral attenuation shown in Fig. 3 is significant even after the filtering described in Sect. 3.1. This data scatter results in a large range of calibrated model parameters. We believe that the scatter of calibrated G,ν in both WS and FS models could be narrowed down by reducing the uncertainty of measured attenuation data. We also believe that for the FS model the calibrated values and spread could be reduced by modifying the objective function (Eq. 6) with an additional constraint on the wave number, which presently shows an extremely large range as shown in Fig. S4. Regardless, the FS model needs much larger G,ν to produce attenuation rates comparable to the observed data. This is the case not only for the pack ice but also for the MIZ. When analyzing data from the Antarctic MIZ reported in Kohout and Williams (2013) and by further constraining ${k}_{\mathrm{r}}={k}_{\mathrm{ow}}/\mathrm{1.7}$, Mosig et al. (2015) obtained $G=\mathrm{4.9}\times {\mathrm{10}}^{\mathrm{12}}$ Pa, $\mathit{\nu}=\mathrm{5}\times {\mathrm{10}}^{\mathrm{7}}$ m^{2} s^{−1}. Though far below the present case, these values are still orders of magnitude above the intrinsic values measured from solid ice (Weeks and Assur, 1967). Meylan et al. (2018) discussed the attenuation behavior among different dissipative models by assuming small k_{r}−k_{ow}. Specifically, they showed that the FS model produced ${k}_{\mathrm{i}}\approx \frac{{\mathit{\rho}}_{\mathrm{ice}}\left(\mathrm{1}+V\right){h}^{\mathrm{3}}}{\mathrm{6}{\mathit{\rho}}_{\mathrm{water}}{g}^{\mathrm{6}}}\mathit{\nu}{\mathit{\sigma}}^{\mathrm{11}}$ and the pure viscous case of the WS model (i.e., the model in Keller, 1998) produced ${k}_{\mathrm{i}}\approx \frac{\mathrm{4}{\mathit{\rho}}_{\mathrm{ice}}h}{{\mathit{\rho}}_{\mathrm{water}}{g}^{\mathrm{4}}}\mathit{\nu}{\mathit{\sigma}}^{\mathrm{7}}$. The higher the power in σ, the more likely ν is to match the measured k_{i} in the highperiod (lowfrequency) wave range. The FS model also naturally leads to large G. As shown in Mosig et al. (2015), inverting the dispersion relation shown in Eqs. (5a) and (5c) gives $Gi\mathit{\sigma}{\mathit{\rho}}_{i}\mathit{\nu}=\mathrm{6}\frac{{\mathit{\rho}}_{w}{\mathit{\sigma}}^{\mathrm{2}}gk{\mathit{\rho}}_{w}hk{\mathit{\rho}}_{i}{\mathit{\sigma}}^{\mathrm{2}}}{{h}^{\mathrm{3}}{k}^{\mathrm{5}}\left(\mathrm{1}+V\right)}$, the leading term of which yields $G\approx O\left({k}_{\mathrm{r}}^{\mathrm{5}}\right)$ in small k_{r}.
5.3 Thoughts on modeling wave–ice interaction
Damping models play a crucial role in spectral attenuation. At present, using any specific model to describe wave attenuation is tentative. There are many identified damping mechanisms. For instance, boundary layer under the ice cover (Liu and MolloChristensen, 1988; Smith and Thomson, 2019), spilling of water over the ice cover (Bennetts and Williams, 2015), interactions between floes (Herman et al., 2019a, b) and jet formation between colliding floes (Rabault et al., 2019) have all been reported in laboratory or field studies. A full waveinice model that considers all important mechanisms is not yet available. While theoretical improvements are needed to better model wave propagation through various types of ice covers, practical applications cannot wait. Calibrated models that are capable of reproducing key observations must be developed parallel to model improvements. The present study provides a viable way to calibrate two such models available in WW3. These models lump all dissipative mechanisms in the ice cover into a viscous term. This type of model calibration study has two obvious utilities. First, with proper calibration, models can capture the attenuation of the most energetic part of the wave spectrum. Second, the discrepancies may be used to motivate the development of models in which there are missing mechanisms, thereby helping future model development. Because different physical processes may play different roles under various ice morphology, collecting wave data to calibrate these models for various ice types is necessary. Finally, more observations with higherquality data will improve the modeling of the wave–ice interaction and the robustness of the models.
In this study, we use the wave spectra retrieved from SAR imagery to examine the spectral attenuation in young pack ice. These images were obtained in the Beaufort Sea on 12 October 2015. According to the analysis of data retrieved over several hundred kilometers, the observed decrease in wave energy and lengthening of dominant waves towards the interior ice are consistent with earlier in situ observations. We investigate wave attenuation of dominant spectral densities per wave number between two arbitrary locations in the region separated by less than 60 km. Similar attenuation rates are observed for all wave numbers from 0.019 to 0.045 m^{−1} (estimated wave period from 9 to 15 s). After isolating the iceinduced attenuation from these data, we calibrate two viscoelastictype waveinice models through an optimization procedure between the measured data and theoretical results. Both models can generally match the observed attenuation corresponding to the energetic portion of the wave spectra, with a large difference in dispersion (see wave number plots in Figs. 6 and S4). For the WS model, the calibrated shear modulus and the viscous parameter in the region beyond the first appearance of leads with thicker ice is slightly larger and smaller, respectively, than in the region closer to the ice edge before the first appearance of leads. The resulting wave number is within 5 % of that from the openwater dispersion. The wave–ice interaction is complicated due to many coexisting physical processes. At present no models have fully integrated all identified processes. The present study demonstrates a method based on measured data to calibrate existing models so that they can be applied to meet the operational needs. As the models improve, further calibration exercises may be performed accordingly. For example, the eddy viscosity model (Liu and MolloChristensen, 1988) attributed the wave attenuation entirely to the water body under the ice cover. This model is also included in WW3 as the switch IC2, with the eddy viscosity as a tuning parameter. The present analysis can also be used to calibrate that model or models that try to combine both dissipation from the ice cover and the eddy viscosity from the water beneath (e.g., Zhao and Shen, 2018). From a physical point of view, dissipative mechanisms may be present simultaneously inside the ice cover and the water body underneath. However, the calibration of complex models with multiple coexisting processes is a difficult task that requires a much more dedicated study.
We conclude by noting that highresolution spatial data from remote sensing provide new opportunities to investigate the wave–ice interaction over a large distance and with different ice types. However, details of ice types and temporal observations are in development. To reach a full understanding and thus a complete waveinice model, collaboration from observation and modeling efforts is required.
The definition of the first appearance of leads was introduced in Stopa et al. (2018b). Here we provide more details of the methodology used. The SAR sea surface roughness images in Fig. 1 of Stopa et al. (2018b) are divided into 5.1 km×7.2 km subimages with a 50 % overlap of adjacent subimages in the range–azimuth domain. Each subimage contains 512 pixels×512 pixels. The FAL location for each range position is defined as the minimum azimuth position where largescale ice features were detected. Detection of largescale ice features is applied to each SAR subimage as the following. We first compute a onedimensional spectrum of the SAR subimage to produce an image modulation spectrum. The spectrum is then normalized by the maximum energy contained in wavelengths from 100 to 300 m (the wavelength range of the dominant sea state for this event). When the ratio of the average of the normalized image spectra with wavelengths in the range of 600–1000 m and the dominant ocean wave wavelength range from 160–220 m exceeds 0.8, we deem there to be a “largescale” feature such as lead within the image. Figure A1 shows two representative examples of detecting ice leads from SAR images captured before and after the FAL. From the criterion above, there are no leads in the top case, but leads are found in the bottom case. Also notice the change in the probability distribution of the roughness: the mean value changes (lower in the nonlead case compared lead case) and the standard deviation (lower in the nonlead case compared to the lead case).
We investigate the decay of the dominant energy component $E\left({k}_{\mathrm{r}},{\mathit{\theta}}_{{k}_{\mathrm{r}}}\right)$ over a long distance along selected tracks with fixed longitude (145, 146, …, 150^{∘} W). Figure B1 shows $E\left({k}_{\mathrm{r}},{\mathit{\theta}}_{{k}_{\mathrm{r}}}\right)$ (blue dots) collected within a given longitude interval 150±0.1^{∘} W starting from 74.5^{∘} N towards the north. E(k_{r}) for k_{r}<0.019 m^{−1} is too scattered to show any attenuation trend. As k_{r} increases, the data become tighter with a clear decay trend. To obtain attenuation coefficient $\widehat{\mathit{\alpha}}\left({k}_{\mathrm{r}}\right)$, we fit the $E\left({k}_{\mathrm{r}},{\mathit{\theta}}_{{k}_{\mathrm{r}}}\right)$ by an exponential curve (black line) in each panel based on an exponential wave decay assumption.
The resulting $\widehat{\mathit{\alpha}}\left({k}_{\mathrm{r}}\right)$ against k_{r} is shown in Fig. B2 as well as another five curves associated with different longitudes, processed in the same way. Figure B2 shows an increasing trend of $\widehat{\mathit{\alpha}}\left({k}_{\mathrm{r}}\right)$ against k_{r} as expected, and $\widehat{\mathit{\alpha}}({k}_{\mathrm{r}}<\mathrm{0.019}$ m^{−1}) are mostly negative. This method masks the variation of ice morphology and peak shifting of the power spectral density; thus it is insufficient to understand the damping mechanism in the water–ice interaction. Hence, as our interest is to determine the attenuation in more consistent ice conditions, a shorter distance between measuring points is adopted in Sect. 3.
The code and data necessary to reproduce the results presented in this paper can be obtained from the corresponding author. Buoys data and the cruise relevant data of SeaState projection is referred to https://drive.google.com/open?id=0B9Au2ZqQBM5YTJPWXBsV2Q0THc (last access: 20 June 2020, Thomson et al., 2018). The Sentinel1A SAR data are processed by ESA from https://scihub.copernicus.eu/dhus/#/home (last access: 12 January 2018). AMSR2 sea ice concentrations are from https://seaice.unibremen.de/data/amsr2 (last access: 21 June 2020, Spreen et al., 2008). SMOS sea ice thickness data are from https://seaice.unibremen.de/data/smos (last access: 21 June 2020, Huntemann et al., 2014). CFSR wind data are from https://rda.ucar.edu/datasets/ds094.0/ (last access: 20 September 2019, Saha et al., 2011). OSCAR currents data are https://doi.org/10.5067/OSCAR03D01 (ESR, 2019; Bonjean and Lagerloef, 2002).
The supplement related to this article is available online at: https://doi.org/10.5194/tc1420532020supplement.
JS and FA retrieved the wave spectra dataset from SAR images. SC and HHS planned the research and prepared the manuscript with contributions from all authors. SC performed the analysis.
The authors declare that they have no conflict of interest.
The present work is supported by the EUFP7 project SWARP under grant agreement 607476 and the Office of Naval Research (grant numbers N000141310294, N000141712862 and N0001416WX01117). Sukun Cheng was partly funded by the US Office of Naval Research project DASIM II (grant no. N000141812493).
This research has been supported by the EUFP7 project SWARP (grant no. 607476) and the Office of Naval Research (grant nos. N000141310294, N000141712862 and N0001416WX01117). Sukun Cheng was partly funded by the US Office of Naval Research project DASIM II (grant no. N000141812493).
This paper was edited by Michel Tsamados and reviewed by Harry Heorton and one anonymous referee.
Ardhuin, F., Stopa, J., Chapron, B., Collard, F., Smith, M., Thomson, J., Doble, M., Blomquist, B., Persson, O., Collins III, C. O., and Wadhams P.: Measuring ocean waves in sea ice using SAR imagery: A quasideterministic approach evaluated with Sentinel1 and in situ data, Remote Sens. Environ., 189, 211–222, 2017.
Ardhuin, F., Boutin, G., Stopa, J., GirardArdhuin, F., Melsheimer, C., Thomson, J., Kohout, A., Doble, M., and Wadhams, P.: Wave attenuation through an Arctic marginal ice zone on 12 October 2015: 2. Numerical modeling of waves and associated ice breakup, J. Geophys. Res.Oceans, 123, 5652–5668, 2018.
Bennetts, L. G. and Squire, V. A.: Model sensitivity analysis of scatteringinduced attenuation of icecoupled waves, Ocean Model., 45, 1–13, 2012.
Bennetts, L. G. and Williams, T. D.: Water wave transmission by an array of floating discs, P. Roy. Soc. A, 471, 20140698, https://doi.org/10.1098/rspa.2014.0698, 2015.
Bonjean, F. and Lagerloef, G. S. E.: Diagnostic model and analysis of the surface currents in the tropical Pacific Ocean, J. Phys. Oceanogr., 32, 2938–2954, https://doi.org/10.1175/15200485(2002)032<2938:DMAAOT>2.0.CO;2, 2002.
Cheng, S., Rogers, W. E., Thomson, J., Smith, M., Doble, M. J., Wadhams, P., Kohout, A. L., Lund, B., Persson, O. P., Collins III, C. O., Ackley, C. F., Montiel F., and Shen H. H.: Calibrating a viscoelastic sea ice model for wave propagation in the Arctic fall marginal ice zone, J. Geophys. Res.Oceans, 122, 8770–8793, 2017.
Collins III, C. O., Rogers, W. E., Marchenko, A., and Babanin, A. V.: In situ measurements of an energetic wave event in the Arctic marginal ice zone, Geophys. Res. Lett., 42, 1863–1870, 2015.
Collins III, C. O., Doble, M., Lund, B., and Smith, M.: Observations of surface wave dispersion in the marginal ice zone, J. Geophys. Res.Oceans, 123, 3336–3354, 2018.
Comiso, J. C., Parkinson, C. L., Gersten, R., and Stock, L.: Accelerated decline in the Arctic sea ice cover, Geophys. Res. Lett., 35, L01703, https://doi.org/10.1029/2007GL031972, 2008.
De Carolis, G., Olla, P., and Pignagnoli, L.: Effective viscosity of grease ice in linearized gravity waves, J. Fluid Mech., 535, 369–381, 2005.
Doble, M. J., De Carolis, G., Meylan, M. H., Bidlot, J. R., and Wadhams, P.: Relating wave attenuation to pancake ice thickness, using field measurements and model results, Geophys. Res. Lett., 42, 4473–4481, 2015.
ESR: OSCAR third degree resolution ocean surface currents, Ver. 1. PO.DAAC, CA, USA, https://doi.org/10.5067/OSCAR03D01, 2009.
Fox, C. and Haskell, T. G.: Ocean wave speed in the Antarctic marginal ice zone, Ann. Glaciol., 33, 350–354, 2001.
Hasselmann, S. and Hasselmann, K.: Computations and parameterizations of the nonlinear energy transfer in a gravitywave spectrum. Part I: A new method for efficient computations of the exact nonlinear transfer integral, J. Phys. Oceanogr., 15, 1369–1377, 1985.
Hasselmann, S., Hasselmann, K., Allender, J. H., and Barnett, T. P.: Computations and parameterizations of the nonlinear energy transfer in a gravitywave spectrum. Part II: Parameterizations of the nonlinear energy transfer for application in wave models, J. Phys. Oceanogr., 15, 1378–1391, 1985.
Herman, A., Cheng, S., and Shen, H. H.: Wave energy attenuation in fields of colliding ice floes – Part 1: Discreteelement modelling of dissipation due to ice–water drag, The Cryosphere, 13, 2887–2900, https://doi.org/10.5194/tc1328872019, 2019a.
Herman, A., Cheng, S., and Shen, H. H.: Wave energy attenuation in fields of colliding ice floes – Part 2: A laboratory case study, The Cryosphere, 13, 2901–2914, https://doi.org/10.5194/tc1329012019, 2019b.
Huntemann, M., Heygster, G., Kaleschke, L., Krumpen, T., Mäkynen, M., and Drusch, M.: Empirical sea ice thickness retrieval during the freezeup period from SMOS high incident angle observations, The Cryosphere, 8, 439–451, https://doi.org/10.5194/tc84392014, 2014.
Johnsen, H. and Collard, F.: Sentinel1 ocean swell wave spectra (OSW) algorithm definition, Sentinel1 IPF Development (Project No.: 355) Report., 2009.
Keller, J. B.: Gravity waves on icecovered water, J. Geophys. Res.Oceans, 103, 7663–7669, 1998.
Kohout, A. L. and Williams, M. J. M.: Wavesinice observations made during the SIPEX II voyage of the Aurora Australis, 2012, Australian Antarctic Data Centre, 2013.
Kohout, A. L., Williams, M. J. M., Dean, S. M., and Meylan, M. H.: Storminduced seaice breakup and the implications for ice extent, Nature, 509, 604–607, https://doi.org/10.1038/nature13262, 2014.
Kohout, A. L., Williams, M. J. M., Toyota, T., Lieser, J., and Hutchings, J.: In situ observations of waveinduced sea ice breakup, DeepSea Res. Pt. II, 131, 22–27, 2016.
Komen, G., Hasselmann, K., and Hasselmann, K.: On the existence of a fully developed windsea spectrum, J. Phys. Oceanogr., 14, 1271–1285, 1984.
Lange, M., Ackley, S., Wadhams, P., Dieckmann, G., and Eicken, H.: Development of sea ice in the Weddell Sea, Ann. Glaciol., 12, 92–96, 1989.
Liu, A. K. and MolloChristensen, E.: Wave propagation in a solid ice pack, J. Phys. Oceanogr., 18, 1702–1712, 1988.
Marko, J. R.: Observations and analyses of an intense wavesinice event in the Sea of Okhotsk, J. Geophys. Res.Oceans, 108, 3296, https://doi.org/10.1029/2001JC001214, 2003.
MATLAB and Global Optimization Toolbox R2016a: The MathWorks Inc., Natick, Massachusetts, United States, 2016.
Meier, W. N.: Losing Arctic sea ice: Observations of the recent decline and the longterm context, in: Sea Ice, Third Edition, edited by: Thomas, D. N., Hoboken, Wiley Blackwell, 290–303, 2017.
Meylan, M. H., Bennetts, L. G., and Kohout, A. L.: In situ measurements and analysis of ocean waves in the Antarctic marginal ice zone, Geophys. Res. Lett., 41, 5046–5051, 2014.
Meylan, M. H., Bennetts, L. G., Mosig, J. E., Rogers, W. E., Doble, M. J., and Peter, M. A.: Dispersion relations, power laws, and energy loss for waves in the marginal ice zone, J. Geophys. Res.Oceans, 123, 3322–3335, 2018.
Montiel, F., Squire, V. A., Doble, M., Thomson, J., and Wadhams, P.: Attenuation and directional spreading of ocean waves during a storm event in the autumn Beaufort Sea marginal ice zone, J. Geophys. Res.Oceans, 123, 5912–5932, 2018.
Monteban, D., Johnsen, H., and Lubbad, R.: Spatiotemporal observations of wave dispersion within sea ice using Sentinel1 SAR TOPS mode, J. Geophys. Res.Oceans, 24, 8522–8537, https://doi.org/10.1029/2019JC015311, 2019.
Mosig, J. E., Montiel, F., and Squire, V. A.: Comparison of viscoelastictype models for ocean wave attenuation in icecovered seas, J. Geophys. Res.Oceans, 120, 6072–6090, 2015.
Newyear, K. and Martin, S.: Comparison of laboratory data with a viscous twolayer model of wave propagation in grease ice, J. Geophys. Res.Oceans, 104, 7837–7840, 1999.
Rabault, J., Sutherland, G., Gundersen, O., and Jensen, A.: Measurements of wave damping by a grease ice slick in Svalbard using offtheshelf sensors and opensource electronics, J. Glaciol., 63, 372–381, 2017.
Rabault, J., Sutherland, G., Jensen, A., Christensen, K. H., and Marchenko, A.: Experiments on wave propagation in grease ice: combined wave gauges and particle image velocimetry measurements, J. Fluid Mech., 864, 876–898, 2019.
Robin, G. D. Q.: Wave propagation through fields of pack ice, Philos. T. R. Soc. Lond. A, 255, 313–339, 1963.
Rogers, W. E., Thomson, J., Shen, H. H., Doble, M. J., Wadhams, P., and Cheng, S.: Dissipation of wind waves by pancake and frazil ice in the autumn Beaufort Sea, J. Geophys. Res.Oceans, 121, 7991–8007, 2016.
Rosenblum, E. and Eisenman, I.: Sea ice trends in climate models only accurate in runs with biased global warming, J. Climate, 30, 6265–6278, 2017.
Saha, S., Moorthi, S., Pan, H. L., Wu, X., Wang, J., Nadiga, S., Tripp, P., Kistler, R., Woollen, J., and Behringer, D.: The NCEP climate forecast system reanalysis, B. Am. Meteorol. Soc., 91, 1015–1058, 2010.
Saha, S., Moorthi, S., Wu, X., Wang, J., Nadiga, S., Tripp, P., Behringer, D., Hou Y. T., Chuang, H. Y., Iredell, M., and Ek, M.: NCEP climate forecast system version 2 (CFSv2) 6hourly products, Research Data Archive at the National Center for Atmospheric Research, Computational and Information Systems Laboratory, https://doi.org/10.5065/D61C1TXF, 2011.
Shen, H., Perrie, W., Hu, Y., and He, Y.: Remote sensing of waves propagating in the marginal ice zone by SAR, J. Geophys. Res.Oceans, 123, 189–200, 2018.
Shen, H. H., Ackley, S. F., and Hopkins, M. A.: A conceptual model for pancakeice formation in a wave field, Ann. Glaciol., 33, 361–367, 2001.
Smith, M. and Thomson, J.: Ocean surface turbulence in newly formed marginal ice zones, J. Geophys. Res.Oceans, 124, 1382–1398, 2019.
Smith, M., Stammerjohn, S., Persson, O., Rainville, L., Liu, G., Perrie, W., Robertson, R., Jackson, J., and Thomson, J.: Episodic reversal of autumn ice advance caused by release of ocean heat in the Beaufort Sea, J. Geophys. Res.Oceans, 123, 3164–3185, 2018.
Snyder, R. L., Dobson, F. W., Elliott, J. A., and Long, R. B.: Array measurements of atmospheric pressure fluctuations above surface gravity waves, J. Fluid Mech., 102, 1–59, 1981.
Spreen, G., Kaleschke, L., and Heygster, G.: Sea ice remote sensing using AMSRE 89 GHz channels, J. Geophys. Res., 113, C02S03, https://doi.org/10.1029/2005JC003384, 2008.
Squire, V. A., Dugan, J. P., Wadhams, P., Rottier, P. J., and Liu, A. K.: Of ocean waves and sea ice, Annu. Rev. Fluid Mech., 27, 115–168, 1995.
Squire, V. A.: Of ocean waves and seaice revisited, Cold Reg. Sci. Technol., 49, 110–133, 2007.
Squire, V. A.: A fresh look at how ocean waves and sea ice interact, Philos. T. R. Soc. A, 376, 20170342, https://doi.org/10.1098/rsta.2017.0342, 2018.
Squire, V. A.: Ocean wave interactions with sea ice: a reappraisal, Annu. Rev. Fluid Mech., 52, 37–60, https://doi.org/10.1146/annurevfluid010719060301. 2020.
Stopa, J. E., Ardhuin, F., Chapron, B., and Collard, F.: Estimating wave orbital velocity through the azimuth cutoff from spaceborne satellites, J. Geophys. Res.Oceans, 120, 7616–7634, 2015.
Stopa, J. E., Sutherland, P., and Ardhuin, F.: Strong and highly variable push of ocean waves on Southern Ocean sea ice, P. Natl. Acad. Sci. USA, 115, 5861–5865, https://doi.org/10.1073/pnas.1802011115, 2018a.
Stopa, J., Ardhuin, F., Thomson, J., Smith, M. M., Kohout, A., Doble, M., and Wadhams, P.: Wave attenuation through an Arctic marginal ice zone on 12 October 2015: 1. Measurement of wave spectra and ice features from Sentinel 1A, J. Geophys. Res.Oceans, 123, 3619–3634, 2018b.
Stroeve, J. and Notz, D.: Changing state of Arctic sea ice across all seasons, Environ. Res. Lett., 13, 10, https://doi.org/10.1088/17489326/aade56, 2018.
Thomson, J., Ackley, S., GirardArdhuin, F., Ardhuin, F., Babanin, A., Boutin, G., Brozena, J., Cheng, S., Collins, C., Doble, M., and Wadhams, P.: Overview of the arctic sea state and boundary layer physics program, J. Geophys. Res.Oceans, 123, 8674–8687, 2018.
Wadhams, P., Squire, V. A., Goodman, D. J., Cowan, A. M., and Moore, S. C.: The attenuation rates of ocean waves in the marginal ice zone, J. Geophys. Res.Oceans, 93, 6799–6818, 1988.
Wang, R. and Shen, H. H.: Experimental study on surface wave propagating through a grease–pancake ice mixture, Cold Reg. Sci. Technol., 61, 90–96, 2010a.
Wang, R. and Shen, H. H.: Gravity waves propagating into an icecovered ocean: A viscoelastic model, J. Geophys. Res.Oceans, 115, C06024, https://doi.org/10.1029/2009JC005591, 2010b.
WAVEWATCH III^{®} Development Group (WW3DG): User manual and system documentation of WAVEWATCH III^{®} version 6.07. Tech. Note 333, NOAA/NWS/NCEP/MMAB, College Park, MD, USA, 465 pp. +Appendices, 2019.
Weeks, W. F. and Assur, A.: The mechanical properties of sea ice, Cold regions research and engineering lab, Hanover, NH, 1967.
Zhao, X. and Shen, H. H.: Wave propagation in frazil/pancake, pancake, and fragmented ice covers, Cold Reg. Sci. Technol., 113, 71–80, 2015.
Zhao, X. and Shen, H. H.: A threelayer viscoelastic model with eddy viscosity effect for flexuralgravity wave propagation through ice covers, Ocean Model., 151, 15–23, 2018.