Prediction of Convergence Zone Parameters Based on Improved Ray Model

1. Introduction

The convergence zone represents a distinctive acoustic propagation phenomenon within the deep-sea environment. This phenomenon is characterized by a periodically focused acoustic energy structure, which facilitates the transmission of acoustic signals over extended distances with reduced attenuation, thereby significantly augmenting the operational efficacy of sonar detection systems [1]. Predicting the characteristic parameters of the deep-sea convergence zone necessitates an examination of its distance, width, and gain. Investigating these parameters enables underwater combat units to precisely understand the influence of convergence zone characteristics on equipment performance. This understanding provides a supplementary basis for decision-making in the formulation of detection maneuver strategies during deep-sea operations.

In contemporary research, the primary approach for predicting the characteristic parameters of convergence zones involves analyzing the propagation paths of acoustic rays through the development of acoustic field models [2,3,4,5,6]. Notably, the Bellhop ray acoustic field model, introduced by Porter and Bueker, serves as a prominent method for modeling acoustic fields in horizontally inhomogeneous environments [7]. This model has gained widespread application in recent years for examining acoustic propagation characteristics. For instance, Weishuai et al. utilized the Douglas-Peucker algorithm to extract oceanographic factors from a database and employed the BELLHOP ray tracing model to derive acoustic propagation characteristics, thereby offering a predictive tool for studying underwater acoustic propagation in oceanic front environments [8] . Similarly, Liu et al. utilized a three-dimensional ocean acoustic framework, integrating the MITgcm and BELLHOP ray model, to investigate the coupling effects of the Luzon cold vortex and tidal waves [9]. Additionally, Ma et al. applied the BELLHOP ray theory model to explore the acoustic field characteristics associated with a pair of cyclonic vortices and a typical anticyclonic vortex [10].

Moreover, there is still a significant discrepancy between predictions using traditional ray models and acoustic experimental observations [12,13]. Therefore, data-driven prediction methods have been widely applied. Yang et al. designed an intelligent convergence zone recognition algorithm based on convolutional neural networks [14]. By inputting sound speed profiles and acoustic propagation loss data, the algorithm automatically extracts the spatiotemporal characteristics of caustic zones and achieved synchronous prediction of convergence zone distance and width in the East China Sea experiment, significantly improving efficiency compared to traditional physical models. Xu et al. developed a convergence zone prediction model based on high-resolution ocean fronts and tested 24 machine learning algorithms using k-fold cross-validation [15]. Li et al. constructed a convergence zone parameter prediction model under the influence of mesoscale eddies by integrating numerical simulations and measured data in multiple dimensions [16]. Ibebuchi and Richman used an autoencoder neural network to encode nonlinear sea surface temperature patterns in the tropical Pacific to predict the Nino 3.4 index [17]. Xu et al.utilized physics-informed machine learning methods to identify and predict acoustic convergence zone characteristics of mesoscale eddies under limited data conditions [18]. However, data-driven prediction methods struggle to cover multi-spatiotemporal scale acoustic-ocean joint datasets and face dual challenges of data quality and model interpretability due to their opaque physical mechanisms.

The actual marine environment is highly spatially heterogeneous and dynamically complex. Sound speed profiles are not uniformly distributed but are influenced by factors such as temperature and salinity, exhibiting significant horizontal gradient characteristics [19]. Particularly in the upper ocean (mixed layer), due to the effects of solar radiation, seasonal temperature differences, and ocean current shearing, sound speed profiles often show strong inhomogeneity [20]. Neither traditional ray models nor data-driven prediction methods have fully considered the horizontal inhomogeneity of sound speed in the deep-sea mixed layer, resulting in cumulative offsets in ray trajectories and directly affecting the prediction accuracy of convergence zone distance and width.

In this article, a horizontal sound speed gradient term is introduced into the ray equation, and a fourth-order Runge-Kutta numerical solution algorithm is designed. This breakthrough overcomes the limitation of traditional ray models that cannot represent the influence of horizontal sound speed gradients, thereby improving the prediction accuracy of convergence zone ray trajectories. The main innovations are as follows.

