An ocean acoustical ray tracing tool based on Fermat’s least time principle in a real environment

Gonçalves, João 1

CINAV (Navy Research Center) / Escola Naval (Navy Academy) 2810-001 Almada, Portugal

ABSTRACT Using an ocean-acoustic parametric ray tracing propagation models in 2D and 3D, based on Fermat’s least time principle, two programs were developed with the ability to study sound propagation in a real ocean environment. REOART 2D and REOART 3D use three di ff erent databases to withdraw information regarding the ocean topography and climatology to compute the acoustic ray trajectories in a real environment.

KEYWORDS: Ray-tracing, Ocean-acoustic propagation, Simulation in a real environment

1. INTRODUCTION

Over the last four decades, underwater acoustics, and more precisely, ocean-acoustic propagation models, have been intensely researched in response to practical requirements [1]. As a result, numerous types of models have been developed, including ray tracing, parabolic equation (PE), and normal mode. The development and application of techniques for forecasting sound propagation trajectories in the ocean is critical for military and societal purposes. It enables, for example, the military to pinpoint places with a low probability of detection. It enables, for instance, the development and application of acoustical technologies for imaging and measuring underwater structures and oceanic parameters in the civil field. The development of underwater acoustic models permit the simulation of sonar performance under laboratory conditions. Therefore, these models are now used to predict sonar performance at sea and designing optimized sonar systems [2]. In this work, two ocean-acoustic propagation programs were developed — REOART 2D and REOART 3D ( R eal E nvironment O cean- A coustic R ay T racing). These programs use an ocean-acoustic parametric ray-tracing propagation models in 2D and 3D, based on Fermat’s least time principle. These models, constituted by di ff erential algebraic system of equations (DAE’s), are very robust and reliable — see [3] and [4]. Through these ray tracing propagation models and utilizing the GEBCO ( General Bathymetric Charts of the Ocean ), WOD18 ( World Ocean Database 2018 ) and TEOS-10 ( Thermodynamic Equation of State 2010 ) databases, the REOART 2D and REOART 3D are capable of studying and forecasting sound propagation in a real ocean environment.

1 joaopdgoncalves@gmail.com

a slaty. inter.noise 21-24 AUGUST SCOTTISH EVENT CAMPUS O ¥, ? GLASGOW

2. RAY TRACING PROPAGATION MODELS

We present below the acoustic propagation models used in this work deduced by Raposo [4] [3]. Here is the 2D model 

x ′ − a = 0 y ′ − b = 0 a ′ b − b ′ a + b v ( x , y ) ∂ v ∂ x − a v ( x , y ) ∂ v ∂ y = 0 aa ′ + bb ′ = 0

, (1)

and here is the 3D model, 

x ′ − a = 0 y ′ − b = 0 z ′ − c = 0 (1 − a ′ ) a ′ + (1 − c 2 ) c ′ − ac ( a ′ + c ′ ) − b ( ab ′ + cb ′ ) + 1

, (2)

v  ∂ v ∂ x + ∂ v

∂ z  dv dx a + dv

dz c  ( a + c )  = 0 (1 − b ′ ) b ′ + (1 − c 2 ) c ′ − bc ( b ′ + c ′ ) − a ( ba ′ + ca ′ ) + 1

dy b + dv

v  ∂ v ∂ y + ∂ v

∂ z  dv dx a + dv

dz c  ( b + c )  = 0 aa ′ + bb ′ + cc ′ = 0

dy b + dv

where x , y , z and v stand, respectively, for the ray positioning and the speed of sound. The numerical resolution of the DAE’s (1) and (2) presents di ffi culties that require the use of appropriate solvers [5], such as the algorithm developed in MATLAB by Tony Roberts [6], which it was used in this work.

3. 2D AND 3D FORMULATION IN REAL ENVIRONMENT

The majority of propagation models assume that sound propagates in a two-dimensional plane. That is not always true. This is not the case when sound speed presents spatial gradients perpendicular to the propagation plane. In this scenario, the propagation takes place in a three-dimensional environment. Additionally, when a ray reflects in the bottom of a 2D plane, it may exit the initial 2D plane [7]. For these reasons, a three-dimensional propagation model must be employed to adequately understand sound propagation. Under the 2D acoustic propagation, hypothesis it is supposed that an acoustic ray propagates in the ocean in a vertical plane that intersects the surface geodesically. Therefore, the acoustic ray’s azimuth along the plane will not be constant along its path. The study of sound propagation will be conducted in this work up to distances of approximately 80 kilometers. To do so, we can utilize the transport formulas for geodesic coordinates adopted from Barahona Fernandes’ Hydrography Manual [8] to determine the variation of latitude ( ∆ ϕ ), the variation of longitude ( ∆ L ), and the meridian convergence ( ∆ α ). Both programs — REOART 2D and REOART 3D — use three di ff erent databases to withdraw information regarding the topography and climatology to compute the acoustic ray trajectories in a real environment. These databases are the GEBCO, the WOD18 and the TEOS-10 software library.

