Extensive, large-scale simulations of laser beam propagation in the atmosphere show that the nonlinear focusing properties of air can be used to deliver relatively high fluence (energy per unit area) to a remote target embedded in turbulence [1–3]. These simulations, while physically accurate, can require considerable computing resources to calculate the evolution of ensemble average quantities, such as the spot size. Here, a reduced model based on the paraxial, enveloped wave equation was employed to describe three-dimensional propagation of a powerful laser beam in a random dielectric, incorporating the effects of diffraction, dissipation, and Kerr nonlinearity. The model led to differential equations that can be readily solved, enabling rapid variation of parameters and the development of theoretical results to guide full-scale simulations and experiments. After outlining the theory, a comparison with field experiments and simulations is presented to illustrate the utility and limitations of the model and to demonstrate a novel effect: chromatic aberration due to self-focusing in the presence of dissipation.

The dielectric can be described by the effective refractive index $n(\omega )=⟨n_0(\omega )⟩+i⟨n_i(\omega )⟩+\delta n_0(\omega )+n_{2}I$, where $n_0(\omega )$ is real linear refractive index evaluated at the frequency $\omega$, $n_i(\omega )$ is the effective imaginary refractive index including the effects of absorption and large angle scattering, $⟨⟩$ represents an ensemble average, $\delta n_0(\omega )$ describes fluctuations in the real refractive index with $⟨\delta n_0⟩=0$, $n_2$ is the real, adiabatic nonlinear refractive index, and $I$ is the beam intensity defined below.

The electric field $\mathbf{E}(\mathbf{r},t)=\frac{1}{2}E{e}^{i\psi }{\mathbf{e}}_x+\text{c.c.}$ is written in terms of complex amplitude $E(\mathbf{r},t)$ and rapidly varying phase $\psi (z,t)=k_{0}z-{\omega }_{0}t$, where $k_0=2\pi /{\lambda }_0$ is the wavenumber, ${\lambda }_0$ is the wavelength, and ${\omega }_0$ is the frequency. Upon transforming $z$, $t\to z$, $\tau =t-z/{\mathrm{v}}_g$, where ${\mathrm{v}}_g$ is the group velocity, and neglecting finite pulse length effects, the propagation equation reduces to [4]:

The method of moments [5–8] with weight functions equal to the intensity $I({\mathbf{r}}_{\perp },z)\equiv c⟨n_0⟩{|E({\mathbf{r}}_{\perp },z)|}^2/8\pi$ was used to obtain equations of motion for the centroid coordinate $\mathbf{R}(z)=\iint {\mathrm{d}}^2{r}_{\perp }{\mathbf{r}}_{\perp }I({\mathbf{r}}_{\perp },z)/N(z)$, momentum $\mathbf{P}(z)=\iint {\mathrm{d}}^2{r}_{\perp }\mathit{S}({\mathbf{r}}_{\perp },z)$, long-term radius squared $R_L^2(z)=\iint {\mathrm{d}}^2{r}_{\perp }{\mathbf{r}}_{\perp }^{2}I({\mathbf{r}}_{\perp },z)/N(z)$, and Hamiltonian $H(z)=c⟨n_0⟩\iint {\mathrm{d}}^2{r}_{\perp }({\nabla }_{\perp }{E}^{*}·{\nabla }_{\perp }E-{\overline{n}}_2{|E|}^4)/16\pi$. Here, $N(z)=\iint {\mathrm{d}}^2{r}_{\perp }I({\mathbf{r}}_{\perp },z)$ and $\mathit{S}({\mathbf{r}}_{\perp },z)=c⟨n_0⟩(E^{*}{\nabla }_{\perp }E-E{\nabla }_{\perp }{E}^{*})/16\pi i$ was the flux density.

Ensemble averaging was affected by employing the Markov approximation in which refractive index fluctuations are $\delta -$ correlated in the direction of propagation; that is, the covariance is expressed as $⟨\eta ({\mathbf{r}}_{\perp },z)\eta ({\mathbf{r}}_{\perp }^{\prime },z^{\prime })⟩=C_{\eta }({\mathbf{r}}_{\perp }-{\mathbf{r}}_{\perp }^{\prime },z)\delta (z-z^{\prime })$, where $C_{\eta }({\mathbf{r}}_{\perp }-{\mathbf{r}}_{\perp }^{\prime },z)$ is the second-order cumulant [9–11]. Then, $⟨N(z)⟩=N(z)$, $d⟨\mathbf{R}(z)⟩/dz=⟨\mathbf{P}(z)⟩/k_{0}N(z)$, and $(d/dz+\gamma )⟨\mathbf{P}(z)⟩=c⟨n_0⟩\iint {\mathrm{d}}^2{r}_{\perp }{\mathrm{\Gamma }}_4({\mathbf{r}}_{\perp },{\mathbf{r}}_{\perp };{\mathbf{r}}_{\perp },{\mathbf{r}}_{\perp },z){\nabla }_{\perp }{\overline{n}}_2({\mathbf{r}}_{\perp },z)/16\pi k_0$, along with

