1.

Introduction

The proliferation of satellites in Earth orbit, increasing in both number and economic value, makes the problem of collision with orbital debris very pressing. A promising solution for this problem is debris removal with the help of a ground-based pulsed laser. In this approach, laser pulses ablate debris material, changing the debris velocity and moving the debris to a lower orbit where natural burn-up takes place.^1,2 One problem with energy delivery to the orbit is self-focusing in the atmosphere. The nonlinear effects typically are not strong enough to destroy the laser beam in the atmosphere, but the acquired phase distortion is sufficient to filament the beam during propagation to debris (see Fig. 1), thus degrading the efficiency of the debris cleaning.³

Fig. 1

The normalized intensity distribution $I(r,z)/I(\mathrm{0,0})$ {with [$I(r,z)={|A(r,z)|}^2$]}. (a) $P/P_{\mathrm{cr}}=50$, (b) $P/P_{\mathrm{cr}}=1000$, (c) $P/P_{\mathrm{cr}}=2500$, and (d) $P/P_{\mathrm{cr}}=5000$. These conditions correspond to a focusing distance $L=1000\mathrm{km}$ for linear propagation.

We start the analysis following⁴ with the basic requirements for the laser pulse on the target. Then, we discuss beam propagation and focusing, completing the requirements for the laser pulse. We will see that in typical situations the pulse power must be hundreds of times larger than the critical power for self-focusing. Then we consider the thin window (TW) model^5,6 for a description of nonlinear propagation through the atmosphere. Within this model, the propagation through the atmosphere produces only phase distortions. Our numerical modeling⁶ demonstrates that the TW model describes surprisingly good beam propagation through the atmosphere even in the case of beam filamentation. The applicability of the TW model yields simple expressions for the average beam description on the target: the Strehl ratio and the beam quality $M^2$.

If propagation through the atmosphere produces only phase distortion, we suggest that preimposing the phase correction on the beam with the help of the phase mask will compensate for the detrimental effect of the atmosphere. Our modeling demonstrates that a complete compensation of the detrimental nonlinear effects is possible.

The phase correction to compensate nonlinear effects depends on the density profile in the atmosphere, which is not well defined, the laser intensity, which can fluctuate, etc. As a result, it is important to present the detailed analysis of the compensation sensitivity. During the long, free propagation to the orbit even small phase disturbances can greatly modify the beam. The study of this sensitivity is the main subject of the paper.

The big problem with sensitivity studies is to find a good figure of merit (FOM). We selected Strehl ratio and the beam quality as FOM. These are the parameters characterized laser intensity and spot size on the target. Also, in Ref. 6, we derived the simple relation that calculated these quantities in terms of laser field existing in the atmosphere, solving the laser propagation problem to the orbit. The FOM approach provides the quantative measure of the sensitivity to phase self-modulation, change of atmospheric density profile, laser altitude, etc.

We conclude that the suggested phase mask method is very robust and that its implementation is an attractive way to increase the laser system performance.

2.

Interaction of Laser Radiation with Debris

High intensity pulsed laser radiation incident on debris vaporizes surface material, creating recoil momentum that changes the debris velocity. It is clear that an optimal laser intensity exists for any specified pulse duration. At low intensity, surface temperature and evaporation are low and recoil momentum is small. At high intensity, a large fraction of laser energy is used to create plasma, which contributes little to momentum change of the debris. A crucial parameter for pulsed laser debris removal is the coupling coefficient $C_m$, the ratio of momentum imparted to the target, to the incident laser energy, $C_m=\mathrm{\Delta }P/E$. A review of data of $C_m$ dependence on intensity for different materials is presented in Refs. 1 and 7. The experimental data from different groups demonstrate that for broad ranges of wavelength, pulse duration, and pulse energy, the coupling coefficient maximum is reached at the intensity

The fluence corresponding to optimal coupling is given as

Let us derive the requirement for laser pulse energy that corresponds to delivering the optimal fluence to debris targets. The energy $E$ delivered by the laser to the target is required to be $E=F\pi r^2$, where $r$ is the radius of the beam in the focal plane. For beam radius that accounts for beam quality and beam diffraction,

The required laser pulse energy $E$ for delivering the pulse fluence for optimal coupling is found by combining Eqs. (1) and (2) giving⁴