3.1. GEBCO bathymetric data The data provided by GEBCO enables to describe the bathymetry of the ocean almost worldwide. In this work it as used the latest version (2021) of the bathymetric data relative to the North Atlantic

Ocean. In Figure 1 it is possible to see a part of the database through the representation of the area surrounding the United Kingdom.

Figure 1: Bathymetric map of the area surrounding the United Kingdom.

3.2. WOD18 climatology data The WOD18 is a NOAA ( National Oceanic and Atmospheric Administration ) database that contains a variety of data types relevant to oceanography. Nonetheless, we calculated the speed of sound using only the practical salinity and in situ temperature. The points with climatological data (practical salinity and in situ temperature) in the area surrounding the United Kingdom are depicted in Figure 2. These points are separated by a quarter of a degree.

3.3. TEOS-10 The TEOS-10 software library contains numerous oceanographic functions for computing relevant ocean variables. To calculate the sound speed, the Gibbs-SeaWater (GSW) Oceanographic Toolbox was used, which is a subset of this library. Nevertheless, there are several irregularities when utilizing the function that calculates the speed of sound. To begin, practical salinity and in situ temperature cannot be employed, as TEOS-10 operates on an absolute salinity and conservative temperature basis, as demonstrated in [9]. Second, the function used to determine the sound speed requires three variables: absolute salinity, conservative temperature, and hydrostatic pressure. This necessitated the usage of the gsw_SA_from_SP and gsw_CT_from_t functions to convert practical salinity to absolute salinity and in situ temperature to conservative temperature, respectively. Additionally, we calculated the hydrostatic pressure using the gsw_p_from_z function. Finally, with the three variables defined, we can use the gsw_sound_speed function to calculate the sound speed and define the sound velocity field.

4. REOART 2D AND REOART 3D

As previously stated, these programs —- REOART 2D and REOART 3D — can simulate in mode 2D or 3D acoustic propagation in a real ocean environment. As already mentioned, these programs use the

Figure 2: Climatological data of the area surrounding the United Kingdom.

ray tracing propagation models developed by Raposo [4] and the databases mentioned before. From geographical coordinates, time and simulation area, REOART 2D and REOART 3D , figure out information about the topography and climatology, and compute the corresponding grid of spatial sound speed and all the necessary gradients used by (1) and (2). These programs are distinguished from others by their capacity to calculate successive ray reflections in the bottom and sea surface, as well as their use of real bathymetry and sound speed.

5. NUMERICAL SIMULATIONS

In this section we will carry out preliminary simulations in order to test the developed programs. All simulations will be performed in the area surrounding the United Kingdom. In Section 5.1, we will illustrate the emission of a single acoustic beam using the REOART 2D program. In Section 5.2, we will simulate the emission of multiple acoustic rays with di ff erent initial directions. In Section 5.3, we will simulate the emission of a single acoustic ray using the REOART 3D program.

5.1. REOART 2D — single ray simulation The first simulation was performed at geographic coordinates 50 ◦ N and 15 ◦ W on azimuth 270 ◦ — note that geographic coordinates to the west have a negative sign. We can observe the bathymetry of the area surrounding the United Kingdom and the geodesic (represented by the red line) in Figure 3. The geodesic is the 2D propagation plane and, as seen, it cannot extend beyond 80 kilometers. The bathymetry of the 2D propagation plane is depicted in Figure 4. The following step is to create the sound velocity field. Similarly to the bathymetry determination, we employed the same 2D propagation plane to define the plane’s climatology. As a result, this enables the establishment of the sound speed grid and the corresponding gradients. Two figures are showed in Figure 5. On the left, a similar figure to Figure 2 can be seen, but with the addition of the geodesic representation and some green points. These are the points with climatological data closest to the trajectory. The second one, on the right, is an approximation of the

Figure 3: Bathymetric map of the area surrounding the United Kingdom and the representation of the geodesic (2D propagation plane).

0 Bathymetry from lat. 50.0 and long. -15.0 on azimuth 270

-500

-1000

-1500

Depth [m]

-2000

-2500

-3000

-3500

-4000

-4500

Range [km] 0 10 20 30 40 50 60 70

Figure 4: Representation of the bathymetry of the geographic coordinates 50 ◦ N and 15 ◦ W on azimuth 270 ◦ .

