Reﬂection of Acoustic Wave through Multilayered Porous Sea Ice Sandwiched between the Water and Air Half-Spaces

: To establish an accurate sea ice model is a tremendous challenge in Arctic acoustic research. Regarding this matter, a multilayered porous sea ice model is proposed based on Biot’s theory in this paper. Assuming that the model is sandwiched between the water and air half-spaces, the reﬂection coefﬁcient of an incident wave from water into ice is deduced and contrasted with the solution calculated by impedance transfer method (ITM) to demonstrate the veriﬁcation of the model. Furthermore, the inﬂuences of frequency, porosity and layering on reﬂection coefﬁcients are analyzed. The results reveal that the reﬂection coefﬁcient is closely associated with layering and porosity. Therefore, it is reasonable and necessary to simultaneously take the layering and porosity of ice into consideration. Different from the existing layered or porous ice model, the presented model synthesizes the layered characteristic and porous structure of ice, which better portrays the real condition of sea ice. It is an improvement of the broadly used stratiﬁed or porous sea ice model, which provides ideas for further sea ice modeling.


Introduction
With the growing exploitation of the Arctic, extensive research into Arctic acoustics has been carried out in contemporary marine science and engineering fields [1]. Due to the seasonal existence of ice in most of the Arctic Ocean, it causes an effect on the propagation of acoustic waves [2]. Consequently, having a further understanding about sea ice models is essential. Arctic ice has a complicated structure because of the changeable environmental condition [3,4], which makes it challenging to establish a plausible model. When it comes to sea ice modeling, scholars conducted a study concentrating mainly on its three characteristics: layered, porous and rough surface. Chen [5] simplified the ice model by ignoring the surface roughness and regarded acoustic reflection from sea ice as from flat plate. What is more, they deduced the reflection coefficient in detail and analyzed the sound field propagation loss using the Bellhop simulation model. McCammon [6] put forward the multilayered ice model and extended the impedance transfer method to calculate the reflection coefficient at a water-ice interface, and then they determined the effects of the physical parameters and layering of ice and snow cover on reflections. Based on the layered model, Liu [7] assumed that the water-ice and ice-air interfaces are rough and used the self-consistent perturbation method proposed by Kuperman [8] to obtain the boundary conditions of rough surfaces to calculate the reflection coefficient. Then, a numerical analysis was conducted about the relations between reflection coefficients and ice thickness, grazing angle and frequency of acoustic wave. Kuperman treated rough surfaces with a Gaussian distribution rather than concentrating on the specific scale of ridges generated by collisions of large packs of floating ice floes [9]. Nevertheless, Diachok [10] calculated the reflectivity from ridges based on the Burke-Twersky model [11], in which ridges are considered as elliptical half-cylinders situated at the bottom surface of ice. Besides the layering and rough surface of ice, scholars took porosity into account in the process of modeling sea ice. Schwarz [12] modeled sea ice as a transversely isotropic brine-saturated porous medium, and the cylindrical pores align themselves vertically. In accordance with the model, Yew [13] studied the reflection and refraction of waves at the interface of water and porous sea ice and the effect of ice thickness and porosity of the skeleton layer on reflection.
Comparing these models, it can be found that most models are established based on one of the three characteristics, deviating from the actual condition of ice. Therefore, in an attempt to obtain precise results during analyses of transmission of acoustic signals in ice, a more comprehensive model is needed to incorporate those characteristics as much as possible.
Based on this situation, a multilayered porous sea ice model which combines the layered characteristic and porous structure is proposed in this paper. Due to the difference of the formed time of ice layers, the density and sound velocity vary in each layer [14], which has a different effect on the propagation waves. Moreover, when the plane wave propagates into pores within ice, it gives rise to vibration of the fluid trapped in pores, and thus its energy is dissipated on account of the friction between the fluid and ice frame [15]. Hence, it is reasonable to model sea ice as multilayered and porous. Figure 1 depicts the model in detail: sea ice is regarded as a multilayered porous medium embedded between the water and air half-spaces. The ice is divided into n layers whose physical properties are different vertically but the same horizontally [16], and lateral dimensions are assumed to be infinite. When an incident wave propagates at the interface between water half-space and the first layer of ice, it is reflected as an acoustic wave in water and gives rise to three transmitted waves [17]: a fast longitudinal wave L f , a transverse wave T and a slow longitudinal wave L s . When these three refracted waves impinge on the next boundary, these waves generate three reflected waves and three refracted waves, and this process continues at the following ice interfaces [18]. Finally, when reaching the ice-air interface, the longitudinal wave is transmitted into air half-space.

