Bidimensional Ray Tracing Model for the Underwater Noise Propagation Prediction

Abstract

An increasing attention has recently been paid to the effect of the underwater noise field generated by ship activities on the marine environment. Although this problem is widely discussed in international treaties and conventions, it has not yet found a consolidated technical-scientific treatment capable of quantifying the level of underwater noise emissions produced by naval systems. As part of a national research collaboration, a novel code has been developed to predict noise propagation according to the Ray Tracing approach. Such optical geometry-based technique allows for calculating the Transmission Loss (TL) trend in its respective contributions: geometrical loss (due to the distance between the source and receiver), dissipation loss (due to the characteristics of the propagation environment), and reflection loss (due to the surfaces that delimit the field). The simulation requires as input parameters the source info as spatial position, frequency, and sound pressure level (SPL) as well as the sea properties like seabed depth, the speed of sound profile, the layers thickness the water column is divided into, the sea salinity, temperature, and pH. The simulation code provides the SPL spatial distribution useful as a fast industrial tool in the future studies addressed to identify the emission limits for the protection of marine wildlife.

1. Introduction

A topic of scientific interest for many years now is the study of noise propagated in the sea. As pioneering contributions in this sense, it is worth recalling the works of Knudsen [1] and Wenz [2], published at the beginning of the twentieth century. The sound disturbances emitted by naval systems are certainly among the main components of noise pollution—both sonic and ultrasonic—negatively affecting the marine ecosystem [3,4,5,6,7]. Awareness perspectives are supported by binding regulations of certain states or of EU Directives on the care of the natural environment, such as to impose limits on underwater noise levels [8]. The United Nations Convention on the Law of the Sea (UNCLOS) established underwater noise produced by both coastal and offshore anthropogenic activities as a source of pollution of the marine environment [9]. Article 1 defines the “pollution of the marine environment”. It is classified on the basis of continuous emissions (such as drilling, ship traffic) and intermittent (such as pile drivers, geoseismic surveys), which are added to the noise normally present in the marine environment. The characteristics of ship noise are a function of various parameters such as the type of vessel, its size, type of propulsion, and speed performance. The noise produced by the propeller in particular could have an emission range of up to hundreds of kilometers. In 2008, underwater noise was defined as a qualitative parameter of the state of health of the marine environment, thus requiring Member States to consider the problem and take precautionary measures aimed at containing cross-border pollution emissions (European Community, through the “Marine Strategy Framework Directive” [10]). In order to assess the extent of the environmental impact, national and international organizations were the first to issue procedures for measuring the noise radiated into the water by ships [11,12]. This, in order to regulate the levels of noise pollution deriving from navigation as the annexes of MAR.POL. (International Convention for the Prevention of Pollution from Ships [13]), only concerns forms of pollution by substances and not by forms of energy. Although there are still substantial gaps, a large growing research community supports the conclusion that noise pollution can severely affect the physiology, communication, behavior, and energy of marine animals with both acute and chronic effects as well [14,15,16,17,18,19,20,21,22]. Many studies on submarine transport noise were provided, highlighting important details on the noise levels of individual ships operating under various operating conditions. Research has always aimed to develop analytical or empirical data-based models for the diagnosis, identification, and classification of underwater acoustic wave sources. The noise produced by commercial ships dominates above all in the low frequency range of deep-ocean, but it remains more difficult to predict in coastal proximity where other local sound components intersect [3,23,24,25,26,27,28]. On the other side, war boats and submarines represent perhaps one of the few cases where noise levels are more limited, even if more for strategic reasons than ecological ones [29,30,31,32,33,34,35]. The typical noise components can be categorized into three main classes [8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32]: contributions due to the components of the hull, i.e., the powertrain and auxiliary systems; noise produced by the ship’s propellers; resonant phenomena of cavitation of the propellers and sustained turbulence around the underwater part of the hull. The ship’s service generator is the main source of underwater noise, especially of a tonal nature, generated by the ship at low speed: they are independent of the ship’s speed. The high-speed noise band is mostly given by some propulsive components with a frequency range above 100 Hz. Experimental field studies have been conducted in order to analyze the formation and propagation of the acoustic waves produced by a small ship in condition of the shallow sea [36]. The authors also evaluated the influence of boundary conditions, such as the seabed, with the noise spectrogram. Another experience concerns the evaluation of the noise radiated by commercial ships in transit on the coast of southern California (Santa Barbara Channel, SBC) [37]. Novel detection techniques of submarine sound levels as well as to identify the ship passages have been studied [38,39,40]. In the course of this article, a fast methodology for the prediction of the noise emitted by a source will be studied by applying the Ray Tracing theory, useful for calculating the wave propagation paths as a function of characteristic parameters: the different distribution speed of sound, absorption characteristics of the region fluid, and surfaces reflection properties. The Matlab^® code developed allows for simulating the mapping of hydrophones positioned at different heights of the seafloor: the focus of authors has been placed on the analysis of noise levels in the Adriatic Sea never presented in literature so far. This work in quality of “Case Report” does not aim to introduce an innovative tool in the study of underwater noise but rather to present its potential in a preliminary project phase between the University of Naples “Federico II” and Fincantieri Company. The results were compared with the Pekeris model [41]: a maximum error of 3 dB was found both at the assigned heights and up to a range of 25 km.