left figure. The requirement for interpolation of climatological data throughout the trajectory is clearly visible. To do this, the concept of the nearest neighbor was applied. Thus, the climatological data were associated with the trajectory’s nearest neighbor. As a result, it was possible to determine the climatology in depth along the trajectory. After defining the climatology in the 2D propagation plane, we can finally compute the sound velocity field using the previously described TEOS-10 functions. Figure 6 shows a simulation of the acoustic ray propagation at the geographic coordinate 50 ◦ N and

Map with the distribution of points with the climatological data of the area surrounding the United Kingdom

50.4

60

50.3

58

50.2

56

50.1

Latitude [º]

Latitude [º]

50

54

49.9

52

49.8

50

49.7

49.6

48

-20 -18 -16 -14 -12 -10 -8 -6 -4 -2 0

-16.6 -16.4 -16.2 -16 -15.8 -15.6 -15.4 -15.2 -15 -14.8 -14.6

Longitude [º]

Longitude [º]

Figure 5: Representation of the map with the distribution of points with the climatological data of the area surrounding the United Kingdom, the geodesic and the nearest neighbors (the green points).

15 ◦ W on azimuth 270 ◦ at a depth of zero meters. We can observe the ray refracting along its path of propagation and reflecting in the sea surface in this figure. As is well known, the sound is described as "lazy". This means that the sound bends towards the regions of slower sound speed. As illustrated in Figure 7, this occurred during the simulation.

Acoustic ray propagation

0

500

1000

1500

2000

Depth [m]

2500

3000

3500

4000

4500

Range [m] × 10 4 0 1 2 3 4 5 6 7 8

Figure 6: Representation of the simulation of a single ray in the geographic coordinates 50 ◦ N and 15 ◦ W on azimuth 270 ◦ at a depth of zero meters.

Acoustic ray propagation

0

1514

1512

500

1510

1000

1508

Sound speed [m/s]

1506

1500

Depth [m]

1504

2000

1502

1500

2500

1498

3000

1496

1494

Range [m] × 10 4 0 1 2 3 4 5 6 7

Figure 7: Representation of the simulation of a single ray in the geographic coordinates 50 ◦ N and 15 ◦ W on azimuth 270 ◦ at a depth of zero meters and the sound velocity field.

5.2. REOART 2D — multiple rays simulation The second simulation was performed at the geographic coordinates 48 ◦ N and 9 ◦ W on azimuth 200 ◦ using multiple rays at zero meters deep. Obtaining the bathymetry and climatology of the 2D propagation plane is exactly the same as in the previous simulation. We will determine the bathymetry and then the climatology through the closest climatological points to the propagation plane. In Figure 8 it is possible to see the bathymetry of the area surrounding the United Kingdom as the 2D propagation plane (the red line). In Figure 9 we can see the real bathymetry of the 2D propagation plane. Figure 10 illustrates how the nearest neighbor concept is used to establish the climatology of the 2D propagation plane (the green points are the neighbors). The multiple ray simulation was performed at geographic coordinates 48 ◦ N and 9 ◦ W on azimuth 200 ◦ (Figure 11). Five rays were used in this simulation, each five degrees apart vertically, commencing at zero (0 ◦ − 5 ◦ − 10 ◦ − 15 ◦ − 20 ◦ ), at a depth of one meter.

5.3. REOART 3D — single ray simulation For the final simulation, we used REOART 3D to do a 3D simulation. This simulation was done at the geographic coordinates 53 ◦ N and 13 ◦ W with a single ray. REOART 3D does not permit degrees when calculating acoustic propagation. It is only e ff ective in meters. As a consequence, the axes will be measured in meters rather than degrees from our starting geographic position, along longitude and latitude. In Figure 12 it is possible to see the 3D simulation area (the red square). Note that geographic position 53 ◦ N and 13 ◦ W is the westernmost lower vertex. This position is our starting position. The bathymetry of the 3D simulation area is shown in Figure 13. As with REOART 2D , the 3D area’s climatology is determined using the nearest neighbor concept (Figure 14). The green points represent the neighbors. The simulation was performed at a depth of 350 m, 8x10 4 m from the initial geographic coordinates (53 ◦ N and 13 ◦ W) measured from longitudes and latitudes. In Figure 15 it is possible to see the acoustic ray refracting and reflecting in the bottom.

Figure 8: Bathymetric map of the area surrounding the United Kingdom and the representation of the geodesic (2D propagation plane — the red line).

Bathymetry from lat. 48.0 and long. -9.0 on azimuth 200

0

-500

-1000

-1500

Depth [m]

-2000

-2500

-3000

-3500

-4000

0 10 20 30 40 50 60 70

Range [km]

Figure 9: Representation of the bathymetry of the geographic coordinates 48 ◦ N and 9 ◦ W on azimuth 200 ◦ .

6. CONCLUSIONS

The aim of this work was to test preliminary REOART 2D and REOART 3D programs and assess the credibility of the simulations carried out. Not only do the trajectories seem credible, but also the reflections on the surface that at the bottom of the ocean seem to satisfy Snell’s law.