Consider the initial field $E({\mathbf{r}}_{\perp },z=0)=E_0\exp (-r_{\perp }^2/w^2+i{k}_0{r}_{\perp }^2/2F)$, where $E_0$ is the amplitude, $w$ is the spot radius, $F$ is the radius of curvature, and $P=c⟨n_0⟩E_0^2{w}^2/16$ is the power. Then $N(z)=N_0\exp (-\gamma z)$, where $N_0=\pi E_0^2{w}^2/2$ and ${\mathrm{\Gamma }}_4$ is computed by assuming the intensity remains Gaussian in form. For ${\nabla }_{\perp }{\overline{n}}_2\to 0$, Eqs. (1b) and (1c) reduce to

A solution requires knowledge of the cumulant $C_{\eta }({\mathbf{r}}_{\perp }-{\mathbf{r}}_{\perp }^{\prime },z)=\int \mathrm{d}{z}^{\prime }⟨\eta ({\mathbf{r}}_{\perp },z)\eta ({\mathbf{r}}_{\perp }^{\prime },z^{\prime })⟩=\iint {\mathrm{d}}^2{\kappa }_{\perp }{\mathrm{\Phi }}_{\eta }({\mathit{\kappa }}_{\perp },{\kappa }_z=0;z)\exp (i{\mathit{\kappa }}_{\perp }·{\mathbf{r}}_{\perp })$, where ${\mathrm{\Phi }}_{\eta }({\mathit{\kappa }}_{\perp },{\kappa }_z;z)$ denotes the spectrum of fluctuations of $\eta$ in $({\mathit{\kappa }}_{\perp },{\kappa }_z)$ wavenumber-space. Henceforth, propagation in an atmosphere with a Kolmogorov-von Kármán spectrum is considered; that is, ${\mathrm{\Phi }}_n(\mathit{\kappa };z)=A{C}_n^2(z){({\kappa }^2+{\kappa }_0^2)}^{-11/6}{e}^{-{\kappa }^2/{\kappa }_m^2}$ [9–11], where $C_n^2(z)$ is the refractive index structure constant, $A=5/18\pi \mathrm{\Gamma }(1/3)$, ${\kappa }_0=2\pi /L_0$, $L_0$ is the outer scale length, ${\kappa }_m=5.92/{\ell }_0$, and ${\ell }_0$ is the inner scale length.

Equation (2), derived using the Markov approximation in the framework of the paraxial envelope equation, is valid if (i) ${\lambda }_0\ll {\ell }_0$, (ii) ${\lambda }_0\beta \ll 1$, and (iii) propagation range $\ll z_m$. Here, $\beta =2{\pi }^{2}A{C}_n^2{k}_0^2{\kappa }_m^{-5/3}\mathrm{\Gamma }(-5/6,{\kappa }_0^2/{\kappa }_m^2)e^{{\kappa }_0^2/{\kappa }_m^2}$ is the scattering coefficient from eddies, and $z_m=({\lambda }_0/4C_n^2{L}_0^{2/3}){({\ell }_0/L_0)}^{1/6}$ is the characteristic distance over which the intensity changes significantly due to backscattering [9].

The cumulant term in Eq. (2b) is expressible in terms of Laguerre polynomials, $2{\pi }^{2}A{C}_n^2(z){\kappa }_m^{-5/3}\mathrm{\Gamma }(-5/6){\nabla }_{\perp }^2[L_{5/6}(-{\kappa }_m^2{r}_{\perp }^2/4)]$ (for ${\kappa }_0\ll {\kappa }_m$). Care must be exercised in evaluating this term. The approach here is informed by the procedure employed to obtain the structure function of the complex phase in [9]. In this procedure, the transverse coherence length $r_0(z)={(4/3M{C}_n^2{k}_0^{2}z)}^{3/5}$, where $M=(6/5){\pi }^{2}A\mathrm{\Gamma }(1/6)/2^{5/3}\mathrm{\Gamma }(11/6)$, plays a central role. When $z$ is sufficiently large, $r_0(z)\ll {\ell }_0$, the Laguerre polynomial can be expanded for ${\kappa }_m{r}_{\perp }\ll 1$, and the long-term mean-square spot radius is