We now consider a specific example, in which $\lambda =1 \hspace{0.17em}\mu \mathrm{m}$, $L=1000\mathrm{km}$, $D=3\mathrm{m}$, and $M^2=2$, a value of $M^2$ that can be achieved for high-energy lasers by using spatial filters and adaptive-optic systems. For this case, $r\sim 34\mathrm{cm}$, the required pulse energy is $E=11\sqrt{\tau (\mathrm{ns})}\mathrm{kJ}$, and pulse power $P=1.1\times {10}^{13}/\sqrt{\tau (\mathrm{ns})}$. For the optimal point of view from laser technology,⁴ pulse duration $\sim 3\mathrm{ns}$ the required power $P=6\mathrm{TW}$. This is much higher than the critical power for self-focusing in atmosphere, $P_{\mathrm{cr}}=4\mathrm{GW}$, which corresponds to $P/P_{\mathrm{cr}}=1500$. The required power $P$ scales as $P\sim M^4{L}^2/D^2$ and the ratio of $P/P_{\mathrm{cr}}$ for practical values of parameters is above 1000.

3.

Laser Beam Propagation Modeling

Specific modeling is needed to describe the long distance propagation between the atmosphere and the orbit, which necessitates high accuracy of nonlinear effects modeling. We will use a code applied and validated for modeling of this problem.³

The propagation of the laser beam is described by the nonlinear Schrodinger equation,^2,7 i.e.,

The resulting normalized equation in dimensionless variables has the form⁶

It is possible to show that the fastest growing perturbations resulting in filamentation are axisymmetric,¹⁰ and that filamentation, at least initially, breaks the beam into ring-like structures. The formation and breakup of the ring structure are the well-documented pathway to Gaussian beam filamentation,¹¹ and beam propagation can be described in the main approximation within the axisymmetric version of Eq. (3).

Let us consider the propagation of an initially Gaussian laser beam. On the surface at $z=0$

Here, $P_{\text{in}}$ is the normalized input power of the laser beam, and the dimensionless parameter $C=k_0·R_0^2/F=L_D/2F$ is the initial beam prefocusing parameter. $I_0=I(\mathrm{0,0})=P_{\text{in}}[W]/(\pi P_0[W])=4P_{\text{in}}/P_{\mathrm{cr}}$. $F$ has the meaning of a focal distance that in this case is the debris height $L$. Therefore, dimensional initial prefocusing is given by $C=L_D/2L$.

We would like to stress that the problem under consideration, though similar in terms of the basic equations to numerous self-focusing studies (see, e.g., Refs. 3, 12 and references therein), is rather different in terms of underlying physics. We consider here light propagation over a finite distance (the thickness of the nonlinear layer), with the focusing point located beyond this region, where the propagation is linear. In this case, the self-focusing effect compresses the beam, but without the catastrophic collapse of all the energy into a small volume. Yet even the small phase disturbances after the atmospheric pass during the long propagation to the target can greatly modify the beam. The TW model developed in Ref. 6 is a powerful tool to evaluate the beam structure.

4.

Thin Window Model

The nonlinear effects are important only in a thin layer of the atmosphere (in comparison with the distance to the target); see Fig. 2.⁶ After passage through the atmosphere, we have free propagation wherein the effects acquired due to nonlinear interaction beam can modify the beam. The TW makes evaluation of this process quite simple.

Fig. 2

Illustration of the TW model.

For a Gaussian initial beam, the laser field after propagation through the atmosphere has the explicit description⁶

The value of $z_1$ is approximately a few times $h$; the exact value must be evaluated after comparison with numerical modeling results. It is clear that there must exist an optimal value, since for small $z_1$ we cut out the part of the atmospheric propagation, and if $z_1$ is too large free propagation will modify the solution (the window is no longer thin). Due to the exponential dependence in Eq. (5), the optimal value is about a few atmospheric thicknesses. The detailed studies presented in Ref. 6 give the optimal value of $z_1$ as $\sim 5h$. Following results will be presented for this $z_1$ value.

After the beam exits the atmosphere, we have free linear propagation (see Fig. 2). From Eq. (5), we see that the phase is not quadratic, and the ensuing beam propagation is not described by the simple formulae available for the focused Gaussian input. We see that the curvature of the phase corresponds to additional focusing, and the atmosphere serves as a focusing astigmatic lens. As a result, the maximal field intensity is reached before the linear focal plane. By varying the prefocusing parameter, we can partially compensate for propagation through the atmosphere.³