Theoretical Modeling
Appl. Sci. 2021, 11, x FOR PEER REVIEW 2 of 11 of ridges generated by collisions of large packs of floating ice floes [9]. Nevertheless, Diachok [10] calculated the reflectivity from ridges based on the Burke-Twersky model [11], in which ridges are considered as elliptical half-cylinders situated at the bottom surface of ice. Besides the layering and rough surface of ice, scholars took porosity into account in the process of modeling sea ice. Schwarz [12] modeled sea ice as a transversely isotropic brine-saturated porous medium, and the cylindrical pores align themselves vertically. In accordance with the model, Yew [13] studied the reflection and refraction of waves at the interface of water and porous sea ice and the effect of ice thickness and porosity of the skeleton layer on reflection.
Comparing these models, it can be found that most models are established based on one of the three characteristics, deviating from the actual condition of ice. Therefore, in an attempt to obtain precise results during analyses of transmission of acoustic signals in ice, a more comprehensive model is needed to incorporate those characteristics as much as possible.
Based on this situation, a multilayered porous sea ice model which combines the layered characteristic and porous structure is proposed in this paper. Due to the difference of the formed time of ice layers, the density and sound velocity vary in each layer [14], which has a different effect on the propagation waves. Moreover, when the plane wave propagates into pores within ice, it gives rise to vibration of the fluid trapped in pores, and thus its energy is dissipated on account of the friction between the fluid and ice frame [15]. Hence, it is reasonable to model sea ice as multilayered and porous. Figure 1 depicts the model in detail: sea ice is regarded as a multilayered porous medium embedded between the water and air half-spaces. The ice is divided into n layers whose physical properties are different vertically but the same horizontally [16], and lateral dimensions are assumed to be infinite. When an incident wave propagates at the interface between water half-space and the first layer of ice, it is reflected as an acoustic wave in water and gives rise to three transmitted waves [17]: a fast longitudinal wave Lf, a transverse wave T and a slow longitudinal wave Ls. When these three refracted waves impinge on the next boundary, these waves generate three reflected waves and three refracted waves, and this process continues at the following ice interfaces [18]. Finally, when reaching the ice-air interface, the longitudinal wave is transmitted into air half-space.  To better illustrate the inner structure of sea ice, the cross-section of the ith ice layer is chosen to portray it. It can be seen from Figure 2 that porous sea ice is comprised of an ice frame and tubular pores in which fluid is trapped. These pores are randomly distributed and approximately along the depth direction, so sea ice is regarded as a transversely isotropic fluid-saturated porous medium. To better illustrate the inner structure of sea ice, the cross-section of the ith ice layer is chosen to portray it. It can be seen from Figure 2 that porous sea ice is comprised of an ice frame and tubular pores in which fluid is trapped. These pores are randomly distributed and approximately along the depth direction, so sea ice is regarded as a transversely isotropic fluid-saturated porous medium.

Calculation of Displacements, Pressures and Stresses
To calculate the reflection coefficient of the plane wave from the model, first, the relations between displacement potentials and displacements, stresses or pressures in different media are deduced. Then, a set of simultaneous equations for the amplitudes of displacement potentials is derived using the boundary conditions of the displacements, stresses or pressures at each interface. Eventually, the reflection coefficient can be calculated by solving these equations.
When a plane wave with the angular frequency ω propagates from water into sea ice, the sound field in water contains the incident and reflected waves. However, in air, only the propagated upward wave exists because of its half-space. Therefore, the displacement potentials in water and air are as follows:  Then, the normal displacement uz0, uza and the pressure p0, pa [19] in water and air are given by the following: where ρ0 and ρa are the density of water and air. Taking the ith layer of sea ice as an example, the sound field contains three reflected and three refracted waves. Thus, the displacement potentials of the fast longitudinal wave, slow longitudinal wave and transverse wave can be written as follows:

Calculation of Displacements, Pressures and Stresses
To calculate the reflection coefficient of the plane wave from the model, first, the relations between displacement potentials and displacements, stresses or pressures in different media are deduced. Then, a set of simultaneous equations for the amplitudes of displacement potentials is derived using the boundary conditions of the displacements, stresses or pressures at each interface. Eventually, the reflection coefficient can be calculated by solving these equations.
When a plane wave with the angular frequency ω propagates from water into sea ice, the sound field in water contains the incident and reflected waves. However, in air, only the propagated upward wave exists because of its half-space. Therefore, the displacement potentials in water and air are as follows: where α 0 = k 2 0 − k 2 x , k x and k 0 refer to the vertical wavenumber, horizontal wavenumber and wavenumber in water, respectively; V is the reflection coefficient of the plane wave; x is the vertical wavenumber in air and k a is the wavenumber in air. Then, the normal displacement u z0 , u za and the pressure p 0 , p a [19] in water and air are given by the following: where ρ 0 and ρ a are the density of water and air. Taking the ith layer of sea ice as an example, the sound field contains three reflected and three refracted waves. Thus, the displacement potentials of the fast longitudinal wave, slow longitudinal wave and transverse wave can be written as follows: where x and k fl,i , k sl,i and k t,i are the vertical wavenumbers and wavenumbers of the three kinds of waves, respectively; i are the amplitudes of the three kinds of waves to be determined and the subscripts + and − signify waves propagating upward (i.e., in the positive z-axis) and downward, respectively. The normal displacements of the ice frame and the fluid in pores can be represented by displacement potentials as follows: where r i1 and r i2 represent the ratio of the velocity of the fluid in pores over the velocity of the frame for the fast longitudinal and slow longitudinal waves, and r i3 is the ratio of the amplitude of displacement of the fluid and of the frame. Then, according to the stress-strain relations in Biot's theory [20], the normal and tangential stresses of the ice frame and the pressure of fluid in pores can be obtained: where P, Q and R are elasticity coefficients and N is the shear modulus of the frame. For detailed derivation of these parameters and wavenumbers of three kinds of waves, the reader can refer to the paper by Allard [21]. As the same factor e j(ωt−k x x) is present in the expressions of displacement potentials in these media, it will be neglected in the following calculation.

Boundary Condition
Continuity relations of normal displacements and stresses exist at the interface between water and porous sea ice. Meanwhile, shear stress in the ice vanishes. Therefore, the boundary conditions [22] are as follows: where β is the porosity of the porous sea ice. When substituting the normal displacement and pressure of water in Equations (3) and (4) and Equations (8)- (13) into Equations (14)- (17), a series of equations about the unknown amplitudes of displacement potentials can be derived, which can be described in the form of matrices: , representing the vector of amplitudes of displacement potentials of water and the first layer of sea ice.
For simplicity of calculation, it is assumed that adjacent ice layers have the same porosity. Nevertheless, it is also applicable to adjacent layers with different porosities as long as the boundary conditions are changed. Additionally, the continuity of displacements, stresses and pressure of the ice frame and the fluid in pores at the ith ice-ice interface can be given as follows: where the subscript i means the ith ice-ice interface, as Figure 3 shows.
(17), a series of equations about the unknown amplitudes of displacement potentials can be derived, which can be described in the form of matrices: where B0(kx) is the 4 × 8 coefficient matrix, For simplicity of calculation, it is assumed that adjacent ice layers have the same porosity. Nevertheless, it is also applicable to adjacent layers with different porosities as long as the boundary conditions are changed. Additionally, the continuity of displacements, stresses and pressure of the ice frame and the fluid in pores at the ith ice-ice interface can be given as follows: ,, ,, where the subscript i means the ith ice-ice interface, as Figure 3 shows. Similarly, when substituting Equations (8)-(13) into Equations (19) and (20), the following can be derived: where Bi(kx) is the 6 × 12 coefficient matrix, ii  +− , representing the vector of amplitudes of displacement potentials of adjacent layers of sea ice. Finally, the boundary condition of the ice-air interface is similar to that of the waterice interface:  Similarly, when substituting Equations (8)-(13) into Equations (19) and (20), the following can be derived: where , representing the vector of amplitudes of displacement potentials of adjacent layers of sea ice.
Finally, the boundary condition of the ice-air interface is similar to that of the waterice interface: (1 − β)u zn + βW zn = u za (22) When substituting the normal displacement and pressure of air in Equations (3) and (4) and Equations (8)-(13) into Equations (22)-(25), the following can be obtained: , signifying the vector of amplitudes of displacement potentials of the nth layer of sea ice and air.
Based on the boundary conditions at all interfaces, a set of matrix equations about Equations (18), (21) and (26) is derived. There are a total of 4 + 6(n − 1) + 4 = 6n + 2 linear equations about amplitudes of displacement potentials. Solving these equations, 6n + 2 solutions are derived: Thus, the reflection coefficient of the plane wave is obtained.