In the absence of self-focusing (i.e., $P/P_{\mathrm{NL}}\to 0$) the limit in Eq. (3a) is given in [9–11,13,14] while the limit in Eq. (3b) has been given in [10,11]. Self-focusing (i.e., $P/P_{\mathrm{NL}}>0$) in the limiting case of Eq. (3a) is discussed in [15]. Equation (3) does not capture the long-term spreading predicted by weak fluctuation theory (Rytov approximation). Nonetheless, as noted in [11], Eq. (3b) is in excellent agreement with experiments (assuming $P/P_{\mathrm{NL}}\to 0$).

(i) Low-power benchmark against experiment at AFRL. Self-focusing is unimportant in the low-power limit. Figure 1 shows the long-term spot radius versus $C_n^2$ for a kW-class, single-mode, $\lambda =1070\mathrm{nm}$ fiber laser with initial spot radius $w=4.5\mathrm{cm}$ focused at a range of 2.2 km by a lens with focal length $|F|\sim 5\mathrm{km}$. For these parameters, the critical power for self-focusing is multiple gigawatts. This propagation test was conducted over a horizontal path at the Air Force Research Laboratory, Kirtland Air Force Base, NM using a Cassegrain telescope with reflective optics. Atmospheric turbulence was driven by the heating of the ground by the sun, and it varied significantly with the time of day. A DIMM (differential image motion monitor) diagnostic was used to measure the transverse coherence length over the 2 km path at 1-minute intervals [16]. The experimental results were recorded by a camera operating at 22 frames/s and each data point (black circle) was obtained by averaging over one minute (1320 frames). Each of the orange circles was the result of averaging ${10}^3$ runs using the full-scale simulation code HELCAP [1]. The dashed green curve was obtained using Eq. (37) of Ref. [10]. Equation (3b) applied and was represented by the blue curve in Fig. 1. Turbulence-induced spreading was observed to be largest in the simulations, followed by Eq. (37) of Ref. [10] and then by Eq. (3b). The closed-form analytical expressions for the spot radius employed a number of approximations and resulted in some arbitrariness in the numerical coefficients of the turbulence term; as a consequence, the green and the blue curves do not overlap. Nevertheless the moment approach does appear to provide a reasonable fit to the measurements.

 

Fig. 1. Long-term spot radius versus refractive index structure constant for a low-power beam (i.e., in the absence of self-focusing). Black circles are data points from field experiments at the Air Force Research Laboratory, NM. Orange circles were obtained using HELCAP full-scale simulations. The analytical results are shown in (i) dashed green, obtained by using Eq. (37) of Ref. [10], and (ii) blue, obtained using Eq. (3b) here.

Download Full Size | PDF

(ii) High-power benchmark against HELCAP simulation. To illustrate the effects of self-focusing and turbulence, Fig. 2 compares Eq. (3b) with simulation results from HELCAP. A $\lambda =1\mathrm{\mu m}$ laser with initial spot radius $w=10\mathrm{cm}$ was focused by a lens with focal length $|F|=8\mathrm{km}$ in air in the absence of extinction ($\gamma =0$). HELCAP results for the spot radius—shown in blue—were plotted at low power (solid curves) and at $P/P_{\mathrm{NL}}=0.8$ (dashed curves), deliberately limiting laser power to be far below the threshold for filamentation. Curves obtained from simulations represented the ensemble average of ${10}^3$ statistically independent realizations of turbulence. The corresponding theory results, Eq. (3b), are shown in red. Figure 2(a) considers the case of strong turbulence, $C_n^2={10}^{-15}{\mathrm{m}}^{-2/3}$. In this case, both theory and simulation showed little difference between high-power and low-power propagation. This behavior can be understood by considering the second-order coherence function ${\mathrm{\Gamma }}_2({\mathbf{r}}_{\perp 1},{\mathbf{r}}_{\perp 2},z)$ and the spatial coherence radius derived therefrom (defined following Eq. (56), chap. 6, Ref. [11]). The HELCAP results—the solid black circles in Fig. 2(a)—show that the coherence radius in this case was smaller than the spot beam radius. The high degree of incoherence associated with the turbulence inhibits nonlinear self-focusing. Alternatively, for weaker turbulence, $C_n^2={10}^{-16}{\mathrm{m}}^{-2/3}$, see Fig. 2(b), the spatial coherence radius was larger than the spot beam radius, and theory and simulation both showed the beam focusing to a smaller spot than for low-power propagation through turbulence. Hence, nonlinear focusing was effective at producing smaller spots and corresponding higher intensities in a turbulent atmosphere provided that the beam remains relatively coherent over the range. As in Fig. 1, for the cases shown in Fig. 2, the moment theory predicts a smaller focal spot than the simulations. (With no turbulence, theory and simulation predict the same spot radius at 8 km).

 

