## Abstract

Light influenced by the turbulent ocean can be fully characterized with the help of the power spectrum of the water’s refractive index fluctuations, resulting from the combined effect of two scalars, temperature and salinity concentration advected by the velocity field. The Nikishovs’ model [ Fluid Mech. Res. **27**, 8298 (2000)] frequently used in the analysis of light evolution through the turbulent ocean channels is the linear combination of the temperature spectrum, the salinity spectrum and their co-spectrum, each being described by an approximate expression developed by Hill [ J. Fluid Mech. **88**, 541562 (1978)] in the first of his four suggested models. The fourth of the Hill’s models provides much more precise power spectrum than the first one expressed via a non-linear differential equation that does not have a closed-form solution. We develop an accurate analytic approximation to the fourth Hill’s model valid for Prandtl/Schmidt numbers in the interval [3, 3000] and use it for the development of a more precise oceanic power spectrum. To illustrate the advantage of our model, we include numerical examples relating to the spherical wave scintillation index evolving in the underwater turbulent channels with different average temperatures, and, hence, different Prandtl numbers for temperature and different Schmidt numbers for salinity. Since our model is valid for a large range of Prandtl number (or/and Schmidt number), it can be readily adjusted to oceanic waters with seasonal or extreme average temperature and/or salinity or any other turbulent fluid with one or several advected quantities.

© 2019 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

## 1. Introduction

The growing interest in the development of optical underwater imaging systems [1], laser communications [2–4] and remote sensing schemes [5] stipulates the request for the accurate predictions of the oceanic turbulence effects on the light wave evolution. Such effects can be characterized with the help of the parabolic equation, the Rytov approximation, the extended Huygens-Fresnel integral and other classic methods as long as the refractive-index power spectrum of turbulent fluctuations is established [6–9]. As applied to oceanic propagation, the theoretical predictions for the major light statistics of the optical waves have been established (see recent review [10]) on the basis of the widely known power spectrum developed by Nikishovs [11]. In particular, the spectral density and the beam spread for coherent and random beams have been analyzed [12], the spectral shifts in random beams were revealed [13], the polarimetric changes of the random electromagnetic beams have been discussed [14], the loss of coherence of initially coherent beams has been addressed [15], the scintillation analysis was provided in Refs. [16] and [17] and the structure functions of various waves were derived in Refs. [18] and [19].

In the situations when there is more than one scalar field advected by turbulence, the total refractive-index power spectrum can be expressed as a linear combination of the individual spectra and their co-spectra. For instance, in the case of the oceanic turbulence with advected temperature and salinity (sodium chloride) concentration, the total spectrum is the linear combination of three spectra: temperature spectrum, salinity spectrum, and their co-spectrum [11] (see also [20]). An important dimensionless quantity determining the turbulent spectrum profile is the Prandtl number, *Pr*_{T}, defined as a ratio of viscous diffusion rate (kinematic viscosity) [m^{2}s] to the thermal diffusion rate (thermal diffusivity) [m^{2}s]. Another important dimensionless quantity is the Schmidt number, *Pr*_{S}, defined as the ratio of kinematic viscosity and mass diffusivity. Since the Prandtl numbers for temperature and the Schmidt numbers for salinity in the water are sufficiently large (greater than one), their power spectra have three regimes: *inertial-convective* regime with the Kolmogorov power law Φ* _{n}*(

*κ*) ∼

*κ*

^{−11/3}at the low wave numbers,

*viscous-convective*regime with the power law Φ

*(*

_{n}*κ*) ∼

*κ*

^{−3}for higher wave numbers and the

*viscous-diffusive*regime at even higher wave-numbers where the spectra decrease very rapidly driven by the diffusion of the advected quantity. Due to two different power laws at the high wave numbers each spectrum forms a characteristic “bump” [21–23].

To achieve a more precise match of the spectra at the boundaries between different regimes, one of the four models proposed by Hill [24] may be applied. Among them, the Hill’s model 1 (H1) that provides the simple analytic approximation but lacks precision in accounting for the bump, and the Hill’s model 4 (H4) that is more precise but involves a non-linear differential equation without known analytic solution [25–29].