2. Materials and Methods

2.1. Theoretical Model of Ray Tracing

The development of predictive tools can also prove extremely useful in order to optimize the noise levels emitted at sea: the mathematical models are based on the phenomenological laws that govern the physics of wave propagation in fluids. The “Ray Tracing theory” represents a possible modeling approach for this purpose. This type of modeling is used in the optical field to analyze the light diffusion, but its applicability can also be extended to the acoustic field. This technique is a high-frequency approximation of the solution of the wave equation. It allows to define the propagation of sound through the emission of rays emitted by a point source in a medium at a constant speed of sound. The rays emitted by the source are perpendicular to the wave fronts, which are instead spherical surfaces radiated by the source itself. The reference theory is due to the Snell–Descartes law, which is used in acoustic geometry to describe the interference phenomena that occur between two homogeneous fluids having different speeds. Let us consider two fluids of which the first has density ρ₁ and speed c₁ and the second which has density ρ₂ and speed c₂ (with c₂ > c₁), as represented in Figure 1.

Snell–Descartes Theorem. The ratio of the cosines of the angles of incidence and refraction is equivalent to the ratio of phase velocities in the two media, or equivalent to the reciprocal of the ratio of the indices of refraction.

The difference in the velocity of the fluids means that the incident ray (inclined by the angle β₁ with respect to the horizontal line) partially reflects upwards symmetrically and partially refracts in the second fluid with an angle β₂ that is obtained from the following relationship (1):

$$\frac{cos{\beta }_1}{c_1}=\frac{cos{\beta }_2}{c_2}$$

(1)

This equation is valid only if β₂ ≤ 1 and therefore cosβ₁ ≤ c₁/c₂.

From this relation, it is possible to introduce another angle β_critic called “critical angle”, equal to Formula (2):

$${\beta }_{critic}=\arccos \left(\frac{c_1}{c_2}\right)$$

(2)

For incidence angles β₁ smaller than β_critic, there is total reflection. The same law can also be extended to study sound propagation phenomena when velocity varies linearly with depth.

Considering the law of variation of speed with depth as per Equation (3) (see Figure 2):

$$c\left(z\right)=c_0+g\left(z-z_0\right)$$

(3)

In which z is the vertical dimension and g is the speed gradient. Substituting the Equation (3) in the Snell–Descartes law Equation (1), the Equation (4) is obtained:

$$cos\beta \left(z\right)=\frac{c\left(z\right)}{c_0}cos{\beta }_0=\left(1+\frac{g}{c_0}\left(z-z_0\right)\right)cos{\beta }_0$$

(4)

