The interaction of two parallel relativistic laser beams in underdense plasmas is investigated by considering the evolution of their wave envelopes. The energy transfer between the two lasers is given by an expression based on the evolution of the total laser power in a regime without beam mixing. It is shown that how the energy is transferred depends nonlinearly on the initial phase difference of the lasers, and the result of the interaction depends on the laser intensity, spot radius, and their separation distance. The results are verified by direct numerical solution of the relativistic nonlinear Schrödinger equations for the laser envelopes as well as particle-in-cell simulation. The study and results should be helpful for understanding the energy transfer behavior of multiple co-propagating laser beams in underdense plasmas.

## I. INTRODUCTION

Intense laser–plasma interactions (LPIs) have attracted much research attention because of their relevance to inertial confinement fusion (ICF),^{1} laser-driven table-top particle acceleration,^{2–4} etc.^{5} The ponderomotive force of an intense laser beam can expel the affected plasma electrons, leading to self-focusing and self-channeling, as well as filamentation of the laser^{6,7} in underdense plasmas. With two or more lasers, more complex nonlinear interaction features can appear.^{8–12} LPI involving two linearly polarized lasers that are orthogonally polarized has been investigated.^{13–16} Depending on the intensity, spot size, and relative phase of the lasers, scenarios, such as attraction, repulsion, and spiraling, of the latter can appear.

In most of the existing works,^{13–16} no account is taken of the initial phase difference, or only discussed cases where the laser beams are exactly in or out of phase. Recent work has shown that the initial phase difference plays an important role in the interaction of lasers,^{17} in particular, the energy transfer becomes intense for a small phase difference. It was noted that the previous case was mainly concentrated on the interaction behavior at the nonlinear stage for a long propagation distance, where the two laser beams are already strongly dissipated and mixed up due to the nonlinear response of plasma. This makes it difficult to calculate the ratio of energy transfer between the lasers. So far, a detailed theoretical analysis of the energy transfer behavior has still not been available. In addition, the length of low density plasma is sometimes only tens of micrometers, and the characteristics of long distance propagation are not yet obvious. Tens of micrometers are also the typical interaction length of two laser beams with beam spot of several micrometers crossing at a small incident angle.^{4} Thus, it is necessary to clarify the interaction dynamics of two lasers at a short propagation distance.

In this paper, we consider the interaction of two parallel circularly polarized laser beams propagating in underdense plasmas using the nonlinear Schrödinger equations (NSEs)^{18–24} for the envelopes of the lasers as well as by two-dimensional (2D) particle-in-cell (PIC) simulation. A relatively simple expression for the evolution of the total laser power in the half space is obtained. It describes how the energy is transferred in an initial interaction stage for a short propagation distance, where the laser energy is not seriously dissipated and the two beams are self-trapped without mixing up. In particular, it is found that there is already significant energy transfer in this regime, and the energy depends sinusoidally on the initial phase difference. The dependencies on the intensities, spot sizes, and separation distance of the lasers are also obtained, and the results are found to agree well with those obtained by numerically solving the NSEs as well as by direct PIC simulation.

## II. ANALYTICAL MODEL

Stationary propagation of two *circularly polarized* laser beams, say lasers 1 and 2, with their electric fields given by $12E1,2(x,y,z)(\xeax+i\xeay)ei(\u2212\omega t+kz+\varphi 1,2)+c.c.$, where *ê*_{x,y} are the unit vectors in the *x* and *y* directions and *ω* and *k* are the light wave frequency and wave vector, respectively, in uniform underdense plasmas along the *z* direction can be described by the coupled NSEs,^{13–24}

where *a*_{1,2}(*x*, *y*, *z*) = *eE*_{1,2}(*x*, *y*, *z*)/*m*_{e}*ωc* are the normalized complex electric fields *E*_{1,2}(*x*, *y*, *z*) of the lasers, *λ* = 2*π*/*k* is the wavelength of the laser light, $kp=n0e2/\u03f50mec2$ is 2*π* times the inverse plasma wavelength, *e*, *m*_{e}, and *n*_{0} are the electron charge, rest mass, and unperturbed plasma density, respectively, *c* and *ϵ*_{0} are the vacuum light speed and permittivity, $\gamma =1+|a|2$, where |*a*|^{2} = |*a*_{1} + *a*_{2}|^{2}, is the relativistic factor of the electrons moving in the total laser field, and $\u2207\u22a52=\u2202x2+\u2202y2$. The coupling of the two lasers is through the nonlinear electron density response to the relativistic ponderomotive force.^{18,20} Here, $n=ne/n0=1+kp\u22122\u2207\u22a521+|a|2$ describes the normalized density perturbations, where we have assumed that the electrons respond to the relativistic laser ponderomotive force and the ions act as homogeneous positive background. With this assumption, one sets *n* = 0 wherever *n* < 0 occurs,^{20} which can happen if electron cavitation (i.e., locally all electrons are driven out by the laser light) occurs. Equation (1) includes the effects of both the laser electric and magnetic fields and is valid for relativistic electron response to the latter, but wave damping and the nonlinear frequency shift $\Delta \omega =\omega 2\u2212k2c2\u2212\omega p2$, where *ω*_{p} = *k*_{p}*c* is the plasma frequency, are ignored.