The widely used Nikishovs’ oceanic turbulence spectrum [11] is developed as a linearized polynomial involving temperature spectrum, salinity spectrum and their co-spectrum, each being based on H1 (see Appendix I for more details). Hence, it has a simple analytic form but lacks precision. It can, in principle, be applied to various Prandtl/Schmidt numbers but, as we will illustrate, substantially deviates from H4, especially for Prandtl/Schmidt numbers much larger than one. Also, a recent model for the oceanic power spectrum [23,30] was obtained by fitting H4 at specific *Pr* values: *Pr* = 7 for temperature spectrum, *Pr* = 700 for salinity spectrum and *Pr* = 13.86 for their co-spectrum. It was revealed that the oceanic power spectrum of [23] has a substantially higher bump at high spatial frequencies as compared with that for the Nikishov’s spectrum, resulting in a greater scintillation index in light waves interacting with such a medium.

Moreover, within the last two decades, several mathematically tractable models based on H4 have also been introduced for atmospheric refractive-index spectrum [31–35], by fitting *Pr* = 0.72 and *Pr* = 0.68 for temperature and humidity Prandtl numbers, respectively. For such low Prandtl numbers in the atmospheric turbulence the viscous-convective regime is not that prevalent but nevertheless the high-frequency bump can still be formed. Among the atmospheric power spectrum models are the Frehlich’s model being a fifth-order Laguerre expansion [32], as well as the Grayshan’s and the Andrews’ models both being products of polynomials and Gaussian functions [33, 34]. All these power spectra provide convenience for theoretical calculations, suggesting that a good analytical fits do not have to be mathematically involved.