• A method for constructing horizontally inhomogeneous sound speed profiles in the deep-sea mixed layer is proposed. By using EOF decomposition technology to extract the spatial modal characteristics of temperature and salinity fields and combining them with the Del Grosso formula, a three-dimensional continuous sound speed field is generated. This approach overcomes the limitation of traditional models that simplify mixed layer sound speed as vertically stratified, achieving dynamic coupled modeling of sound speed profiles with horizontal distance and vertical depth.
• An improved ray tracing algorithm based on non-uniformly distributed sound speed profiles is developed. By introducing the horizontal sound speed gradient term into the ray equation and designing a fourth-order Runge-Kutta numerical solution algorithm, the bottleneck of traditional ray models in representing the influence of horizontal sound speed gradients is overcome, thereby improving the prediction accuracy of convergence zone ray trajectories.
• The enhanced algorithm for forecasting convergence zone characteristic parameters was validated using environmental parameters from a representative marine area. The findings indicate that the refined ray model successfully accounted for the forward displacement of the convergence zone induced by the horizontal gradient of sound speed within the mixed layer, thereby enhancing the forecast’s accuracy.

The rest of the article is organized as follows. Section 2 is the analysis of the non-uniform sound speed profile in the deep-sea mixed layer. Section 3 is the forecast of convergence zone characteristic parameters based on an improved ray model. Section 4 is the simulation. Section 5 is the conclusion.

2. Analysis of the Non-uniform Sound Speed Profile in the Deep-sea Mixed Layer

This article constructs a non-uniform sound speed profile model for the mixed layer by measuring temperature, salinity, and pressure data and substituting them into empirical formulas [21,22]. Firstly, the driving factors of sound speed in the mixed layer are analyzed from the perspective of physical mechanisms. It is clarified that the horizontal gradients of the temperature and salinity fields are the main sources of spatial heterogeneity of sound speed, while the pressure field dominates the vertical stratification structure through the hydrostatic effect. Secondly, within the framework of the classical Del Grosso sound speed formula [23], the temperature and salinity parameters are decomposed into coupled functions of horizontal distance and vertical depth using gridded measured temperature, salinity, and pressure data, thereby establishing the non-uniform distribution relationship of sound speed with horizontal distance and vertical depth within the mixed layer of the sea area. Finally, the mixed layer depth is introduced as a boundary threshold to construct a piecewise continuous sound speed profile model. The steps are shown in Figure 1.

2.1. Analysis of Factors Affecting Sound Speed in the Mixed Layer

This section employs the empirical formulae for ocean sound velocity calculation [23,24,25], to analyze the impact mechanisms of the mixed layer sound velocity from the perspective of physical mechanisms, which are shown in equation (1), (2) and (3):

$$\begin{array}{cc}\hfill C& =1449.2+4.6T-0.055T^2+0.000029T^3\hfill \\ \hfill & +\left(1.34-0.01T\right)\left(S-35\right)+0.016Z.\hfill \end{array}$$

(1)