Fig. 2. Long-term spot radius versus propagation distance for a beam focused at 8 km. Blue curves are the results of HELCAP simulations, and the red curves were generated using Eq. (3b). Dashed curves represent cases where $P/P_{\mathrm{NL}}=0.8$, while solid curves represent linear propagation; that is, $P/P_{\mathrm{NL}}\to 0$. Black dots denote the spatial coherence radius computed in HELCAP. In panel (a) $C_n^2={10}^{-15}{\mathrm{m}}^{-2/3}$, and in panel (b) $C_n^2={10}^{-16}{\mathrm{m}}^{-2/3}$.

Download Full Size | PDF

Because of the $z$-dependence, the turbulence term in Eq. (3) eventually dominates all others and the beam spreads without bound. For $P/P_{\mathrm{NL}}>1$, however, the radius can initially decrease. In fact, based on Eq. (3b) for a collimated beam ($F\to \infty$) there is a critical level of turbulence $C_{n,\text{crit}}^2\approx 0.75{[{⟨n_0⟩}^5{(P/P_{\mathrm{NL}}-1)}^4/M^{23/6}]}^{1/3}/k_0^3{w}^{11/3}$ below which the beam completely pinches (completely collapses) in a finite distance before turning round and eventually spreading. This peculiar behavior is one manifestation of the limitations associated with analytical approaches [17]. Of course, simulations are also bedeviled by gridding issues as one approaches the collapsed stage. As discussed in Ref. [6], when self-focusing occurs, the near-axis “rays” converge to a focus before the rays further out in the wings of the intensity distribution do. The moment approach prediction for the collapse distance combines these varied behaviors, and generally overpredicts the collapse distance [6]. Furthermore, the moment approach predicts the evolution of the average spot size, a positive definite quantity. If the average spot size collapses to zero, then the spot size in any particular instance must be identical to zero. As demonstrated in [1], this is not the case. Indeed, the beam can be highly filamented for $P/P_{\mathrm{NL}}\gg 1$ [1], and then an averaged spot size description is not meaningful.

(iii) Atmospheric extinction and chromatic aberration. Figure 3 demonstrates one of the effects of dissipation obtained by solving Eq. (2) for $\lambda =1$, 3, 5 μm, with (a) $\gamma =0$ and (b) $\gamma =1{\mathrm{km}}^{-1}$ and all other parameters being the same as in Fig. 2. For $\gamma =0$ the radii for $\lambda =1$, 3 μm are comparable because the $\lambda =3\mathrm{\mu m}$ beam experiences more diffractive spreading, but it is more resistant to turbulent spreading (i.e., it has a larger coherence length). Comparison of Figs. 3(a) and 3(b) shows that the fundamental is least affected by extinction, because while extinction does not modify turbulent spreading, damping the beam power reduces self-focusing and thus increases diffractive spreading—the more so, the longer the wavelength. Note that in practice, $\gamma$ can vary significantly with wavelength [18]. The minimum spot radii and focal points, comparing low-power and high-power cases in the presence of damping, are listed in Table 1.

 

Fig. 3. Long-term spot radius versus propagation distance for $P/P_{\mathrm{NL}}=0.8$ and wavelengths $\lambda =1\mathrm{\mu m}$ (blue, solid), $\lambda =3\mathrm{\mu m}$ (green, short dashed), and $\lambda =5\mathrm{\mu m}$ (red, long dashed). (a) Extinction rate $\gamma =0$ and (b) extinction rate $\gamma =1{\mathrm{km}}^{-1}$. A beam with initial spot radius $w=10\mathrm{cm}$ was focused by a lens with focal length $|F|=8\mathrm{km}$ in air with $C_n^2={10}^{-15}{\mathrm{m}}^{-2/3}$.

Download Full Size | PDF

Table 1. Minimum Spot Radius at Focal Point for Three Wavelengths, Comparing the Low-Power Case with the High-Power Case, in the Presence of Extinction $\gamma =1{\mathrm{km}}^{-1}$

View Table

The moment approach is used to analyze high-power laser beam propagation in a turbulent, dissipative medium. Comparisons with low-power field experiments as well as with simulations in the high-power regime are presented. It is shown that dissipation can lead to chromatic aberration.