The aim of this paper is to employ the H4 model for developing a power spectrum that can be readily adapted to the wide range of Prandtl number and Schmidt number but, at the same rate, remain numerically accurate and mathematically tractable. Such new spectrum can be applied for classic homogeneous oceanic turbulence in a variety of thermodynamic states. Moreover, it can also be employed for optical characterization of other turbulent media, with one or two advected quantities, for example, fresh turbulent water contaminated by a chemical. The ratio of kinematic viscosity to diffusivity is orders of magnitude larger for temperature and salinity fluctuations in water than for temperature and humidity fluctuations in air [22], resulting in much larger Prandtl numbers and Schmidt number (averaging at *Pr*_{T} = 7 for temperature and *Pr*_{S} = 700 for salinity [10,11,21,23,30]. More importantly, these numbers can exhibit variation. For example, at 1 Bar pressure (water surface), at fixed salt concentration at 34.9ppt, as the temperature increases from 0°C to 30°C, *Pr*_{T} decreases from 13.349 to 5.596 and *Pr*_{S} decreases from 2393.2 to 456.1 (more details in Appendix II). On the other hand certain oil-like substances can lead to Prandtl numbers up to 10^{5}.

The paper is organized as follows. In section 2.1 we introduce an approximate fit to H4 for *Pr* ∈ [3, 3000] for handling situations in which a single scalar is advected by turbulence. Section 2.2 provides the validity estimation of the fit by means of the dissipation constraint. Section 3.1 introduces an H4-based power spectrum of oceanic turbulence by applying the fit of Section 2.1 to three spectra: temperature-based, salinity-based and the their co-spectrum. In Section 3.2 we apply the new spectrum model for analytic calculation of the scintillation index of the spherical wave for a variety of Prandtl numbers and Schmidt numbers and illustrate its agreement with numerical calculation based on the H4 model as well as its discrepancy with the result based on the Nikishovs’ spectrum. In Section 4 all the results are briefly summarized.

## 2. Analytic approximation of the Hill’s model 4 for one advected quantity

#### 2.1. The polynomial-Gaussian fit of the Hill’s model 4

Let us begin by developing the three-dimensional (3D) scalar power spectrum of the refractive index fluctuations produced by a single quantity advected by a fluid, as an analytical fit to the H4 model. According to [29] the 3D scalar spectrum Φ* _{n}*(

*κ*) is related to a universal dimensionless function

*g*(

*κη*) by expression

*β*is the Obukhov-Corrsin constant;

*ε*is the dissipation rate of the kinetic energy [m

^{2}s

^{−3}];

*η*is the Kolmogorov microscale [m];

*χ*is the ensemble-averaged variance dissipation rate. The units of

*χ*depend on the advected quantity, for instance, for temperature, salinity and refractive index they are K

^{2}s

^{−1}, g

^{2}s

^{−1}and s

^{−1}, respectively.

According to model H4, *g*(*κη*) is the solution of the nonlinear ordinary differential equation [24,29]

*z*=

*κηa*

^{−1},

*c*=

*βa*

^{4/3}

*Pr*

^{−1}. Generally, as recommended by Hill [21] the following fixed values of parameters

*a*= 0.072,

*b*= 1.9, and

*β*= 0.72 [36] will be used throughout the paper. The numerical solution

*g*(

*κη*) of Eq. (2) exhibits an obvious “bump” with height

*g*

_{max}gradually increasing with growing values of

*Pr*, for

*Pr*> 0.2 [29].

Considering the mathematical convenience given by Gaussian function in related calculation [6,37], we keep the form — a product of a polynomials and a Gaussian function [38] — from Andrews’ model [33, 39]. In addition, to include the significant effect of *Pr* on *g*(*κη*) [29], we separate (*κη*) and *c*, and give them different values of power [40]:

*h*

_{2},

*h*

_{5}and

*h*

_{7}are limited to positive numbers and then fitted the numerical solution

*g*(

*κη*) of Eq. (2) with Eq. (4) by iterative least squares method. This resulted in the following expression for

*g*(

*κη*):

*Pr*∈ [3, 3000]. On substituting from Eq. (5) into Eq. (1), we find that the fitted power spectrum based on H4 takes form

To examine the accuracy of the fitted solution for *g*(*κη*), we compare it with analytic formula H1 (see Appendix I) and H4 (calculated numerically) for several values of the Prandtl number. Figure 1 includes comparison among the H1 model, H4 model, our analytic fit and model in [23] at *Pr* = 7 and *Pr* = 700. At these values our fit is in a very good agreement with the H4 model and the model in [23]. Figures 2(a)–2(d) illustrate the ability of our analytic fit to match H4 model within the wide range [3,3000] of the Prandtl numbers. Figure 2(a) produced for *Pr* = 3 implies that the fitted function *g*(*κη*) calculated from Eq. (5) and that obtained from H1 model agree with that based on H4 model reasonably well while Eq. (5) provides with a better bump-shape reconstruction. Further, Figs. 2(b)–2(d), obtained for *Pr* = 30, 300 and 3000, respectively, illustrate that *g*(*κη*) in Eq. (5) is very consistent with that based on H4, but it is not the case for *g*(*κη*) based on H1, the latter leading to drastic underestimation of the bump’s strength. Hence, if one uses H4 as a benchmark, the proposed *g*(*κη*) in Eq. (5) and, hence, the Φ(*κ*) in Eq. (6) show advantages in describing the spectral bump for all *Pr* ∈ [3, 3000], and especially so for *Pr* ⩾ 30.

We note that *g*(*κη*) in Eq. (5) has the same functional form as those in the Andrews’ [33] and in the Grayshan’s [34] power spectra being the product of a polynomial and a Gaussian function. On the other hand, it also relies on *κ* and *Pr* independently, making the bump height *g*_{max} increase with growing *Pr* values, which is of importance in achieving the precision in the spectrum model valid for the wide range of *Pr*.

Thus, most existing H4-based models such as Yi’s model [23] and Andrews’ model [33], use fixed Prandtl number and/or fixed Schmidt number. These models are very accurate. For analysis of power spectra with *Pr*_{T} = 7 and Schmidt number *Pr*_{S} = 700, model in Ref. [23] can be used. However, as we have clarified, the diversity of oceanic environment leads to drastic deviation from values *Pr*_{T} = 7 and *Pr*_{S} = 700. In such cases model in Eq. (6) covering range [3, 3000] becomes very convenient.

#### 2.2. Dissipation constraint analysis

Although a good numerical fit for the power spectrum was obtained in the previous section for a wide range of *Pr*, it appears of importance to estimate whether the new model agrees with the basic laws of fluid mechanics. In this section we will verify to which degree the dissipation constraint [32,35] being the consequence of the scalar transport equation is satisfied by Eq. (5). The dissipation constraint can be expressed via integral [32]:

*g*(

*x*) agrees with the dissipation constraint.

The values of X calculated on substituting Eq. (5) into Eq. (7) are shown in Table 1 and Fig. 3 for different values of the Prandtl number. The table entries imply that Eq. (5) underpredicts X by 8.3% and 22.2% when *Pr* = 3 and 3000, respectively. This discrepancy is intensified with the increasing Prandtl number which, in its turn, corresponds to the increasing bump height. This indicates that a more precise fitting of the spectral bump occurs at the expense of deviation from the dissipation constraint.

For comparison, X corresponding to other models were calculated in [35], and are shown in Table 2. These models fit H4 at *Pr* = 0.72. Apart from the Frehlich model, the values of X for the rest of the models do exhibit certain deviation from value 1. On comparing entries of Tables 1 and 2 it appears that the model in Eq. (5) meets the dissipation constraint in a satisfactory manner.

The power spectrum model given by Eq. (6) is the first of the two main results of our paper. It can be used for any turbulent medium in which a single scalar is advected by the velocity field and is applicable for the wide range of the Prandtl numbers, *Pr* ∈ [3, 3000]. Moreover, it is confirmed that model (6) is in the reasonable agreement with the dissipation constraint.

## 3. The two-scalar oceanic refractive-index spectrum and optical scintillation

#### 3.1. Practical example: oceanic refractive-index spectrum

In this section we will demonstrate how to extend the power spectrum model given by Eq. (6) from one to two advected scalars with different Prandtl numbers and Schmidt numbers. The best illustration of such a medium is the oceanic turbulence, hence we will adopt the corresponding notations to maintain practicality.

On following the Nikishov’s linearized polynomial approach [11] we assume that the fluctuating portion of the refractive index *n* is given by linear combination

*T′*=

*T*− 〈

*T*〉 and the salinity concentration fluctuation

*S′*=

*S*− 〈

*S*〉,

*T*and

*S*representing the instantaneous temperature and salinity values and the angular brackets standing for the statistical average. In Eq. (8)

*A*is the thermal expansion coefficient and

*B*is the saline contraction coefficient [11].

Therefore, the power spectrum of the refractive index fluctuations can be expressed as the following linear combination of the temperature spectrum Φ_{T}(*κ*), the salinity spectrum Φ_{S}(*κ*), and the co-spectrum Φ_{TS}(*κ*) [11]:

*Pr*

_{T}and

*Pr*

_{S}are the temperature Prandtl number and salinity Schmidt number, respectively,

*Pr*

_{TS}= 2

*Pr*

_{T}

*Pr*

_{S}(

*Pr*

_{T}+

*Pr*

_{S})

^{−1}is the coupled Prandtl-Schmidt number [20,23], and

*χ*are the ensemble-averaged variance dissipation rates which are related by expressions [41]

_{i}*d*being the eddy diffusivity ratio and

_{r}*ω*being the relative strength of temperature-salinity fluctuations. In most practical cases,

*d*and

_{r}*ω*are related as [41]

The power spectrum model given by Eqs. (9) and (10) is the second of the two main results of our paper. It gives the analytic description of optical turbulence with two advected quantities while each of them may have Prandtl or Schmidt numbers in interval [3, 3000]. In its present form the new model is specifically applied for the oceanic waters in which the turbulent velocity field advects temperature and salinity concentration. However, it can be readily adapted for a variety of other turbulent media based on two scalars, for instance a fresh water contaminated by a chemical substance or any exotic turbulent fluid. Furthermore, a straightforward linearization of three or more advected quantities can be worked out in a similar manner.

#### 3.2. The scintillation index of a spherical wave

We will now illustrate how the analytical fit for the oceanic power spectrum given by Eqs. (9) and (10) compares with the Nikishovs’ spectrum (based on analytical model H1) and the H4 model (handled via numerical calculations) as it manifests itself in the evolution of the scintillation index of light propagating though an extended turbulent channel. The scintillation index of an optical wave at distance *L* from the source plane is generally defined by expression [6]

*I*is the instantaneous intensity of the wave, and the angular brackets stand for the ensemble average.

We will restrict ourselves to the scintillation index of a spherical wave being one of the most important parameters of light-turbulence interaction. In addition to the Rytov variance, i.e., the scintillation index of the plane wave, it can serve as an indicator of the global turbulence regime (weak/focusing/strong) [6]. This quantity was extensively studied for the atmospheric propagation in dependence from a number of turbulence parameters [6]. Based on the Nikishovs’ spectrum the scintillation index of the spherical wave was also previously evaluated on propagation in the oceanic water (c.f. [13]).

Let us now analytically derive the expression for the scintillation index of the spherical wave in the water channel described by spectrum in Eqs. (9) and (10). According to the Rytov perturbation theory it takes form

*i*also used in the paper),