The initial Gaussian laser beams are given by $a1,2(x,y,0)=a10,20\u2061exp\u2212(x\u2213d/2)2+y2w10,202+i\varphi 1,2$, where $a10,20$ are the initial amplitudes, $w10,20$ are the initial spot radii, *d* is the initial separation distance, and *ϕ*_{1,2} are the initial phases of the two lasers. For simplicity, in the following, we shall set $w10=w20=w0$. It is well known that the NSE without dissipation and source has a constant of motion, i.e., *∂*_{z}∫∫|*a*_{1,2}|^{2}*dxdy* = 0, where the integration is over the entire *x*, *y* space. This conservation law for propagating light wave envelopes is independent of the details of LPI. Thus, we have *∂*_{z}∫∫(|*a*_{1}|^{2} + |*a*_{2}|^{2})*dxdy* = 0, i.e., in the absence of dissipation and source, the total light power is conserved as the two lasers propagate and interact in the plasma. We can use this property to qualitatively deduce the behavior of the interaction by looking at the interaction of the lasers at the very beginning.

For simplicity, we introduce the half space laser powers *P*_{1,2} = ∫∫_{x>0,x<0}*Idxdy*, where *I* = *I*_{0}|*a*|^{2} is the total light intensity, *I*_{0} = 1.37(*λ* [*μ*m])^{−2} × 10^{18} W/cm^{2}, and $a=a1+a2ei(\varphi 2\u2212\varphi 1)$ describes the total electric field strength. The argument is based on the physically reasonable *ad hoc* assumption, to be verified *a posteriori* by numerical solution of the two NSEs (1) as well as by PIC simulation, that the initial evolution of the total power within a half space can already determine the subsequent propagation characteristics of the two lasers. Without the loss of generality, we set *ϕ*_{1} = 0 and *ϕ*_{2} = Δ*ϕ*, and *a* = *a*_{1} + *a*_{2}*e*^{iΔϕ} also satisfies Eq. (1). From Eq. (1), one finds that^{20}

Substituting the profile $a=a10\u2061exp\u2212(x\u2212d/2)2+y2w02+a20\u2061exp\u2212(x+d/2)2+y2w02+i\Delta \varphi $ and after some straightforward algebra, one obtains

where we have assumed that the profile of laser pulse changes slowly and the rate of energy transfer remains unchanged. This equation shows the effect of the laser parameters *d*, $w0$, Δ*ϕ*, etc., on the behavior of the power for *z* ≳ 0. When Δ*ϕ* = 0 or *π*, the right-hand side of Eq. (3) vanishes so that no energy transfer occurs. In particular, it is shown that the energy transfer depends nonlinearly on the initial phase difference, and for 0 < Δ*ϕ* ≤ *π*/2, the energy transfer increases with the phase difference, which is opposite to the previous conclusions for a later interaction stage.^{17} Such a discrepancy can be understood as follows. For a long propagation distance, Δ*ϕ* will increase and finally exceed *π*/2. As a result, the two laser beams repel each other, and the energy transfer process is gradually terminated. Thus, the interaction length becomes larger for a smaller Δ*ϕ* and the energy transfer increases as Δ*ϕ* decreases at a late interaction stage. However, for a short propagation distance considered here, the interaction length can be regarded as the same and the energy transfer is only dependent on the transfer ratio determined by Δ*ϕ*, and thus, the energy transfer increases as Δ*ϕ* increases.