2unding the United Kingdom

Figure 10: Representation of the map with the distribution of points with the climatological data of the area surrounding the United Kingdom, the geodesic and the nearest neighbors.

Acoustic ray propagation

0

Latitude [) 4644 ~«422«40~=«Ssti( KD Longitude Latitude () 478. ara

500

1000

1500

Depth [m]

2000

2500

3000

3500

4000

4500

Range [m] × 10 4 0 1 2 3 4 5 6 7 8

Longitude [')

Figure 11: Representation of the simulation of the multiple rays, with five degrees apart vertically, in the geographic coordinates 48 ◦ N and 9 ◦ W on azimuth 200 ◦ at a depth of one meter.

7. FUTURE WORK

The preliminary satisfactory results suggest the following developments:

– Validate the ray trajectories predicted by REOART 2D using classic ray tracing software under the same real sound speed grid;

58. 56 Map with the distribution of points with the climatological data of the area surr 43.4

Figure 12: Bathymetric map of the area surrounding the United Kingdom and the representation of the 3D simulation area (the red square).

3D Bathymetry

-200

-300

-400

-500

Depth [m]

-600

-700

-800

-900

-1000

53

-11.6 -11.8 -12 -12.2

52.4 52.5 52.6 52.7 52.8 52.9

-12.4 -12.6 -12.8

-13 52.1 52.2 52.3

Latitude [º]

Longitude [º]

Figure 13: Bathymetry of the 3D simulation area.

– Improve the e ffi ciency of the programs developed, namely the REOART 3D program.

ACKNOWLEDGEMENTS

Firstly, I would like to express my very great appreciation to my thesis advisor Miguel Moreira, for his valuable and constructive suggestions during the planning and development of this research work.

Map with the distribution of points with the climatological data of the area surrounding the United Kingdom

53.6

58

53.4

53.2

56

53

52.8

Latitude [º]

Latitude [º]

54

52.6

52.4

52

52.2

52

50

51.8

48

51.6

-16 -14 -12 -10 -8 -6 -4 -2 0 2

-14.5 -14 -13.5 -13 -12.5 -12 -11.5 -11 -10.5

Longitude [º]

Longitude [º]

Figure 14: Representation of the map with the distribution of points with the climatological data of the area surrounding the United Kingdom, the 3D simulation area and the nearest neighbors.

Acoustic ray propagation

-200

-300

-400

-500

Depth [m]

-600

-700

-800

-900

-1000

10

8

6

× 10 4

10 9

4

8

7 6 5 4

2

× 10 4

3 2 1

0 0

Distance to the initial geographic point [m] - Longitude

Distance to the initial geographic point [m] - Latitude

Figure 15: Representation of the simulation of a single ray at 8x10 4 m from the initial geographic coordinates (53 ◦ N and 13 ◦ W) measured from longitudes and latitudes at a depth of 350 m.

To my co-advisors, Carlos Costa and Mário Gatta, I would also like to extend my thanks for all the knowledge passed in the area of oceanography and their valuable technical support on this work. Last but not least, I acknowledge CINAV the Portuguese Navy’s Research Center for the support in carrying out this work.

REFERENCES

[1] Daljit S Ahluwalia and Joseph B Keller. Wave propagation and underwater acoustics . Springer, 1977.

[2] Paul C Etter. Underwater acoustic modeling: principles, techniques and applications . CRC Press, 1995. [3] P. Raposo, M. Gatta, and M. Moreira. Ocean-acoustic raytracing propagation models based on fermat’s least time principle. Inter-Noise , 2019. [4] P. Raposo. Estudo da propagação acústica 3D em meios não homogéneos – Raytracing baseado no príncipio de Fermat . Escola Naval, 2019. [5] Uri M Ascher and Linda R Petzold. Computer methods for ordinary di ff erential equations and di ff erential-algebraic equations , volume 61. Siam, 1998. [6] A. J. Roberts. Di ff erential algebraic equations solvers. Technical report, The University of Adelaide, 1998. [7] Tiago CA Oliveira, Ying-Tsong Lin, and Michael B Porter. Underwater sound propagation modeling in a complex shallow water environment. Frontiers in Marine Science , page 1464, 2021. [8] J.A. Barahona Fernandes. Manual de Hidrografia . Instituto Hidrográfico, 1967. [9] Feistel, R TJ, McDougall, FJ Millero, DR Jackett, DG Wright, BA King, GM Marion, C Chen, P Spitzer, and S Seitz. The international thermodynamic equation of seawater 2010 (teos-10): Calculation and use of thermodynamic properties. Global Ship-based Repeat Hydrography Manual, IOCCP Report No , 14, 2009.