*k*= 2

*πn*

_{0}/

*λ*

_{0},

*λ*

_{0}is the wavelength in vacuum,

*n*

_{0}(ca. 1.33) is the averaged refractive-index of the ocean water.

On substituting from Eq. (9) into Eqs. (14) and (15) we arrive at the scintillation index in of the form

_{3}

*F*

_{2}is the generalized Hypergeometric function [42]. Equations (16)–(18) combined with the Eqs. (11) and (12) imply that ${\sigma}_{I}^{2}$ is proportional to

*χ*

_{T}and

*ε*

^{−1/3}but show a somewhat complex dependence on

*ω*,

*Pr*

_{i}and

*η*.

In order to compare the values of ${\sigma}_{I}^{2}(L)$ corresponding to different oceanic power spectrum models, we set *n*_{0} = 1.33, *λ*_{0} = 532nm, *β* = 0.72, *A* = 2.56 × 10^{−4}deg^{−1}l, *B* = 7.44 × 10^{−4}g^{−1}l, *χ*_{T} = 1 × 10^{−5}K^{2}s^{−1}, *ε* = 1 × 10^{−2}m^{2}s^{−3}, but vary *ω* and average temperature 〈*T*〉. Varying 〈*T*〉 leads to the changes of *Pr _{i}* and