It is shown that the above equation of the *ad hoc* model does not include plasma density since we have assumed that the laser profiles are unchanged and there is no dissipation in the initial stage of the propagation of the two laser beams in low density plasma. For two parallel relativistic laser beams with spot size of several to dozens of micrometers, propagating in underdense plasmas, the initial interaction stage can be estimated at the propagation distance *z* < 50 *μ*m and the propagation time *t* < 200 fs.^{17} In this case, the plasma effect can be neglected except that it provides the electrostatic charge-separation force to balance the laser ponderomotive force such that the two laser beams can propagate in plasmas without mixing up. For laser pulse undergoing intense profile modulation, which is induced by the plasma, for example, intense self-focusing, the laser parameters after modulation should be used. In addition, the model breaks down if the distance between the two lasers is smaller than $2ln\u20612w0$ or $d<2ln\u20612w0$. In this case, the peak field in the interference region exceeds that of a laser, and the two laser-driven channels can merge. In the following, we shall verify the above conclusions and assumptions by numerically solving Eq. (1) as well as by PIC simulation.

## III. SOLUTION OF THE NSEs

We now numerically solve the coupled NSEs (1) for *λ* = 1 *μ*m and the initial plasma density *n*_{0} = 0.1*n*_{c}, where *n*_{c} = *ω*^{2}*ϵ*_{0}*m*/*e*^{2} is the critical density. The initial laser profiles are $a1,2=12(\xeax+i\xeay)a10,20\u2061exp\u2212(x\u2213d2)2+y2w02exp(0,i\Delta \varphi )$.

In Figs. 1(a)–1(c), the two lasers are identical except for their initial position and phase difference Δ*ϕ*. These differences determine how the LPI starts since the beam with delayed phase will encounter a plasma that is slightly disturbed and that can affect the later propagation of both beams. Panel (a) shows that for Δ*ϕ* = *π*/2, the intensity of the upper (lower) laser decreases (increases). That is, energy is transferred between the lasers, but their separation distance remains unchanged. Panel (b) for Δ*ϕ* = 0 (no phase shift) shows that the two lasers merge into one. Panel (c) for Δ*ϕ* = 3*π*/4 shows that energy is transferred from the upper to the lower laser, and the two lasers repel each other. That is, the two lasers tend to attract each other if |Δ*ϕ*| < *π*/2 and repel each other if |Δ*ϕ*| > *π*/2.^{17} Panels (d)–(f) are for Δ*ϕ* = *π*/2. Self-focusing of the laser beams can be clearly observed. The self-focusing distance (say, from *z* = 0 to the small red region) of about 16 *μ*m in panel (f) is consistent with that of the well-known expression for thin laser beams (where the paraxial approximation is valid)^{25} $Zf=\pi w02/(\lambda Pin/Pcr\u22121)=12.9$ *µ*m, where *P*_{in} is the input laser power and $Pcr\u223c17.1\omega 2/\omega p2$ GW is the critical power for self-focusing. Figure 1(d) for $a10=2$ and $a20=22$ shows that laser 2 becomes more intense after propagating dozens of micrometers. Panel (e) for a large separation distance (*d* = 6 *μ*m) shows that the energy transfer is weaker, as expected. Panel (f) for $w0=4$ *μ*m and *d* = 8 *μ*m shows that lasers at large separation distance do not interact even if their spot radius is large. This can be expected since the overlapping region involves only very weak fields, especially after self-focusing of each of the lasers.

The results obtained here can be interpreted in terms of competition among the effects of the ponderomotive forces of the two laser pulses on the background plasma and its response. As can be seen in Fig. 2(a) for the laser intensity and electron density profiles along *x* at *z* = 6 *μ*m for Δ*ϕ* = *π*/2, the two lasers each expels the local electrons and creates two electron cavities, which in turn trap and guide the laser beams in the LPI. In fact, we see that both lasers can fully expel the plasma electrons in their center regions. The expelled electrons are compressed around the cavity edges, especially in the region (*x* ∼ 0) between them. Clearly, how the laser beams propagate depends on the balance of the light and plasma pressures in the overlapping region. However, since the electrons are expelled and compressed by the laser light, higher local electron pressure corresponds to higher light pressure, as can be deduced from the asymmetry of the electron density profile near *x* = ±0. Accordingly, how the lasers attract or repel each other can be qualitatively determined by examining the local forces in this region. In particular, the lasers will merge (repel) if the plasma pressure at *x* = ±0 is less (higher) than the light pressures. For example, in Fig. 1(b) for Δ*ϕ* = 0, the two identical laser beams merge since they both start at *x* = 0 with their weakest field region and the ponderomotively expelled electrons could not build up a sufficiently high pressured plasma wall between them to prevent their merging.