The free propagation with complex enough phase distribution is not easy to evaluate. But for average beam characteristics in the focal plane (the peak intensity and the average beam spot size), it is possible to get compact expressions using the fact that the field in the focal plane is proportional to the Fourier transform of the field exiting the atmosphere.⁶

It is convenient to use the Strehl ratio, the ratio of the intensity in the center after propagation of the field, to that produced by the linear propagation of the Gaussian beam [Eq. (5)]. The electric field in the center of the focal spot is proportional to $\int U(r)r\mathrm{d}r$ giving for $S$⁶

Recalling that $U(r)=\sqrt{I_0}\exp (-\frac{r^2}{2}+i\mathrm{\varphi })$, $U_G(r)=\sqrt{I_0}\exp (-\frac{r^2}{2})$, we can calculate the Strehl ratio explicitly getting the analytical expression for $S$⁶

The integral estimate of the square of the spot size on the target can be calculated in a similar way, using the Parseval formula⁶ as follows:

Similarly, it is convenient to calculate the beam quality parameter $M^2$, the ratio spot size squared [Eq. (7)] to the value calculated for the Gaussian beam. For the exit from the atmosphere field distribution [Eq. (5)], $M^2$ can be evaluated analytically, as follows:⁶

The above-presented expressions for $S$ and $M^2$ give a quantitative description of the beam degradation due to the self-phase modulations. We got this advance in beam evaluation because, different from laser system studies, we are interested in the intensity distribution in the focal plane only.

The derived equations can have multiple applications. One can estimate quantitatively the detrimental effect of self-focusing in a case without compensation, or evaluate the effect of the finite aperture.

The fact that all effects of nonlinear interaction with atmosphere produce only phase distortion suggests that if we impose the initial phase equal $\mathrm{\varphi }$ [Eq. (5)] with the opposite sign, the propagation through the atmosphere will compensate the nonlinear phase distortion. During the long propagation to the orbit, even a small deviation from the TW model can produce a large deviation. We proceed to evaluate in detail the applicability of the model.

5.

Results and Discussion

Figures 3Fig. 4Fig. 5Fig. 6–7 show that the results obtained by using the simple TW model are in good agreement with the direct nonlinear Schrödinger equation (NLSE) simulations. Here $h=0.000506$, $C=k_0{R}_0^2/L=5.93$, $L=1000\mathrm{km}$, $R_0=1\mathrm{m}$, ${\lambda }_0=1.06\mu \mathrm{m}$, and $P_{\mathrm{in}}/P_{\mathrm{cr}}$ varies from 100 to 7500, which corresponds to $b=0.2$ and $b=15$, respectively. Let us note that the nonuniformity is important: even for $P_{\mathrm{in}}/P_{\mathrm{cr}}$ few thousands, the $b$ values are about 10 and we do not expect the massive beam filamentation.

Fig. 3

(a, b) Comparison of Strehl ratio $S$ and beam quality parameter $M^2$ computed by the NLSE [Eq. (3) (black line)] and TW model [Eqs. (6), (8) (red line)]. Green points correspond to the values computed by NLSE [Eq. (3)] with laser elevation at the height $z_0=3\mathrm{km}$. (c, d) Errors calculations, with $\text{error}(S)=|S_{\mathrm{NLSE}}-S_{\mathrm{TW}}|$ and $\text{error}(M^2)=|M_{\mathrm{NLSE}}^2-M_{\mathrm{TW}}^2|/M_{\mathrm{TW}}^2$.

Fig. 4

(a, b) $P_{\mathrm{in}}/P_{\mathrm{cr}}=2500$—conditions correspond to Fig. 1(c); (c, d) $P_{\mathrm{in}}/P_{\mathrm{cr}}=5000$—conditions correspond to Fig. 1(d). (a, c) Intensity distribution along $z$. (b, d) Radial distribution at the initial focal point $L=1000\mathrm{km}$, $C=5.93$. Black line: solution of NLSE [Eq. (3)] with initial condition [Eq. (4)]. Green line: solution of TW model [Eq. (5)]. Red line: NLSE with preimposed phase [Eq. (9)]. Blue line: propagation of Gaussian beam in linear case.

Fig. 5