*η*(more details in the Appendix II). In Fig. 4 the scintillation index ${\sigma}_{I}^{2}(L)$ calculated from Eq. (16) is plotted for different values of

*ω*and 〈

*T*〉 (solid black curve or curve 2) and is compared with that calculated from the H4 model (dashed red curve or curve 1) and also with that given by the Nikishovs’ H1-based model (dotted blue curve or curve 3). Figure 4 suggests that if H4 is considered as the most precise model, a standard, then the scintillation index given in Eq. (16) provides a very accurate approximation to that based on H4, which is not the case for that based on the Nikishovs’ model. Such an advantage stems from the fact that the information about the bumps occurring in all three spectra, Φ

_{T}, Φ

_{S}and Φ

_{TS}, at high spatial wave numbers has been taken into the account by Eq. (16) at its best.

As Fig. 4 indicates, the average temperature is not in the direct relation with the scintillation index. Indeed, the scintillation index is directly affected by the variance of the fluctuating temperature rather than by its average value.

## 4. Summary

The power spectrum of the refractive index fluctuations is an effective tool for optical characterization of the statistics of turbulent medium (up to the second order) in its homogeneous regime. The widely accepted model H4 for the power spectrum of a single scalar advected by the turbulent fluid is based on a nonlinear differential equation and can always be used for numerical calculation of the spectrum. However, for light propagation problems it is desirable to have an accurate yet tractable analytical model for the power spectrum. For the turbulent media with low Prandtl numbers, as is the case for the atmospheric turbulence, such analytical models have been developed long time ago by fitting an analytic curve to the H4 profile. However, for turbulent media in which the scalar fields are advected at high and possibly variable Prandtl numbers, the accurate analytical fits have not yet been introduced.

In this paper we have introduced a very simple yet accurate analytical model for the power spectrum of turbulence with either one or two advected scalars. We have illustrated our procedure for the best example of such a medium - oceanic turbulence - and have obtained a very good agreement between such power spectrum and that based on the H4 numerical curve. We have also shown that the Nikishovs’ model widely used to characterize oceanic optical turbulence is much less accurate. As a consequence, our model’s close fit to the H4-based spectrum result in the close agreement for the statistics of the propagating light. In particular, we have illustrated that the scintillation index of a spherical wave as calculated from our power spectrum model is practically indistinguishable from that based on the H4 model. At the same rate, the Nikishovs’ spectrum based scintillation index is largely underestimated.