Model Validation
In this section, the computed values of the reflection coefficient are compared with the results calculated by ITM [6] to demonstrate the validation of the model. These two methods can be set with the same parameters for comparison, which overcomes the diffi-culties of obtaining precise data from the changeable environment of the Arctic. The ITM is widely utilized to calculate the reflection coefficient of acoustic waves from multilayered materials. Therefore, the ITM was employed to verify the model. The principle of the ITM is to establish the relations of displacements and stresses between the first and the last layer. As for the displacement and stress field in the ith interface of ice, it can be described by the vector U(i): Furthermore, U(i) can be expressed by the amplitudes of displacement potentials: where A i is the matrix of order 6, signifying the coefficient matrix of the ith interface of ice, and its elements can be deduced by Equations (8)- (13). In the paper, double-layered porous sea ice is chosen as the model to simplify the calculation. Then, from the bottom surface of ice to the top surface, the transfer matrix G, which relates U(3) and U(1), is derived: where A 2d and A 2u are the coefficient matrix of the second interface in the first layer and the second layer of ice. Based on the boundary condition at the ice-air half-space interface along the condition p f 1 = −βp 0 and τ zx1 = 0, the transfer matrix H from the air to the upper ice layer can be obtained: where U(a) = [u a , p a ] T . Therefore, the transfer matrix from the bottom surface of the ice to the ice-air half-space interface is deduced: Furthermore, using the two remaining boundary conditions at the water-ice interface, the following can be derived: where U(0) = [u 0 , p 0 ] T and K is the transfer matrix. Solving Equation (34), the reflection can be calculated. For the double-layered model, the acoustic velocity and density are 1500 m/s and 1000 kg/m 3 in water and 343 m/s and 1.2 kg/m 3 in air, respectively. The density of the ice frame and the fluid in pores is 1.2 and 1000 kg/m 3 , respectively. The porosity of the two layers is 0.1, and the other parameters of porous ice are given in Table 1. For long-range propagation of acoustic waves in the Arctic, the transmission loss increases dramatically with the increase in frequency. Therefore, we focused on the low-frequency region below 2000 Hz. Using the parameters in Table 1, reflection coefficients as a function of incident angle for sound frequencies of 500, 1000, 1500 and 2000 Hz were calculated by solving the system of equations method (EM) and the ITM. To better distinguish our method from the ITM in the following figures, the former method is described as the EM. Figure 4 shows the comparisons between the EM and the ITM on the basis of the same double-layered porous ice model.
Curves of the reflection coefficient through EM and ITM nearly coincide with each other at different frequencies in Figure 4. The minor differences in those peaks are caused by the difference in accuracy of calculation between the two methods. Nevertheless, the similarity between the EM and ITM results is enough to confirm the validity of the model. Consequently, a conclusion can be drawn that the proposed double-layered porous model is valid, which simultaneously demonstrates the correctness of the multilayered and porous model. the ITM in the following figures, the former method is described as the EM. Figure 4 shows the comparisons between the EM and the ITM on the basis of the same doublelayered porous ice model. Curves of the reflection coefficient through EM and ITM nearly coincide with each other at different frequencies in Figure 4. The minor differences in those peaks are caused by the difference in accuracy of calculation between the two methods. Nevertheless, the similarity between the EM and ITM results is enough to confirm the validity of the model. Consequently, a conclusion can be drawn that the proposed double-layered porous model is valid, which simultaneously demonstrates the correctness of the multilayered and porous model.

Analysis of Factors Affecting the Reflection Coefficient
In Section 3, the double-layered porous ice model was verified to be correct, on the basis of which a further study was made on the factors influencing the reflection coefficient, such as frequency, porosity and ice thickness.