Figure 2(b) shows that the laser power *P*_{1} (in *x* > 0) depends monotonously on *z* for all phase differences Δ*ϕ* of the two lasers. It also suggests that the laser with delayed phase will continue to gain energy. For showing the extent of the energy transfer between the two lasers, it is convenient to introduce a parameter *η* for the averaged initial (*z* < 30 *μ*m, where the plasma effect has not significantly changed the beam profile) slopes of the *P*_{1}(*z*) curves in (c). Quantitatively, *η* is defined as the slope of the fitting line of *P*_{1} and *z* with the least square method. Figure 2(c) is for *η* vs the initial phase difference Δ*ϕ*. One can see that *η* behaves as sin(Δ*ϕ*), in agreement with Eq. (3) of our theory. However, the peak of *η* occurs at Δ*ϕ* slightly less than *π*/2, which can be attributed to the repulsion of lasers when |Δ*ϕ*| is very close to *π*/2, such that the energy transfer is terminated earlier than other cases with a phase difference smaller than *π*/2. For completeness, we show in Fig. 2(d) the phase angle arctan(*a*_{i}/*a*_{r}), where *a*_{r} and *a*_{i} are the real and imaginary parts of the total laser field *a*, at *z* = 6 *μ*m for Δ*ϕ* = *π*/2. One can see the mixing of the phases of the two lasers as they merge, and that in the laser overlapping region, the total laser field is of mixed polarization. Overall, the distribution of the total light intensity is modified by the coherent superposition of laser electric fields, and the self-focusing effect keeps each beam propagating in a limited region, which leads to the fact that the two laser beams keep transferring energy without mixing up.

Figure 3 shows the effect of the laser parameters on *P*_{1} (left column) and its average starting slope *η* (right column). One can see in panels (a), (c), and (e) that *P*_{1} always increases monotonically. That is, the light energy is transferred from the *x* < 0 region to the *x* > 0 region and the degree of transfer increases with *a*_{2}, as can also be seen in panel (b) for *η* vs *a*_{2}. In fact, *η* depends almost linearly on *a*_{2}, in agreement with Eq. (3) from the theory. Panel (c) shows that the slope of *P*_{1} decreases with an increase in *d* [but note the discussion on panel (e) below], and panel (d) for *η* vs *d* shows that *η* decreases exponentially for *d* > *w*_{0}, which is also in agreement with Eq. (3). Panel (e) for the dependence of *P*_{1} on *w*_{0} at fixed *d*/*w*_{0} (=2) shows that the slope of *P*_{1} does not depend on *w*_{0}, but *P*_{1} increases with *w*_{0}, also in agreement with Eq. (3). However, panel (f) for *η* vs *w*_{0} shows that although the magnitude is consistent with that from Eq. (3), the curve (blue asterisks) from the numerical solution of Eq. (1) is not a straight line as that from Eq. (3). This can be attributed to the fact that the laser beams with different spot radii would undergo different profile modulation, which slightly influences the interaction. In particular, the amount of energy transferred should be estimated by using the laser parameters after modulation if the laser beams would undergo intense and rapid profile modulation.

## IV. PIC SIMULATION RESULTS

In the above, we have shown that the results from our model Eq. (3) for the evolution of the LPI agree with that of the solution of the NSEs (1). In order to see the interaction scenario in more realistic situations as well as at longer propagation distances, we perform 2D PIC simulations of the two-laser interaction in plasma using the EPOCH code.^{26} The simulation box is 151 *μ*m × 60 *μ*m with 3020 × 1200 grids in the *x*, *z* plane (−1 *μ*m < z < 150 *μ*m, −30 *μ*m < x < 30 *μ*m), with 40 ions and 40 electrons per grid. The laser wavelength is 1.06 *μ*m, the plasma is located in *z* > 0, and its density is *n*_{e} = 0.1*n*_{c}, where *n*_{c} is the critical density. Two circularly polarized Gaussian laser beams are incident from the left boundary. Their electric fields are given by $a1=a10,20\u2061exp\u2212x\u2213d/22w02ei\varphi +(0,\Delta \varphi )$, where $a10,20$ are the peak electric fields, *d* is the separation distance between their axes, *w*_{0} is their spot radius, *ϕ* = −*iωt* + *ikz* is the phase of the light wave with frequency *ω* and wave vector *k*, and Δ*ϕ* is the phase difference between the lasers.

Figure 4 shows the laser intensity and plasma density at *t* = 50*T*_{0}. In panel (a), the parameters of the two lasers are the same except that the phase of the upper laser lags by *π*/2. We see that at larger *z*, the lower laser loses energy to the upper laser, which confirms the energy transfer scenario predicted by the NSEs (1) and the initial power evolution Eq. (3). However, at *z* > 30 *μ*m, the two beams start to repulse each other.^{4} This can be attributed to the change in the relative phase difference. The phase velocity of a laser beam propagating in underdense plasmas can be written as $vp=c/1\u2212ne/\gamma nc$, where $\gamma =1+|a|2$. Accordingly, as the upper laser gains energy, its phase velocity decreases and that of the lower laser increases. Thus, the relative phase difference increases and, when $\Delta \varphi >\pi 2$, repulsion occurs. The electrons in the plasma between the lasers are piled up in the region *x* ∼ 0, as can be seen in panel (b). Since the plasma refractive index decreases when the electron density increases, the laser light bends away from the latter region. Panel (c) shows the laser intensity for the case *I*_{2} = 2*I*_{1}. In panel (d), for the corresponding electron density, we see that the electrons between the two lasers are driven to the *x* > 0 region as the lasers enter the plasma, but are then driven to the *x* < 0 region after the lasers have propagated for about 10 *μ*m. This is due to the fact that the more intense lower laser loses energy to the upper one as they propagate. Panel (e) for *d* = 6 *μ*m and $w0$ = 2 *μ*m and panel (f) for *d* = 8 *μ*m and *w*_{0} = 4 *μ*m show that the interaction between the lasers is weak if the lasers are too far apart, even if the spot radius is larger. In particular, when the distance of two laser beams becomes larger than 4*w*_{0}, i.e., *d* ≥ 4*w*_{0}, the value of Eq. (3) is on the order of 10^{−3}, which indicates that the interference field becomes negligible compared to the field of each laser beam and in this case no interaction occurs. The light enhancement noticeable in both lasers shown in panel (f) can be attributed to self-focusing.

We now consider the validity of our *ad hoc* model [Eq. (3)] involving the energy of the two regions. It is convenient to introduce the parameter *μ* = *ε*_{1}/(*ε*_{1} + *ε*_{2}), where *ε*_{1,2} = (1/2)∫∫_{x>0,x<0}*ϵ*_{0}*E*^{2}*dxdz* is the electric field energy in the upper and lower half spaces. That is, *μ* is the ratio of the energies (of the two lasers) in *x* > 0. Assuming that dissipation is negligible, one can write

where $P10,20$ are the incident laser powers in the half spaces *x* > 0 and *x* < 0, and *∂*_{z}*P*_{1} is given by Eq. (3). We introduce the parameter *η*′ of the slope of the linear fitting of *μ*(*t*) for *t* < 30*T*_{0}. Figure 5(a) for the evolution of *μ* shows that *μ* increases linearly with *t*, in agreement with Eq. (4). Panel (b) for *η*′ vs Δ*ϕ* clearly shows a sin(Δ*ϕ*) dependence, also in agreement with the theory. Panel (c) for *η*′ vs *d* shows that *η*′ decreases exponentially with *d*. This result also agrees with that from Eq. (3), except when the laser spots are initially too close ($d<2ln\u20612w0$) together. Panel (d) for *η*′ vs *w*_{0} also shows agreement with Eq. (3) except when $w0$ is too large.

## V. SUMMARY

In summary, we have investigated the interaction of two co-propagating laser beams in underdense plasmas. It is shown that if the two lasers are initially not too close together, the interaction features and energy transfer characteristics can be described by a relatively simple expression, obtained from the coupled NSEs governing the two lasers and their interaction, for the evolution of the total laser power in the half space of (initially) one of the lasers. Based on this model, the ratio of energy transfer can be calculated for different laser intensities, spot sizes, separation distance of the lasers, and the phase difference, in an initial interaction stage. The predictions from the simple theory agree well with that from the solutions of the coupled NSEs as well as the corresponding 2D PIC simulations.

## ACKNOWLEDGMENTS

This work was supported by the National Key R&D Program of China, Grant No. 2016YFA0401100; the National Natural Science Foundation of China (NNSFC), Grant Nos. 11575031, 11875092, and 11705120; the National ICF Committee of China; the Featured Innovation Project of Educational Commission of Guangdong Province of China, Grant No. 2018KTSCX352; the Natural Science Foundation of Top Talent of SZTU, Grant No. 2019010801001; and the Open Fund of the State Key Laboratory of High-Field Laser Physics at SIOM-CAS. R. X. Bai would like to thank C. Z. Xiao, K. Jiang, T. Y. Long, R. Li, Y. C. Yang, C. N. Wu, and M. J. Jiang for useful discussions and help.