The key difference between our new power spectrum and the previously introduced ones stems from the fact that in our model the Prandtl/Schmidt numbers of the advected scalars are given as parameters that may vary independently in the interval *Pr* ∈ [3, 3000]. Such a large range produces a slight deviation of the new spectrum from the H4 model for all the *Pr* values which can be substantially reduced if a smaller range is set. However, this very large range of the *Pr* values leads to our model’s main advantage over all other models developed so far: the intrinsic flexibility in accounting for various thermodynamic states. Since, in general, the Prandtl/Schmidt numbers are very sensitive to a number of factors, such as the average temperature, pressure, concentration, etc., it is of utmost importance to have a power spectrum being flexible enough to account for these variations. Moreover, this very feature becomes of fundamental convenience when one may want to apply our model for homogeneous turbulence developed in a fluid other than the oceanic water. Hence we envision that our development may find uses in ecology, meteorology, food processing, medicine and some other technologies that use light propagation in a variety of turbulent media.

## Appendix I: Nikishovs’ H1-based model

Here we summarize the H1 model and the Nikishov’s model based on it (see [10,11,24,41] for more details).

H1 describes the three-dimensional scalar power spectrum of a single advected quantity as

*g*(

*κη*) is set to

## Appendix II: Oceanic turbulence parameters depending on average temperature

Table 3 gives the values of several average temperature-based parameters in oceanic refractive-index spectrum when salinity 〈*S*〉 = 34.9ppt. The values of kinematic viscosity *υ* and temperature Prandtl number *Pr*_{T} are calculated with the codes available in http://web.mit.edu/seawater/ [43,44]. The values of salinity diffusivity *α*_{S} for 〈*T*〉 = 15°C and 30°C are obtained from the diffusivity of Cl^{−} [45]; the *α*_{S} for 〈*T*〉 = 0°C is calculated from that for 15°C, 20°C, 25°C, 27°C and 30°C by the Stokes–Einstein law [45]. We calculate the Schmidt number by *Pr*_{S} = *ν*/*α*_{S} and the Kolmogorov microscale by *η* = *ν*^{3/4}*ε*^{−1/4}.

## Acknowledgments

The authors thank Mohammed Elamassie, Ming-Yuan Ren and Jian-Dong Zhang for their helpful comments, Wei-Qiang Ding for interesting discussions on numerical methods, and Andreas Muschinski for providing several manuscript documents. Thanks are also due to John Lienhard and his co-workers for developing the website http://web.mit.edu/seawater. We also thank Han-Tao Wang, Ying Sun, Zeng-Quan Yan, Han-Bing Wang and Fan Zhang for their early helpful work in this topic.

## References

**1. **W. Hou, “A simple underwater imaging model,” Opt. Lett. **34**, 2688–2690 (2009). [CrossRef] [PubMed]

**2. **X. Yi, Z. Li, and Z. Liu, “Underwater optical communication performance for laser beam propagation through weak oceanic turbulence,” Appl. Opt. **54**, 1273–1278 (2015). [CrossRef] [PubMed]

**3. **Y. Baykal, “Bit error rate of pulse position modulated optical wireless communication links in oceanic turbulence,” J. Opt. Soc. Am. A **35**, 1627–1632 (2018). [CrossRef]

**4. **Z. Cui, P. Yue, X. Yi, and J. Li, “Scintillation of a partially coherent beam with pointing errors resulting from a slightly skewed underwater platform in oceanic turbulence,” Appl. Opt. **58**, 4443–4449 (2019). [CrossRef] [PubMed]

**5. **O. Korotkova, “Enhanced backscatter in lidar systems with retro-reflectors operating through a turbulent ocean,” J. Opt. Soc. Am. A **35**, 1797–1804 (2018). [CrossRef]

**6. **L. C. Andrews and R. L. Phillips, *Laser Beam Propagation through Random Media, Second Edition* (SPIE, 2005). [CrossRef]

**7. **O. Korotkova, *Random Light Beams: Theory and Applications* (CRC, 2013).