(a, c) Intensity distribution along $z$. (b, d) Radial distribution at the initial focal point $L=1000\mathrm{km}$, $\delta =\mathrm{\Delta }/b$. (a, b) $b=2.01$ ($P_{\text{in}}/P_{\mathrm{cr}}=1000$)—corresponds to Fig. 1(b). (c, d) $b=10.05$ ($P_{\text{in}}/P_{\mathrm{cr}}=5000$)—corresponds to Fig. 1(d).

Fig. 6

The normalized intensity distribution $I(r,z)/I(\mathrm{0,0})$ {with [$I(r,z)={|A(r,z)|}^2$]}. (a) Propagation of Gaussian beam in linear case, (b) $\delta =0$, (c) $\delta =-0.1$, and (d) $\delta =0.1$. $\delta =\mathrm{\Delta }/b$. $C=5.93$ and $b=10.02$—the conditions correspond to Fig. 1(d).

Fig. 7

(a) Strehl ratio as a function $b$ and $\delta =\mathrm{\Delta }/b$.Strehl ratio varies from 0.4 (blue) to 1 (red). (b) Beam quality $M^2$ as a function $b$ and $\delta =\mathrm{\Delta }/b$. $M^2$ varies from 0.997 (red) to 6.5 (blue). Red color means that the Strehl ratio or $M^2$ is close to 1, i.e., the beam is close to Gaussian.

Let us discuss the calculations of the beam average parameters in the focal plane, i.e., Strehl ratio $S$ and beam quality $M^2$. For the TW model both quantities can be calculated analytically, as in Eqs. (6) and (8). The comparison of the Strehl ratio computed with the NLSE solution and the TW model is shown in Fig. 3(a). We see that the TW model is very close to the NLSE solution, reproducing even nonmonotonic $S$ behavior. Let us stress again that the TW model is accurate for the calculation of intensity in spot center even at large $b$ when the beam is far from Gaussian. The calculations of beam quality $M^2$ are presented in Fig. 3(b), together with analytical expression in Eq. (8). It is important to note that this result is also valid for beams quite different from Gaussian. We see that the TW model provides an excellent description at modest $b$, and slightly overestimates beam quality for large $b$, when the beam is already completely destroyed.

Let us mention that in the TW model the atmosphere characterized only by one parameter $b$. Within the NLSE, we can start with different values of intensities and atmospheric heights but the same values of $b$. The coincidence of these runs additionally validates the TW model. In Fig. 3, we present the comparison of the propagation starting from the sea level and from 3-km elevation. It is a practically interesting situation. To minimize the effect of atmospheric turbulence, it is attractive to place the laser as high as possible.¹ One can see the very good TW applicability. With increase of laser elevation, $b$ value decreases, decreasing the self-phase modulation effects.

In Fig. 4, we present the radial intensity distribution and the peak intensity distribution along $z$ for the TW model of Eq. (5) and the corresponding NLSE solution. Here $P_{\mathrm{in}}/P_{\mathrm{cr}}=1500$ and 5000, with the corresponding factor nonlinear phase shift $b=3.02$ or $b=10.05$, respectively.

Figure 4 shows an excellent agreement between the TW model and the full numerical simulations based on the NLSE, both for evolution of the peak intensity with distance and for the radial beam intensity distribution in the focal plane. One can see that the TW model approximates well the exact solution of the NLSE even in the situation with well-developed filamentation [see Figs. 1(c), 1(d), and 4] with the field distribution being very far from the Gaussian beam. Note that in the TW model, the solution depends on the dimensionless parameter $b=I_{0}h[1-\exp (-z_1/h)]\hspace{0.17em}$, only, which simplifies the system optimization.

The excellent results of the TW model give hope that the initial phase predistortion can be used to compensate the nonlinear phase changes. As a result, one can have an almost perfect Gaussian beam at the atmospheric exit and the detrimental effects of self-focusing can be eliminated to a great extent. The new initial condition with the corrected phase is

We compare the solution of the NLSE with the initial condition [Eq. (4)] and chirp $C=5.93$, which corresponds to linear focusing at $L=1000\mathrm{km}$; the solution of the NLSE [Eq. (3)] with a preimposed phase [Eq. (9)]; and the solution of the linear problem with the initial condition [Eq. (4)]. The result is presented in Fig. 4. We see that the initial phase modification compensates nonlinear effects and the solution of the NLSE is very close to the linear one. Note, that as expected, the solution of the NLSE with preimposed phase correction [Eq. (9)] preserves the Gaussian shape in the transversal direction.