$$\left\{\begin{array}{cc}\hfill & V=1449.14+V_T+V_P+V_S+V_{STP},\hfill \\ \hfill & V_T=4.5721T-4.4532\times {10}^{-2}{T}^2\hfill \\ \hfill & -2.604\times {10}^{-4}{T}^3+7.9851\times {10}^{-6}{T}^4,\hfill \\ \hfill & V_P=1.60272\times {10}^{-1}P+1.0268\times {10}^{-5}{P}^2\hfill \\ \hfill & +3.5216\times {10}^{-9}{P}^3-3.3603\times {10}^{-12}{P}^4,\hfill \\ \hfill & V_S=1.39799(S-35)+1.69202\times {10}^{-3}{(S-35)}^2,\hfill \\ \hfill & V_{STP}=(S-35)(-1.1244\times {10}^{-2}T\hfill \\ \hfill & +7.7711\times {10}^{-7}{T}^2+7.7016\times {10}^{-5}P\hfill \\ \hfill & -1.2943\times {10}^{-7}{P}^2+3.1580\times {10}^{-8}PT\hfill \\ \hfill & +1.5790\times {10}^{-9}P{T}^2)+\mathrm{P}(-1.8607\times {10}^{-4}T\hfill \\ \hfill & +7.4812\times {10}^{-6}T+4.5283\times {10}^{-8}{T}^3)\hfill \\ \hfill & +P^2(-2.5294\times {10}^{-7}T+1.8563\times {10}^{-9}{T}^2)\hfill \\ \hfill & +P^3(-1.9646\times {10}^{-10}T).\hfill \end{array}\right.$$

(2)

$$\begin{array}{cc}\hfill & C=(1.389-1.262\times {10}^{-2}t+7.164\times {10}^{-5}{t}^2\hfill \\ \hfill & +2.006\times {10}^{-6}{t}^3-3.21\times {10}^{-6}{t}^4)S\left(s_0\right)\hfill \\ \hfill & +(-1.922\times {10}^{-2}-4.42\times {10}^{-5}t)S{\left(s_0\right)}^{3/2}\hfill \\ \hfill & +1.727\times {10}^{-3}S{\left(s_0\right)}^2,\hfill \end{array}$$

(3)

where C and V represent the sound speed, T is the temperature of sea, S is the salinity of sea, Z is depth and P is pressure. According to the theory of the wave equation, sound speed refers to the phase velocity of a plane wave, which is a longitudinal wave and is related to density and compressibility. It can be expressed as:

$$C={\displaystyle \frac{1}{\sqrt{\rho \beta }}},$$

(4)

where $\rho$ and $\beta$ represent The density and adiabatic compressibility of seawater.

2.1.1. The Effection of Temperature on Sound Speed

Based on the functional relationship between the compressibility of seawater and the temperature field, an empirical relationship model was established using marine observational data [26]:

$$\beta =481\times {10}^{-13}(1-0.00707T).$$

(5)

Let ${\beta }_0=481\times {10}^{-13}$, $S_0=0.00707$, then we have

$$\beta ={\beta }_0(1-S_{0}T).$$

(6)

Within the typical fluctuation range of conventional marine environmental parameters, the variation in density is minimal. By adopting the assumption of approximate constancy, we have

$$C_e={\displaystyle \frac{1}{\sqrt{\rho {\beta }_0(1-S_{0}T)}}}={\displaystyle \frac{C_0}{\sqrt{(1-S_{0}T)}}},$$

(7)

where $C_0$ represents the sound speed at a temperature of $0^{\circ }$C.

By applying the Taylor series expansion method to linearize the above equation and retaining the first and second-order terms to construct a simplified model, then we have

$$C_e=C_0\left(1+{\displaystyle \frac{1}{2}}{S}_{0}T\right).$$

(8)

The expression for the change in sound speed derived through this model is

$$\Delta \mathrm{C}=C_e-C_0={\displaystyle \frac{1}{2}}{S}_{0}T=C_0\times 0.00354T.$$

(9)

Equation(9) indicates that when the temperature changes by $1^{\circ }$C, the relative increase in sound speed is approximately 0.354% of the original value.

2.1.2. The Effection of Salinity on Sound Speed

According to the empirical formula between seawater density and salinity [26], we have

$$\rho ={\rho }_0(1+0.0008S)={\rho }_0(1+{\alpha }_{p}S),$$

(10)

where S represents salinity. The relationship between the compressibility coefficient of seawater and salinity can be expressed as

$$\beta ={\beta }_0(1-0.0024S)={\beta }_0(1-{\alpha }_{k}S).$$

(11)

From equation (10) and equation (11), it can be seen that When the salinity increases by 1, the compressibility decreases by 0.24%, and the sound speed increases; when the density increases by 0.08%, the sound speed decreases. Substituting equation (10) and equation (11) into equation (4), we have

$$C={\displaystyle \frac{C_0}{\sqrt{(1+{\alpha }_{p}S)(1-{\alpha }_{k}S)}}}\cong C_0+{\displaystyle \frac{1}{2}}{C}_0({\alpha }_k-{\alpha }_p)S.$$

(12)

Substituting ${\alpha }_k=0.0024$ and ${\alpha }_p=0.0008$ into equation (12) , we have

$$\Delta C=C-C_0=C_0\times 0.000825S.$$

(13)

From equation (13), it can be concluded that when salinity increases by 1, the relative increase in sound speed is approximately 0.0825% of the original value.

2.1.3. The Effect of Pressure on the Speed of Sound

The change in sound speed due to pressure mainly depends on the variation of the compressibility coefficient $\beta$. For seawater, the greater the pressure, the more difficult it is to compress, and the smaller the compressibility coefficient becomes, which ultimately leads to an increase in sound speed. Therefore, pressure is positively correlated with sound speed. The following is the expression for the relationship between the static pressure of seawater and the compressibility coefficient:

$$\beta ={\beta }_0(1-0.00044P)={\beta }_0(1-{\mu }_{k}P),$$

(14)

where the pressure is measured in units of one standard atmosphere (101325 Pa). Substituting equation (14) into equation (4), the change in sound speed is finally obtained as

$$\Delta C={\displaystyle \frac{C_0}{2}}{\mu }_{k}P=0.00022C_{0}P.$$

(15)

As can be seen from equation (15), when the pressure increases by 1 Pa, the relative increase rate of sound speed is approximately 0.022% of its original value.

2.2. Construction of Non-uniform Sound Speed Profiles

This section elaborates on the data-driven method for constructing non-uniform sound speed profiles, including data acquisition and processing, modeling of the horizontal-vertical distribution of temperature and salinity, and Del Grosso sound speed coupling analysis.

2.2.1. Data Acquisition and Preprocessing

The data source used in this article is the ocean space temperature, salinity, and pressure data, which is derived from the global ocean Argo gridded data set (BOA_Argo) provided by the Argo Real-Time Data Center [27]. The data applied include longitude, latitude, pressure, salinity, temperature, and mixed layer depth. Taking temperature as an example, its surface distribution is shown in Figure 2(a). After obtaining the measured data, the data are preprocessed by regional clipping according to the research sea area. The index corresponding to longitude $i_0$ is obtained through the following minimization formula, and the same applies to latitude j. The final longitude and latitude index ranges are determined as $i\in \left[i_{min},i_{max}\right]$ and $j\in \left[j_{min},j_{max}\right]$, and the clipped sea area is shown in Figure 2(b).

2.2.2. Modeling the horizontal and vertical distribution of temperature and salinity

To accurately characterize the continuous spatial distribution of temperature, salinity, and pressure, this article employs a two-dimensional Empirical Orthogonal Function (EOF) decomposition method, combined with measured data and function interpolation techniques, to construct a temperature field parameterization model [28]. The original temperature and salinity data have a horizontal resolution of $1^{\circ }\times 1^{\circ }$, and are vertically divided at 10m intervals within the mixed layer. The data are preprocessed to eliminate the dimensional differences in the vertical stratification and to highlight the spatial characteristics of the temperature and salinity fields. Each depth layer of the temperature and salinity fields is normalized.

Based on the BOA_Argo data, the horizontal-vertical distribution model of temperature and salinity established using EOF and the comparison with the original data are shown in Figure 3.

2.2.3. Del Grosso acoustic velocity coupling

Considering the horizontal inhomogeneity of the sound speed profile in the deep-sea mixed layer, the sound speed profile is described in segments. The entire sound speed profile can be represented by the following equation:

$$C(r,z)=\left\{\begin{array}{cc}C(r,z),& 0<z<z_M,\\ C\left(z\right),& z_M<z<H,\end{array}\right.$$

(16)

where, r represents the horizontal distance, z represents the vertical depth, $z_M$ represents the depth of the mixed layer, and H represents the sea depth. When $0<z<z_M$, the sound speed profile will fluctuate with changes in horizontal distance and vertical depth, exhibiting non-uniform distribution characteristics. When $z_M<z<H$, the sound speed is located below the mixed layer and varies only with depth. The schematic diagram of the sound speed profile is shown in Figure 4, with the mixed layer above the stratification line. The solid line represents the simplified equivalent sound speed profile that ignores the non-uniform changes in the marine environment, and the shaded area represents the fluctuation range caused by the non-uniform distribution of sound speed.

The Del Grosso algorithm, as shown in equation (1). To ensure the continuity of the sound speed profile at $z_M$, the function C in equation (16) must satisfies $C\in C^1$.

The temperature and salinity at the base of the mixed layer must smoothly transition to the deep-layer constant values, which is specifically realized by introducing the transition function $\eta \left(z\right)$:

$$\begin{array}{cc}\hfill T^{\prime }(r,z)=& \eta \left(z\right)T(r,z)+(1-\eta (z\left)\right)T,\hfill \\ \hfill & {\mathrm{z}}_M-\Delta z\le z\le {\mathrm{z}}_M+\Delta z.\hfill \end{array}$$

(17)

where $\eta \left(z\right)$ is a cubic polynomial weighting function:

$$\eta \left(z\right)=\left\{\begin{array}{cc}& 1,z\le {\mathrm{z}}_M-\Delta z,\\ & 2{\left({\displaystyle \frac{{\mathrm{z}}_M+\Delta z-z}{2\Delta z}}\right)}^3-3{\left({\displaystyle \frac{{\mathrm{z}}_M+\Delta z-z}{2\Delta z}}\right)}^2+1,\\ & {\mathrm{z}}_M-\Delta z<z\le {\mathrm{z}}_M+\Delta z,\\ & 0,z>{\mathrm{z}}_M+\Delta z.\end{array}\right.$$

(18)

To avoid a significant impact on the sound speed, the thickness of the transition layer, $\Delta z$, is set to only 2-5m, to ensure the continuity of the first-order derivative of the sound speed. In this case, the temperature gradient at the base of the mixed layer, ${\displaystyle \frac{\partial T}{\partial z}}$, and the salinity gradient, ${\displaystyle \frac{\partial S}{\partial z}}$, will smoothly decay to zero, avoiding the occurrence of a sawtooth discontinuity in the sound speed at the boundary point.

3. Forecast of Convergence Zone Characteristic Parameters Based on an Improved Ray Model

3.1. Ray model improvement

Starting from the classical ray equation, this section derives the improved ray tracing equation for non-uniform media, designs a numerical solution algorithm, and clarifies the implementation process of the improved model in the prediction of parameters in the convergence zone. According to the eikonal equation, the ray direction is defined as the normal direction of the wavefront, where the wavefront is the set of points with the same acoustic pressure vibration phase. The normal direction at any point on the wavefront is the direction of the phase gradient at that point, representing the direction of sound wave propagation, that is:

$${\displaystyle \frac{d\mathbf{x}}{ds}}={\displaystyle \frac{\nabla \tau }{|\nabla \tau |}}=C\nabla \tau .$$

(19)

Combining the eikonal equation $|\nabla \tau |=1/C$, we have

$${\displaystyle \frac{d\mathbf{x}}{ds}}=C\nabla \tau \Rightarrow \nabla \tau ={\displaystyle \frac{1}{C}}{\displaystyle \frac{d\mathbf{x}}{ds}}.$$

(20)

Differentiate with respect to s, using the chain rule, we have

$${\displaystyle \frac{d}{ds}}\left({\displaystyle \frac{1}{C}}{\displaystyle \frac{d\mathbf{x}}{ds}}\right)=\nabla \left({\displaystyle \frac{1}{C}}\right),$$

(21)

after expanding equation(20), we obtain that

$${\displaystyle \frac{d^2\mathbf{x}}{d{s}^2}}={\displaystyle \frac{1}{C}}\nabla C-{\displaystyle \frac{2}{C^2}}\left({\displaystyle \frac{d\mathbf{x}}{ds}}·\nabla C\right){\displaystyle \frac{d\mathbf{x}}{ds}}.$$

(22)

In the two-dimensional vertical plane (horizontal x, vertical z), let the angle between the sound ray and the horizontal plane be $\theta$, then we have

$${\displaystyle \frac{dx}{ds}}=sin\theta {\displaystyle \frac{dz}{ds}}=cos\theta .$$

(23)

Differentiate with respect to $\theta$, using the geometric relationship, we have

$${\displaystyle \frac{d\theta }{ds}}=-{\displaystyle \frac{1}{C}}{\displaystyle \frac{dC}{dz}}cos\theta .$$

(24)

Finally, we obtain the differential equation system describing the sound ray trajectory in the traditional ray model by equation (19) – (24):

$$\left\{\begin{array}{cc}& {\displaystyle \frac{dr}{ds}}=sin\theta ,\hfill \\ & {\displaystyle \frac{dz}{ds}}=cos\theta ,\hfill \\ & {\displaystyle \frac{d\theta }{ds}}=-{\displaystyle \frac{1}{C\left(z\right)}}{\displaystyle \frac{dc}{dz}}cos\theta .\hfill \end{array}\right.$$

(25)

In the equation, s is the arc length of the sound ray, $\theta$ is the angle between the sound ray and the horizontal plane, and $C\left(z\right)$ is the vertically stratified sound speed. This model does not consider the effect of the horizontal gradient term of the sound speed ${\displaystyle \frac{\partial C}{\partial r}}$.

Considering the non-uniform sound speed profile $C(r,z)$ that depends on both the horizontal position r and depth z, the ray model is improved starting from the eikonal equation, and the model must include the horizontal gradient term ${\displaystyle \frac{\partial C}{\partial r}}$. The propagation of sound waves can be regarded as particles moving along ray trajectories. Therefore, starting from the Hamiltonian mechanics framework, the ray tracing problem is transformed into a particle trajectory problem in classical mechanics. The position coordinates are defined as $mathbfq=(r,z)$ and the momentum as $\mathbf{p}=\nabla \tau$, with the Hamiltonian [29]

$$H={\displaystyle \frac{1}{2}}\left(p_r^2+p_z^2-{\displaystyle \frac{1}{C^2(r,z)}}\right).$$

(26)

The equations of motion in Hamiltonian mechanics are:

$${\displaystyle \frac{d\mathbf{q}}{ds}}={\displaystyle \frac{\partial H}{\partial \mathbf{p}}}{\displaystyle \frac{d\mathbf{p}}{ds}}=-{\displaystyle \frac{\partial H}{\partial \mathbf{q}}}.$$

(27)

After expanding equation (27) and using equation (26), the canonical equation that the ray trajectory satisfies:

$$\left\{\begin{array}{cc}& {\displaystyle \frac{dr}{ds}}={\displaystyle \frac{\partial H}{\partial p_r}}=p_r,\hfill \\ & {\displaystyle \frac{dz}{ds}}={\displaystyle \frac{\partial H}{\partial p_z}}=p_z,\hfill \\ & {\displaystyle \frac{d{p}_r}{ds}}=-{\displaystyle \frac{\partial H}{\partial r}}={\displaystyle \frac{1}{C^3}}{\displaystyle \frac{\partial C}{\partial r}},\hfill \\ & {\displaystyle \frac{d{p}_z}{ds}}=-{\displaystyle \frac{\partial H}{\partial z}}={\displaystyle \frac{1}{C^3}}{\displaystyle \frac{\partial C}{\partial z}}.\hfill \end{array}\right.$$

(28)

From the eikonal equation $|\nabla \tau |={\displaystyle \frac{1}{C}}$, the physical meaning of the momentum $\mathbf{p}=\nabla \tau$ is the direction of the wave vector. Combining this with the sound ray direction angle $\theta$, we can write $p_r={\displaystyle \frac{sin\theta }{C}}$ and $p_z={\displaystyle \frac{cos\theta }{C}}$. Through variable substitution, we obtain the improved ray equation system:

$$\left\{\begin{array}{cc}& {\displaystyle \frac{dr}{ds}}={\displaystyle \frac{sin\theta }{C(r,z)}},\hfill \\ & {\displaystyle \frac{dz}{ds}}={\displaystyle \frac{cos\theta }{C(r,z)}},\hfill \\ & {\displaystyle \frac{d\theta }{ds}}=-{\displaystyle \frac{1}{C(r,z)}}\left({\displaystyle \frac{\partial C}{\partial z}}cos\theta -{\displaystyle \frac{\partial C}{\partial r}}sin\theta \right).\hfill \end{array}\right.$$

(29)

Compared with the classical equation (25), the added term $\left({\displaystyle \frac{\partial C}{\partial r}}\right)sin\theta$ in equation (29) quantifies the contribution of the horizontal gradient to the curvature of the sound ray. When ${\displaystyle \frac{\partial C}{\partial r}}\ne 0$, the sound ray will deflect towards the direction of increasing sound speed, resulting in a lateral displacement of the propagation path, as shown in Figure 5.

The complexity of the improved ray-tracing algorithm makes the equations analytically intractable, necessitating the use of numerical methods. This article employes the fourth-order Runge-Kutta method (RK4) to numerically integrate the improved ray equations. The pseudocode for its algorithm implementation is shown as Algorithm 1.

3.2. Forecast of Convergence Zone Characteristic Parameters

The construction of sound rays in the acoustic field is achieved using the Bellhop acoustic toolbox in Matlab, which further enables the calculation of the distance, width, and gain of the first convergence zone. The subsequent rows contain the actual sound speed values, in meters per second, with each row corresponding to a fixed depth. The actual depths used are taken from the depth data below the fifth line in the environmental file.

When analyzing specific sea areas, bathymetric data can be obtained from the General Bathymetric Chart of the Oceans (GEBCO) database [30]. The plotssp2d command is used to draw the horizontal inhomogeneous distribution of sound speed profiles in the mixed layer. The resulting sound speed profile curves associated with horizontal distance and vertical depth, as well as the color-mapped sound speed profile images, are shown in Figure 6(a) and Figure 6(b), respectively. The plotray command is used to track sound rays and generate the sound field model of the convergence zone, as shown in Figure 6(c). It can be seen that the sound rays in the convergence zone, considering the horizontally inhomogeneous sound speed profile, exhibit significant fluctuations compared to those under the typical Munk profile during propagation. The overall shape of the sound rays is more complex and clearly conforms more closely to the actual sound signal propagation patterns. Based on this, the prediction of the distance and width of the first convergence zone is consistent with the method under the typical Munk profile and can be easily calculated by measuring the upper bending point of the sound ray with a $0^{\circ }$ launch angle and the critical depth bending sound ray. The prediction of the gain in the first convergence zone is also calculated by constructing a transmission loss model, as shown in Fig. Figure 6(d).

4. Simulation

In this section, the acoustic propagation characteristic experiment conducted in the South China Sea in 2018 is selected as the benchmark case data.

4.1. Characteristic Features of Typical Marine Environment and Data Preparation

The typical sea area selected for the simulation experiment of the convergence zone characteristic parameter forecasting algorithm is covered by the sea area of the 2018 spring South China Sea acoustic propagation experiment. The selected area is located at $17.5^{\circ }$N latitude and $114.5^{\circ }$E-$115.5^{\circ }$E longitude in the South China Sea, with a depth of approximately 3500 meters and a horizontal distance of about 110 kilometers. The environmental characteristic data include temperature, salinity, pressure, and mixed layer depth. The temperature and salinity profiles are shown in Figure 7.

Algorithm 1 Improved sound ray tracing algorithm for the ray model.

In the South China Sea experiment, the source depth was 80 m and the source frequency was 300 Hz. The convergence zone characteristic parameters obtained using the first-order internal wave normal mode basis function model are shown in Table 1, which will be used to verify the accuracy of the Bellhop ray model and the improved ray model in predicting convergence zone characteristics.

4.2. Simulation of Convergence Zone Characteristic Parameter Forecasting and Comparative Analysis of Results

The temperature, salinity, and depth data of the typical sea area are substituted into the empirical sound speed formula of Del Grosso. The sound speed profile is then written into the environmental file to generate a simplified equivalent sound speed profile for the entire sea area under the traditional ray model, as shown in Figure 8.

Under the Bellhop ray model, the simple equivalent sound speed profile is relatively smooth. The mixed layer does not take into account the horizontal distance. The sound channel axis is located at a depth of about 800 meters. The minimum sound speed is 1496.74 m/s, the sea surface sound speed is 1515.9 m/s, and the sea bottom sound speed is 1518.67 m/s. Under this sound speed profile, there is a convergence zone effect.

Combining the temperature, salinity, and depth data with the latitude and longitude data of typical sea areas, a horizontal-vertical distribution model of the temperature-salinity field was constructed using the two-dimensional empirical orthogonal functions. By coupling the Del Grosso empirical formula for sound speed, the horizontally non-uniform sound speed profiles in the mixed layer were obtained, as shown in Figure 9 and Figure 10.

Figure 9 shows that the sound speed profile changes with horizontal distance due to the mixed layer’s temperature and salinity inhomogeneity. At the sea surface, sound speed is 1515.86 m/s at $114.5^{\circ }$E, 1517.65 m/s at ${115}^{\circ }$E, and 1519.67 m/s at $115.5^{\circ }$E, while it is nearly constant below the mixed layer. This is because the mixed layer’s temperature and salinity increase with distance, but remain constant below it. The color mapping of sound speed also shows that the mixed layer’s sound speed rises with distance. Although the sound speed considering the mixed layer’s non-uniformity differs from the simple equivalent sound speed, their overall trends are similar, with the sound channel axis depth around 800 m. This confirms the effectiveness of the improved sound speed profile.

Based on the construction of the sound speed profile, the acoustic propagation paths under the Bellhop ray model and the improved ray model are simulated. According to the experimental scenario, the source depth is set to 80 m, the source frequency is set to 300 Hz, and the emission angle range is set to $\pm {10}^{\circ }$. The sound rays of the sound field drawn by the Bellhop ray model and the improved ray model are shown in Figure 11.

The improved ray model sound field shifts significantly forward compared to the Bellhop ray model sound field. This is due to the consideration of the horizontally inhomogeneous sound speed distribution, which changes the turning points of the sound rays. Specifically, the turning point depth of rays with smaller emission angles becomes shallower, while that of rays with larger emission angles becomes deeper. The upper turning points of the $0^{\circ }$ emission angle ray and the critical depth bending ray are shown in Figure 12.

It can be seen that, according to the Bellhop ray model, the distance to the first convergence zone is 63.6742 km, and the width of the first convergence zone is 20.1763 km; according to the improved ray model, the distance to the first convergence zone is 60.7144 km, and the width of the first convergence zone is 20.7463 km.

The first convergence zone gain is obtained by plotting the propagation loss model. The magnitude of the propagation loss varies with the receiver location. Therefore, according to the experimental setup, the receiver is placed at a depth of 200 meters within the convergence zone distance to measure the magnitude of the convergence zone gain. The propagation loss models of the Bellhop ray model and the improved ray model are shown in Figure 13.

The Bellhop ray model shows a propagation loss of about 93 dB outside the convergence zone and a minimum loss of 77.8827 dB inside it, giving a convergence zone gain of around 15.12 dB. The improved ray model has a propagation loss of about 90 dB outside the convergence zone and a minimum loss of 77.334 dB inside, resulting in a gain of about 12.67 dB. Thus, the convergence zone gain from the Bellhop model is slightly higher than that from the improved model. Compare the prediction results of the characteristics parameters of the first convergence zone by the Bellhop ray model and the improved ray model with the experimental results, as shown in Table 2.

The comparison shows that the improved ray model has enhanced the accuracy of the first convergence zone distance prediction by 2.26%, the first convergence zone width prediction by 2.66%, and the first convergence zone gain prediction by approximately 5.85% compared with the Bellhop ray model, which demonstrates the effectiveness of the improvements made to the ray model in this article.

5. Conclusion

This article proposes an improved ray model framework that constructs a sound speed profile model with horizontal gradient corrections by integrating real measured ocean environmental data, and verifies its effectiveness. The model significantly enhances the prediction accuracy of deep-sea convergence zone characteristic parameters, effectively addressing the forecast deviation caused by traditional models ignoring the horizontal inhomogeneity of sound speed in the mixed layer. Future research will employ machinoee learning methods to further improve the prediction capability of convergence zone acoustic field effects and apply it to practical mission scenarios, thereby promoting the development of underwater detection and communication technologies.