**8. **O. O. Chumak and R. A. Baskov, “Strong enhancing effect of correlations of photon trajectories on laser beam scintillations,” Phys. Rev. A **93**, 033821 (2016). [CrossRef]

**9. **R. A. Baskov and O. O. Chumak, “Laser-beam scintillations for weak and moderate turbulence,” Phys. Rev. A **97**, 043817 (2018). [CrossRef]

**10. **O. Korotkova, “Light propagation in a turbulent ocean,” in *Progress in Optics*, vol. 64T. D. Visser, ed. (Elsevier, 2019), pp. 1–43. [CrossRef]

**11. **V. V. Nikishov and V. I. Nikishov, “Spectrum of turbulent fluctuations of the sea-water refraction index,” Fluid Mech. Res. **27**, 82–98 (2000). [CrossRef]

**12. **L. Wei, L. Liu, and J. Sun, “Influence of temperature and salinity fluctuations on propagation behaviour of partially coherent beams in oceanic turbulence,” J. Opt. A-Pure Appl. Opt. **8**, 1052 (2006). [CrossRef]

**13. **E. Shchepakina, N. Farwell, and O. Korotkova, “Spectral changes in stochastic light beams propagating in turbulent ocean,” Appl. Phys. B **105**, 415 (2011). [CrossRef]

**14. **O. Korotkova and N. Farwell, “Effect of oceanic turbulence on polarization of stochastic beams,” Opt. Commun. **284**, 1740–1746 (2011). [CrossRef]

**15. **N. Farwell and O. Korotkova, “Intensity and coherence properties of light in oceanic turbulence,” Opt. Commun. **285**, 872–875 (2012). [CrossRef]

**16. **O. Korotkova, N. Farwell, and E. Shchepakina, “Light scintillation in oceanic turbulence,” Waves Random Complex Media **22**, 260–266 (2012). [CrossRef]

**17. **Y. Baykal, “Scintillation index in strong oceanic turbulence,” Opt. Commun. **375**, 15–18 (2016). [CrossRef]

**18. **Y. Ata and Y. Baykal, “Structure functions for optical wave propagation in underwater medium,” Waves Random Complex Media **24**, 164–173 (2014). [CrossRef]

**19. **L. Lu, X. Ji, and Y. Baykal, “Wave structure function and spatial coherence radius of plane and spherical waves propagating through oceanic turbulence,” Opt. Express **22**, 27112–27122 (2014). [CrossRef] [PubMed]

**20. **J. Yao, Y. Zhang, R. Wang, Y. Wang, and X. Wang, “Practical approximation of the oceanic refractive index spectrum,” Opt. Express **25**, 23283–23292 (2017). [CrossRef] [PubMed]

**21. **R. J. Hill, “Spectra of fluctuations in refractivity, temperature, humidity, and the temperature-humidity cospectrum in the inertial and dissipation ranges,” Radio Sci. **13**, 953–961 (1978). [CrossRef]

**22. **R. J. Hill, “Optical propagation in turbulent water,” J. Opt. Soc. Am. **68**, 1067–1072 (1978). [CrossRef]

**23. **X. Yi and I. B. Djordjevic, “Power spectrum of refractive-index fluctuations in turbulent ocean and its effect on optical scintillation,” Opt. Express **26**, 10188–10202 (2018). [CrossRef] [PubMed]

**24. **R. J. Hill, “Models of the scalar spectrum for turbulent advection,” J. Fluid Mech. **88**, 541–562 (1978). [CrossRef]

**25. **Further, H4 agrees with the Kraichnan spectrum at high wave numbers, while H1 evidently deviates from it. The precision of the Kraichnan spectrum and, hence, of H4 is well established [26–29].

**26. **A. Muschinski and S. M. de Bruyn Kops, “Investigation of hill’s optical turbulence model by means of direct numerical simulation,” J. Opt. Soc. Am. A **32**, 2423–2430 (2015). [CrossRef]

**27. **D. Bogucki, J. A. Domaradzki, and P. K. Yeung, “Direct numerical simulations of passive scalars with pr>1 advected by turbulent flow,” J. Fluid Mech. **343**, 111–130 (1997). [CrossRef]