Analysis of Factors Affecting the Reflection Coefficient
In Section 3, the double-layered porous ice model was verified to be correct, on the basis of which a further study was made on the factors influencing the reflection coefficient, such as frequency, porosity and ice thickness.

Effects of Frequency on Reflection Coefficients
A map of the reflection coefficients is shown in Figure 5, mapped versus incident angle and frequency. The other parameters are the same as the model in Section 3. Obviously, the reflection coefficients approach 1 at 1-100 Hz. With the increase in frequency, there is a more rapid change in reflection coefficient. It is noted that frequency has an apparent effect on reflection coefficients, which decrease at first and then increase at small incident angles. As for larger incident angles, the reflection coefficient decreases with the increase in frequency and even reaches 1 between 82 • and 90 • . It can be explained that total reflection easily happens at those larger incident angles, in which case frequency is not the dominant factor. Moreover, the acoustic wavelength at low frequency is larger than the thickness of ice, so the wave is not trapped in pores without energy loss. Nevertheless, the wavelength is small at higher frequency and, thus, trapped in pores, which causes the loss of the incoming energy. Regarding the area in Figure 5 where the reflection coefficients approach 1 at 1000 Hz, it might be caused by other dominant factors, which needs further exploration.

Effects of Frequency on Reflection Coefficients
A map of the reflection coefficients is shown in Figure 5, mapped versus incident angle and frequency. The other parameters are the same as the model in Section 3. Obviously, the reflection coefficients approach 1 at 1-100 Hz. With the increase in frequency, there is a more rapid change in reflection coefficient. It is noted that frequency has an apparent effect on reflection coefficients, which decrease at first and then increase at small incident angles. As for larger incident angles, the reflection coefficient decreases with the increase in frequency and even reaches 1 between 82° and 90°. It can be explained that total reflection easily happens at those larger incident angles, in which case frequency is not the dominant factor. Moreover, the acoustic wavelength at low frequency is larger than the thickness of ice, so the wave is not trapped in pores without energy loss. Nevertheless, the wavelength is small at higher frequency and, thus, trapped in pores, which causes the loss of the incoming energy. Regarding the area in Figure 5 where the reflection coefficients approach 1 at 1000 Hz, it might be caused by other dominant factors, which needs further exploration.

Effects of Ice Thickness on Reflection Coefficients
Figure 6a depicts the reflection coefficients at different thicknesses of the first layer while keeping the second layer invariable. On the contrary, the second layer's ice thickness varies and the thickness of the first is constant in Figure 6b. An apparent phenomenon in these two figures is that the reflection coefficients decrease with the increase in thickness in most cases. This is because the length of tubular pores increases with the thickness, and the waves trapped in pores could be reflected more times, thus causing more energy losses. Furthermore, it can be seen in Figure 6a that the change in the first layer's ice thickness at small values has a minor effect on the reflection coefficient. Contrarily, Figure 6b shows that the reflection coefficient is deeply influenced by the second layer's ice thickness. To further illustrate this matter, a comparison was made between the double-layered porous model and the single-layered porous model at the same thickness in Figure 6c. The solid line represents the 1.5-m thickness of the single-layered model and the other two lines signify the 1-and 0.5-m thickness of the two layers in the double-layered model, in which the dashed line represents a first layer of 1 m and a second layer of 0.5 m and the dotted line is the opposite. As Figure 6c shows, these lines have the same trend to some degree but differences in detail. Consequently, the layering and the distribution of thickness of sea ice have an effect on the reflection coefficient, and it is necessary to establish the layered model.

Effects of Ice Thickness on Reflection Coefficients
Figure 6a depicts the reflection coefficients at different thicknesses of the first layer while keeping the second layer invariable. On the contrary, the second layer's ice thickness varies and the thickness of the first is constant in Figure 6b. An apparent phenomenon in these two figures is that the reflection coefficients decrease with the increase in thickness in most cases. This is because the length of tubular pores increases with the thickness, and the waves trapped in pores could be reflected more times, thus causing more energy losses. Furthermore, it can be seen in Figure 6a that the change in the first layer's ice thickness at small values has a minor effect on the reflection coefficient. Contrarily, Figure 6b shows that the reflection coefficient is deeply influenced by the second layer's ice thickness. To further illustrate this matter, a comparison was made between the double-layered porous model and the single-layered porous model at the same thickness in Figure 6c. The solid line represents the 1.5-m thickness of the single-layered model and the other two lines signify the 1-and 0.5-m thickness of the two layers in the double-layered model, in which the dashed line represents a first layer of 1 m and a second layer of 0.5 m and the dotted line is the opposite. As Figure 6c shows, these lines have the same trend to some degree but differences in detail. Consequently, the layering and the distribution of thickness of sea ice have an effect on the reflection coefficient, and it is necessary to establish the layered model.