Given the radius R_C = c₀/gcosβ₀, the relationship (4) is reduced as follows:

$$cos\beta \left(z\right)=cos{\beta }_0+\frac{z-z_0}{R_C}$$

(5)

It is observed, therefore, that in the case in which the speed varies with depth, the rays do not maintain straight trajectories but bend defining curved sections, unlike the case in which the speed remains constant and the rays propagate according to rectilinear sections. The method therefore consists in determining of the propagation of sound through linear segments and arcs that are connected along the water column, Figure 2. Relying upon this physical principle, the numerical prediction has been developed in this paper.

2.2. Numerical Prediction Code: Objectives and Functionality of the Program

The simulation program developed and described within this paper implements the basic foundations of this theory. The main purpose was to evaluate the Transmission Loss (TL) in water, and therefore the calculation of the sound pressure level (SPL) of a given point-like source that radiates a sound wave. The phases necessary to achieve this goal will be explained below, giving the reader the opportunity to get the program steps and then describing the applicative example about the noise analysis of the Adriatic Sea (Italy). From a general standpoint, the simulation code is structured in two distinct parts in order to guarantee greater speed and fluidity. The first part is dedicated to the input of the characteristic parameters of the source (the minimum angle, the maximum angle, and the step of emission of the rays) and of the surrounding environment (depth, temperature, salinity, pH). These data are fundamental to evaluate the water absorption coefficient α (Appendix A). Subsequently, the program, through a series of iterative cycles, will provide the different propagation paths of the rays emitted by the source. In the second part of the program, the user has the possibility to set up a matrix of virtual receivers through which to map the trend of the Transmission Loss (TL) and the SPL sound pressure level. The global layout is represented in Figure 3.

2.2.1. Propagation Paths

In evaluating the propagation paths, it was assumed that the velocity varies linearly with the depth through a law of this type (6):

$$c=c_1+g_1\times z_n$$

(6)

where c₁ is the surface speed, g₁ the speed gradient, and z_n the various levels in which the water column is subdivided. The variation of the incident ray along the water column is evaluated as a function of the velocity along the different levels.

According to the Snell–Descartes law at each level:

-

If the angle of the i-th ray is greater than the critical angle, the latter will continue to propagate according to its direction until it reaches the bottom or the surface of the sea and then will reverse its direction.

-

If the angle of the i-th ray is greater than the critical angle, the latter will obviously reverse its direction of propagation before reaching the bottom.

Finally, one set of rays will reverse their direction when they interfere with the sea floor and surface while other rays will reverse their direction before reaching the bottom, tracing arch-shaped trajectories, Figure 4a. Once the different propagation paths have been obtained, the second part of the program requires the definition of the map of receivers such as to discretize the fluid domain, Figure 4b.

2.2.2. Calculation of the Loss of Propagation (TL)

When a wave propagates from a source, a reduction in intensity due both to the absorption of acoustic energy and to a loss of geometric diffusion (Appendix B) could be considered. Absorption was taken into account for the initial stage of the program, aiming to calculate the transmission loss ratio as follows in the Equation (7):

$$TL=20logR+ \alpha R$$

(7)

In the case of rays incident on the sea floor, an absorption coefficient equal to 0.5 has been assumed, determining, therefore, an increase of 3 dB in the propagation loss according to Equation (8):

$$TL=20logR+ \alpha R+3n$$

(8)

where n is the number of reflections on the sea floor.

2.2.3. Calculation of the Sound Pressure Level (SPL)

The SPL magnitude with respect to the acting pressure p is defined as the following quantity (9):

$$L_P=10{\log }_{10}{\left(\frac{p}{p_{ref}}\right)}^2 \left[dB\right]$$

(9)

in which p_ref represents the reference pressure evaluated at 1 m from the source and equal to 1 μPa. In order to evaluate the SPL due to the various sound fields, the linear sum of the pressures is automatically performed in Matab^® (The MathWorks, Inc. Natick, MA, USA), assuming that:

$${\left(\frac{p_i}{p_{ref}}\right)}^2={10}^{SP{L}_i/10} \left[Pa\right]$$

(10)

where p_i is the i-th sound field pressure calculated at the receiver. The total sound field is obtained from the sum of the contributions of all pressures as per Equation (11):

$$SP{L}_{TOT}=10{\log }_{10}\left(\frac{{{\displaystyle \sum }}_i{p}_i}{p_{ref}}\right) \left[dB\right]$$

(11)

2.3. Case Study Overview

As part of a collaboration between the University of Naples “Federico II” and Fincantieri technical department, the following program was developed to evaluate the propagation of sound in the Adriatic Sea, from a source located at a depth of 10 m from the sea surface. Five case studies were carried out in the same hydrogeological site, considering different depths with the purpose of making a final comparison of the different modes of transmission. The parameters used during the application and the expected results are outlined below.

• SOURCE:

A point source located at –10m from the sea surface and with an SPL equal to 100 dB was assumed as a source of disturbance. With a range of rays varying between −90° and +90° with a step of 0.5°.

• SEA AND CHARACTERISTICS:

The Adriatic Sea, as mentioned, was the reference marine environment for this research. This sea develops for a length of 800 km and a width of 200 km. This sea can be divided into three characteristic areas, Figure 5a:

-

A northern area characterized by a low percentage of salt due to the fact that the rivers flow into the basin. This area, which extends up to the cross of Ancona, is characterized by a shallow seabed.

-

A central area that goes from Ancona to the cross of the island of Pianosa, has a higher depth, and in particular, this area is characterized by warmer waters exhaled by the presence of the Ionian currents that penetrate through the Otranto Canal.

-

A southern area that extends along the coasts of Puglia where the maximum depth foreseen for this sea is reached, which is 1233 m in correspondence with the Gargano. This depression is also known as the “Fossa of the lower Adriatic”.

Based on this consideration, five case studies were elaborated by changing the depth of the seabed: from an average value of 250 m up to the maximum depth of 1250 m, with steps of 250 m assumed as the spatial range of investigation. So, the following cases were contemplated: 250 m, 500 m, 750 m, 1000 m, 1250 m. The characteristics of the sea are those reported in the MARCOAST (Marine and Coastal Environmental Information Services) project with which ISPRA (Higher Institute for Environmental Protection and Research) collaborated [42]. The project carried out a monitoring of the characteristics of the Adriatic Sea from 2007 to 2011. For the selection of the parameters, the average surface temperature measured in July 2008 (T = 24 °C) was considered, Figure 5b, while for the salinity, an average value based on ten-year surveys was chosen. The salinity chosen for the application is 37.62 psu (Practical Salinity Units), Figure 5c.

• MATRIX OF RECEIVERS:

For the second part of the program, a matrix of 6 × 6 receivers distributed at intervals of equal width was set up as reported in Figure 4b.

3. Results

Following the logic of the previous paragraphs and exploiting the input data shown, the program will be able to graphically provide the different sound fields. At the end of processing, the program will graphically provide a series of sound fields. The first to be provided are: the direct sound field, the sound field after reflection on the free surface, the sound field after reflection on the seabed. Since sound can be composed, it is observed that the compound sound fields are obtained from the sum of direct and reflected sound fields. In fact, the program provides: the sound field composed between direct and reflected on the free surface, the sound field composed between direct and reflected on the bottom. Finally, the total one is obtained from the composition of the three sound fields, that is, total sound field composed of direct, reflected on the free surface, and reflected on the bottom. This reasoning was carried out for each case study and therefore at each depth. The simulations have been performed using an Intel^® Core™ i7-6500U processor @2.50 GHz with 2 cores. These sound fields (direct, reflected, and combined) obtained for the particular case of maximum height of 250 m are represented in Figure 6. The composed acoustic contours representative of the other heights from the seabed, i.e., 500 m, 750 m, 1000 m, and 1250 m, are reported in Figure 7. The simulations show the importance of considering the reflected fields for the global propagation loss quantification. The global sound field in Figure 6f actually follows almost the same wave pattern as the direct field in Figure 6a but higher in amplitude: the resulting field reaches up to 55 dB against the 50 dB of the direct field only. As the coverage range increases, the SPL decays with a logarithmic law and as a function of the attenuation parameters illustrated above in paragraph 2: the acoustic level reaches a maximum of 40 dB at the distance of 1250 m, Figure 7d.