**28. **X. Sanchez, E. Roget, J. Planella, and F. Forcat, “Small-scale spectrum of a scalar field in water: The batchelor and kraichnan models,” J. Phys. Oceanogr. **41**, 2155–2167 (2011). [CrossRef]

**29. **D. J. Bogucki, H. Luo, and J. A. Domaradzki, “Experimental evidence of the kraichnan scalar spectrum at high reynolds numbers,” J. Phys. Oceanogr. **42**, 1717–1728 (2012). [CrossRef]

**30. **Y. Li, Y. Zhang, and Y. Zhu, “Oceanic spectrum of unstable stratification turbulence with outer scale and scintillation index of gaussian-beam wave,” Opt. Express **27**, 7656–7672 (2019). [CrossRef] [PubMed]

**31. **J. Churnside, “A spectrum of refractive turbulence in the turbulent atmosphere,” Opt. Acta Int. J. Opt. **37**, 13–16 (1990).

**32. **R. Frehlich, “Laser scintillation measurements of the temperature spectrum in the atmospheric surface layer,” J. Atmospheric Sci. **49**, 1494–1509 (1992). [CrossRef]

**33. **L. C. Andrews, “An analytical model for the refractive index power spectrum and its application to optical scintillations in the atmosphere,” Opt. Acta Int. J. Opt. **39**, 1849–1853 (1992).

**34. **K. J. Grayshan, F. S. Vetelino, and C. Y. Young, “A marine atmospheric spectrum for laser propagation,” Waves Random Complex Media **18**, 173–184 (2007). [CrossRef]

**35. **A. Muschinski, “Temperature variance dissipation equation and its relevance for optical turbulence modeling,” J. Opt. Soc. Am. A **32**, 2195–2200 (2015). [CrossRef]

**36. **The value of β may change with environment, and has always been used as an experimental constant. Several experiments gave its range 0.22 ∼ 0.78 [39]. The value 0.72 has been widely used in oceanic optical research, so we use the value 0.72 in Table 1 and Figure 2, but keep it as a parameter in all formulae.

**37. **P. Yue, X. Luan, X. Yi, Z. Cui, and M. Wu, “Beam-wander analysis in turbulent ocean with the effect of the eddy diffusivity ratio and the outer scale,” J. Opt. Soc. Am. A **36**, 556–562 (2019). [CrossRef]

**38. **The Gaussian tail is used here for historical reasons and mathematical convenience. Some well-known power spectrum models developed by Tatarskii [6] and Andrews [33] rely on the Gaussian profile of the high frequency cut-off. Such Gaussian profile is crucial for convergence of integrals involved in evaluation in statistics of light field propagation.

**39. **H. L. Grant, B. A. Hughes, W. M. Vogel, and A. Moilliet, “The spectrum of temperature fluctuations in turbulent flow,” J. Fluid Mech. **34**, 423–442 (1968). [CrossRef]

**40. **We separate (κη) and c in order to raise them to different powers. The applied mathematical recipe is to make the maximum value of g(κη) vary with c, and thus with Pr. Such variation is the key feature of H4 model [29]. If not separated, c and Pr only scale the curve of g(κη) horizontally but do not change its maximum value.

**41. **M. Elamassie, M. Uysal, Y. Baykal, M. Abdallah, and K. Qaraqe, “Effect of eddy diffusivity ratio on underwater optical scintillation index,” J. Opt. Soc. Am. A **34**, 1969–1973 (2017). [CrossRef]

**42. **Y. L. Keke, *The special functions and their approximations* (Academic, 1969).

**43. **K. G. Nayar, M. H. Sharqawy, L. D. Banchik, and J. H. Lienhard V, “Thermophysical properties of seawater: A review and new correlations that include pressure dependence,” Desalination **390**, 1–24 (2016). [CrossRef]

**44. **M. H. Sharqawy, J. H. Lienhard V, and S. M. Zubair, “Thermophysical properties of seawater: a review of existing correlations and data,” Desalination Water Treat. **16**, 354–380 (2010). [CrossRef]

**45. **A. Poisson and A. Papaud, “Diffusion coefficients of major ions in seawater,” Mar. Chem. **13**, 265–280 (1983). [CrossRef]