The parameter $b$ in the phase is not well defined; we used the simple exponential model of the atmosphere, while the real one is more complex and depends on many parameters. The phase mask is design for specific intensity which can vary from shot to shot. This uncertainty we model with the help of parameter $\delta$.

Hence, it is very important to understand how much our results change, if the $b\to b+\delta b$ parameter in the preimposed phase will be different from the $b$ parameter in the acquired phase [Eq. (9)]. We did an extensive investigation of the sensitivity of the results to the value of $\delta b$. The results of the modeling are presented in Figs. 5 and 6. One can see that, even for high intensity pulse $P/P_{\mathrm{cr}}=5000$, a variation of $b$ of about 10% (both positive and negative) does not affect the smooth radial distribution in the focal plane and slightly reduces the peak intensity in the focal point.

The result of the incomplete compensation is a function of both $\delta b$ and $b$. In Fig. 7, we present the Strehl ratio and beam quality $M^2$ as a function of both parameters. We see that the suggested scheme is robust, and nonsensitive to $\delta b$, especially for smaller beam powers. The reason is that we do not need the complete phase compensation; it is sufficient to have the efficient $B$ integral $b$ small enough. We see for an example that even for $b=10$ variation of $b$ by 10% changes Strehl number by 5% only. We see also that detrimental effects of atmosphere reduces with laser elevation (smaller $b$).

The required phase pattern can be preimposed with the help of a phase mask. The spatial phase distribution is independent of the intensity, which simplifies the mask production. The amplitude of the phase variation (mask thickness) is proportional to the intensity. The above studies of $\delta b$ indicate that the results are not very sensitive to the phase amplitude. In practice, a few masks will be sufficient to cover the range of intensities, making the proposed method more practical. The above results can be useful for the material processing applications.

Femtosecond material processing has become more common due to the development of new, compact lasers, high precision processing, and small collateral damage. In many cases, the processing takes place in vacuum in order to eliminate laser breakdown in air. The laser light leaving the entrance window can be distorted due to self-focusing and can attain peak intensity away from the initial focal point. This focal point displacement can be important for three-dimensional memory writing or volume Bragg gratings. The phase distortion in the window distorts the beam, increasing the spot size—an effect detrimental for the precise machining.

Let us demonstrate that this situation is similar to that considered above. For fused silica, the critical power is about 1.5 MW. For typical fs processing (pulse energy 1 mJ and pulse duration 1 ps), the pulse power is about $1000P_{\mathrm{cr}}$ and the situation is similar to that discussed above. Thus we can directly use the results from the previous sections, for example to estimate the spot size and the Strehl ratio. The use of a phase mask to compensate the detrimental windows effects also looks promising.

We have demonstrated that the nonlinear effect of self-focusing in the atmosphere for space debris cleaning can be described with good accuracy within the TW model. Within this model, the nonlinearity produces phase front distortion, serving as a high aberration focusing lens. Optical phase distortion results in displacement of the focusing point and beam filamentation, degrading the system performance.

We demonstrated that the TW model describes beam focusing with strong aberrations, with field distributions in the focal plane far from Gaussian. The pattern of laser field is determined by a single dimensionless parameter $b$ similar to the $B$ integral used in laser design to control the self-phase modulation. The dependence on only a single parameter greatly simplifies the optimization of the beam prefocusing arrangements.

The description of linear propagation after exit from the atmosphere can be simplified using the fact that the field in the focal plane is proportional to the Fourier transform of the field exiting atmosphere. As a result, we obtained simple expressions for the peak intensity (Strehl ratio) and beam quality $M^2$ which can be calculated in terms of the exiting field. The use of these quantities as FOM provides a convenient quantative description of the beam structure on the target.

Because of the high accuracy of the TW model, one can compensate for self-focusing in the atmosphere by a preimposed phase distribution, which will cancel the nonlinear phase acquired during propagation in the atmosphere. Our modeling demonstrated that the detrimental effects of self-focusing even for $P/P_{\mathrm{cr}}\sim 5000$ can be almost completely eliminated by a preimposed phase calculated within the TW model.