Effect of Porosity on Reflection Coefficients
The relation between the reflection coefficient and the porosity of sea ice is shown in Figures 7 and 8. The porosities were chosen as 0.1, 0.2, 0.3 and 0.4. It is shown that there is a trend that reflection coefficients decrease with the increase in porosity, except at some special large angles. It is because sea ice contains more pores when the porosity increases. Then, waves are more easily trapped in many pores and absorbed during the process of propagating in the ice. Therefore, porosity should be included in sea ice modeling. Appl. Sci. 2021, 11, x FOR PEER REVIEW 9 of 11

Effect of Porosity on Reflection Coefficients
The relation between the reflection coefficient and the porosity of sea ice is shown in Figures 7 and 8. The porosities were chosen as 0.1, 0.2, 0.3 and 0.4. It is shown that there is a trend that reflection coefficients decrease with the increase in porosity, except at some special large angles. It is because sea ice contains more pores when the porosity increases. Then, waves are more easily trapped in many pores and absorbed during the process of propagating in the ice. Therefore, porosity should be included in sea ice modeling.
Based on these analyses, a conclusion can be drawn that frequency, layering and porosity have profound effects on the reflection coefficient. Therefore, a comprehensive consideration about the frequency of acoustic waves and ice properties is needed in Arctic acoustic research. Furthermore, it is necessary to take the layering and porosity characteristics into account in sea ice modeling.

Effect of Porosity on Reflection Coefficients
The relation between the reflection coefficient and the porosity of sea ice is shown in Figures 7 and 8. The porosities were chosen as 0.1, 0.2, 0.3 and 0.4. It is shown that there is a trend that reflection coefficients decrease with the increase in porosity, except at some special large angles. It is because sea ice contains more pores when the porosity increases. Then, waves are more easily trapped in many pores and absorbed during the process of propagating in the ice. Therefore, porosity should be included in sea ice modeling.
Based on these analyses, a conclusion can be drawn that frequency, layering and porosity have profound effects on the reflection coefficient. Therefore, a comprehensive consideration about the frequency of acoustic waves and ice properties is needed in Arctic acoustic research. Furthermore, it is necessary to take the layering and porosity characteristics into account in sea ice modeling.

Conclusions
A multilayered porous sea ice model has been proposed in this paper. The reflection coefficient of an incident wave from water into ice was deduced. Then, a simplified double-layered model was established for numerical calculation and verification. On the basis of the model, a further study was made on the factors influencing the reflection coefficients. It was found that the frequency, layering and porosity have a comprehensive effect on the reflection coefficient. Therefore, it is essential to simultaneously take layered and Based on these analyses, a conclusion can be drawn that frequency, layering and porosity have profound effects on the reflection coefficient. Therefore, a comprehensive consideration about the frequency of acoustic waves and ice properties is needed in Arctic acoustic research. Furthermore, it is necessary to take the layering and porosity characteristics into account in sea ice modeling.

Conclusions
A multilayered porous sea ice model has been proposed in this paper. The reflection coefficient of an incident wave from water into ice was deduced. Then, a simplified doublelayered model was established for numerical calculation and verification. On the basis of the model, a further study was made on the factors influencing the reflection coefficients. It was found that the frequency, layering and porosity have a comprehensive effect on the reflection coefficient. Therefore, it is essential to simultaneously take layered and porous structures into account in sea ice modeling. Different from the layered or porous models, the multilayered porous model fully combines the layered and porous characteristics of ice, which more closely approaches the real condition of ice, thus improving the accuracy of sea ice modeling to a great degree. In this model, the roughness of the water-ice interface is omitted. Nevertheless, the sea ice surface roughness affects the interactions between ice and acoustic waves to some extent [23]. Consequently, the roughness of the water-ice interface should also be included when modeling sea ice in subsequent research.