4. Discussion

In order to assess the Ray Tracing trend accuracy, the transmission loss has been compared with the Pekeris waveguide-based model; the purpose of such investigation was to deepen the basic understanding of the propagating mechanisms of underwater noise with respect to well-known theoretical schemes, especially for far-field cases. A maximum range of 25 km has been assumed. Below are the graphs representing a final comparison among the case studies, Figure 8.

The Pekeris model represents the simplest one to characterize shallow water acoustics, considering calm sea, constant depth, homogeneous medium non-absorbent seabed. The analyses showed an acceptable level of correlation between the Ray Tracing model, whose decay is a logarithmic function, and the transcended Pekeris model based instead on the propagation of normal modes; an average error of about 3 dB was estimated between the two approaches, which is acceptable in a preliminary design phase. For Pekeris, the rate of dispersion and waves attenuation in the waveguide depends on the iterative search of mode numbers. The wave equation is solved to determine the normal mode solutions which prevail at large distances from the source, according to Lamb procedure. For the present calculation, four normal modes have been supposed: the pressure field is calculated according to a two-dimensional modal extension (Kirchhoff equation). The discrepancy between two approaches could be found, therefore, in the different schematization of reflection and attenuation mechanisms, with the TL increasing with mode number [41]. Due to the broadband nature of the problem, the model is very computationally demanding. The calculation times essentially depend on the pitch of the rays and the distance to be covered. Performance can certainly be optimized by consolidating the simulation code structure. However, for a preliminary phase, these times are still advantageous with respect to finite element models, which required much more computational cost and mesh modeling steps as well. At the end of the processing, the time required for each simulation was evaluated: it is the function of two main parameters:

• HEIGHT OF THE SEABED

The higher the height of the bottom, the more the processing time required for the simulation increases.

• STEP BETWEEN RADIATION ANGLES

The code allows for inserting the interval between the emission rays of the source; the smaller the step, the more detailed the noise map will be. The computational cost varies quadratically with the depth of the sea, Figure 9. In

Table 1. Computational cost estimation as function of sea depth.

SIMULATIONS	DEPTH	TIME REQUIRED
CASE I	250 m	9.86 min
CASE II	500 m	30.47 min
CASE III	750 m	55.30 min
CASE IV	1000 m	92.05 min
CASE V	1250 m	137.54 min

, the times required by the program in the five cases.

5. Conclusions

This paper is part of the macro-topic inherent in the current need to study forecasting systems for optimizing the sources of underwater noise pollution generated by ship traffic. In particular, the activity arises from a collaboration between the University of Naples “Federico II” (industrial engineering department) and the technical office of Fincantieri Group with the purpose of creating a fast predictive code for a preliminary evaluation of the acoustic propagation fields in submarine environments, with special attention to the basin of the Adriatic Sea. The simulation program in Matlab^® is essentially based on the optical-acoustic technique of Ray Tracing, which allows to characterize the sound propagation path also as a function of the characteristic parameters of the source and of the diffusion medium. Methodology is very suitable for preliminary calculations for two reasons:

• For the identification of the sound field radiated even at a great distance and, therefore, of the impact this may have on the marine environment.
• For the interpretation of the possible experimental values detected by the hydrophones in a field relatively close to the ship for an effective characterization of the disturbance source.

Future developments conducted by the research group aim at expanding the functionality of the program, taking into account also the third spatial dimension, while always limiting the computational cost, which, for applications of this kind, can be excessively expensive. Experimental campaigns of field measurements are also underway to validate this first part of the program.