Generalized nonlinear Schrödinger equations describing the Second Harmonic Generation of femtosecond pulse, containing a few cycles, and their integrals of motion.

An interaction of laser pulse, containing a few cycles, with substance is a modern problem, attracting attention of many researches. The frequency conversion is a key problem for a generation of such pulses in various ranges of frequencies. Adequate description of such pulse interaction with a medium is based on a slowly evolving wave approximation (SEWA), which has been proposed earlier for a description of propagation of the laser pulse, containing a few cycles, in a medium with cubic nonlinear response. Despite widely applicability of the frequency conversion for various nonlinear optics problems solutions, SEWA has not been applied and developed for a theoretical investigation of the frequency doubling process until present time. In this study the set of generalized nonlinear Schrödinger equations describing a second harmonic generation of the super-short femtosecond pulse is derived. The equations set contains terms, describing the pulses self-steepening, and the second order dispersion (SOD) of the pulse, a diffraction of the beam as well as mixed derivatives. We propose the transform of the equations set to a type, which does not contain both the mixed derivatives and time derivatives of the nonlinear terms. This transform allows us to derive the integrals of motion of the problem: energy, spectral invariants and Hamiltonian. We show the existence of two specific frequencies (singularities in the Fourier space) inherent to the problem. They may cause an appearance of non-physical absolute instability of the problem solution if the spectral invariants are not taken into account. Moreover, we claim that the energy preservation at the laser pulses propagation may not occur if these invariants do not preserve. Developed conservation laws, in particular, have to be used for developing of the conservative finite-difference schemes, preserving the conservation laws difference analogues, and for developing of adequate theory of the modulation instability of the laser pulses, containing a few cycles.

Introduction

A well-known process of second harmonic generation (SHG) by the laser beam was the first nonlinear optical effect, which was observed by P. Franken et al. [1] in 1961. One year later, J. Armstrong, N. Bloembergen et al. [2] published the fundamental paper, dealing with the optical frequency conversion, in which the coupled nonlinear equations were developed and many explicit solutions of these equations were obtained in the framework of the plane wave approximation. During passed decades the SHG was widely observed: in plasma optics [3-5] and nonlinear optics, including high-harmonic generation [6-12], and in a medium containing nanoparticles [13, 14], and in semiconductor [15]. In [16] the SHG at the boundary between dielectric media is analyzed. Obviously, this list of papers can be easily increased.

Many researchers pay an attention to an investigation of various factors, which restrict the frequency doubling efficiency. For example, in [17] the effect of a group-velocity dispersion was considered. In [18] the authors have considered an influence of interplaying between two types polarizations of the pulse at the fundamental frequency on the frequency doubling efficiency. With increasing of the laser pulse intensity it is necessary to take into account a cubic nonlinear response of a medium under the SHG. An influence of the self-modulation as well as cross-modulation of the laser pulse on the frequency doubling efficiency was investigated experimentally and theoretically (see, for example [2, 19–28]). It should be stressed that the most general analytical theory of this process has been developed recently in [28], taking into account the time-dependent pulse shape and without using the basic wave energy non-depletion approximation on the base of using the integrals of motion (conservation laws) of the problem. We notice that the authors of [22] have obtained the most (but not all) solutions of this problem in the framework of long pulse duration approximation. In [29] an effective way for increasing of SHG efficiency in bulk medium was demonstrated by using of the incident beam with a ring (tubular) profile. This leads to re-profiling of the beams due to their diffraction and, therefore, to decreasing the beam phase distortions.

New opportunities appeared at SHG under the big phase-mismatching, so-called cascading SHG [30-33]. This process allows us to achieve an effective cubic nonlinear response in a medium with a quadratic nonlinear response. This effect was firstly predicted by Yu. N. Karamzin and A. P. Sukhorukov [34] in 1974. Later, it was used for a compression of the soliton with femtosecond duration [30, 34–50], and for the beam self-focusing [51-53], and for using the cascading SHG for suppression of the optical pulse intensity fluctuations, occurring in a sequence of pulses with randomly variations of their maximal intensities [54-57]. Such sequence of the laser pulses is produced by a laser, operating at the free generation mode. The total duration of the pulses sequences is equal to milliseconds. Therefore, the cascading SHG allows realizing the laser generation mode that is similar to the Kerr-locked mode applying for the femtosecond laser system.

There are a lot of other schemes of SHG, using the modern substances (see, for example, [58-84]) and at a present time the problem is still actual and important for various applications. Among these applications, we stress using of SHG for a measurement of pulse parameters and of parameters of a medium as well as for visualization of various fast processes, occurring in the substances. The most famous approaches of parameters measuring are SPIDER (spectral interferometry for direct electric-field reconstruction) [85] and FROG (frequency-resolved optical gating) [86, 87]. Also, in [88], a technique, based on combination of FROG spectra measurements was proposed to completely characterize the amplitude and phase of an ultrashort pulse in space and time. Mix FROG method was implemented in [89]. Thus, the frequency doubling of the optical pulse, containing a few cycles, is of practical interest until present time. However, in this case it is necessary to take into account the self-steepening of the laser pulses on both fundamental and doubled frequencies. According to our knowledge, the first analysis of this problem was made numerically in [90] at the group velocity synchronism and phase matching, without taking into account the beam diffraction and the second order dispersion influence on the pulses interaction. For writing the corresponding equations, the authors have applied an approach, which was developed first by N. Tzoar and M. Jain [91] under the investigation of the single pulse propagation in a medium with a cubic nonlinear response. Experimental observation of the self-steepening influence on the pulse spectrum deformation was demonstrated in [92]. Then, this approach was successfully used for describing the processes of a laser pulse propagation in the optical fiber [93]. The corresponding nonlinear Schrödinger equation (NLSE) contains the first time derivative of a nonlinear medium response. This equation was derived using the slowly varying envelope approximation (SVEA) and named as generalized NLSE.

Taking into account a laser beam diffraction, T. Brabec and F. Krausz proposed in 1997 a new approach, so-called slowly evolving wave approximation (SEWA). Derived equation was also named as the generalized nonlinear Schrödinger equation (GNLSE) [94]. Let us note that the SEWA requires not only an envelope of a wave packet, but also the phase changing must slowly vary as the pulse covers a distance, which is equal to the wavelength. The main feature of GNLSE consists in its containing of mixed derivatives on time and spatial coordinates, and of the second-order time derivative of a medium nonlinear response (which describes the pulse self-steepening).

It should be emphasized that a differential operator of the GNLSE, written in [91, 93], coincides with a differential operator of the GNLSE, derived in [94] if the optical beam diffraction does not take into account. They differ only by the factor two at the coefficient characterizing the self-steepening term of the equation. Obviously, these equations written in the dimensionless units look the same.

Finding of conservation laws attracts an attention of various authors (see, for example, recently published papers [95-99]). In papers [100, 101] we derived the conservation laws for the frequency doubling process if the self-steepening of the pulses occurs. However, until present time, the set of GNLSEs describing a SHG of super-short femtosecond pulse with taking into account a diffraction of the laser beam (it means that this process is described in the framework of the SEWA) as well as their conservation laws were not derived. We do this in the present paper.

As is well-known, the invariants are very useful for computer simulation of the laser pulse propagation because the equations describing various problems of nonlinear optics are nonlinear ones and as a rule they do not allow us to find their analytical solution. As a consequence, the finite-difference time-domain (FDTD) method [102, 103] for computer simulation of a propagation of the optical pulse, containing a few cycles, is usually employed. For example, in [104] the authors develop the time-transformation method based on FDTD. An alternative method, which permits applying fast numerical algorithms, is a generalized source method (GSM) [105]. Linear GSM was adapted to describe the SHG in diffraction gratings, containing non-centrosymmetric materials [106]. Obviously, there are other papers, in which a computer simulation is applied for an investigation of the SHG problem. As a rule, a split-step method was used at computer simulation.

Another approach for computer simulation is based on developing of the conservative finite-difference schemes [107]. They allow preserving the difference analogues of the corresponding conservation laws. Obviously, it is necessary firstly to write these conservation laws. With this aim, we propose a new transform for the derived GNLSEs, which reduces the equations to the form containing neither mixed derivatives on time and spatial coordinates nor time derivatives of the nonlinear response of a medium. Using a new set of GNLSEs, we obtain the energy invariant, spectral invariants and Hamiltonian. We add some specific requirements to the problem statement (excluding of the singularities in frequency space) to avoid a development of non-physical instability of the SHG process at taking into account the pulses self-steepening. These requirements, together with the spectral invariants, are also important for a validity of the energy conservation law.

Generalized nonlinear Schrödinger equations derivation

Derivation of the equations for SHG process in the framework of SEWA

Derivation of the equations describing the optical frequency doubling in the framework of SEWA for an electric field strength E(z, x, t) starts from the well-known wave equation: written above in 2D case, for example. The laser pulse propagates along z coordinate, t is a time, c is the light velocity in a vacuum, x is a transverse coordinate to the optical pulse propagation direction, an operator is the transversal Laplace operator. We choose only x-coordinate as a transverse one for simplicity. This does not restrict our consideration. The function D(z, x, t) describes an electric field induction and is written in the form: the function P(z, x, t) is a medium polarization, containing its linear and nonlinear responses:

For linear non-instantaneous medium response P(z, x, t) the following representation is widely used, which yields in the relation for a dielectric permittivity: here is a linear electric susceptibility of a medium at the frequency ω.

Function P(z, x, t) describes a nonlinear response of a medium. If we consider a medium with the quadratic nonlinearity, then this function is defined as P = χ(2)E2, and χ(2) is a quadratic electric susceptibility.

The next step for the wave equation reduction consists in representation of an electric field strength E(z, x, t) in a following manner: where E1(z, x, t), E2(z, x, t) are the electric field strengths of waves, propagating with carrier frequency ω, and with doubled frequency 2ω and with wave-numbers k1, k2, correspondingly. In (6) we introduce the slowly varying envelopes for both wave packets and denote them as A1(z, x, t) and A2(z, x, t). Symbols c.c. denote a complex conjugation. The linear polarization of a medium at these frequencies can be written as

The corresponding nonlinear polarization of a medium at chosen frequencies is written as follows

Let us remind that the following relations for wave-numbers are valid: where n, n2 are the medium refractive indexes at the frequencies ω and 2ω, correspondingly.

Taking into account the expressions (6) and (8) and (1) transforms to the following form: where D(z, x, t), j = 1, 2 is a linear part of the electric field induction at the frequencies ω and 2ω, correspondingly:

Let us multiply (9) by exp(− ik1z) or exp(− ik2z). Then, we integrate each of these expressions over the corresponding wave period, taking into account the orthogonal property of the sine and cosine functions. Then, (9) can be reduced to the following equations:

Let us note the Fourier transform property:

Therefore, (11) can be written as:

Taking into account the dispersion relation: one can write (13) as follows

Normally, we expand k(ω) in a series near the frequency ω in (15): and near the frequency 2ω in (16): here β and s are the n–th or m–th order derivatives of the wave number at the frequencies ω and 2ω, correspondingly.

Let us neglect the terms β, s, n, m = 3, 4, … referring to the third order dispersion and other higher terms in (17) and (18). Also, we neglect the term (ω − ω)4. Thus, (17) and (18) take the forms:

After Inverse Fourier transform, we obtain the following equations set for describing the frequency doubling process:

Above Δk = k2 − 2k1 is the phase mismatch. (21) and (22) are called as the set of generalized nonlinear Schrödinger equations (GNLSEs).

Coordinate system transform

First type of coordinate system transform

In a coordinate system moving with the first wave packet (21) and (22) could be re-written as follows:

Let us raise to the second power the expressions in brackets and then neglect the terms with the second order derivative on ξ − coordinate:

By grouping the terms in (26) and (27) and dividing them by 2ik, j = 1, 2, we obtain:

The set of Eqs (28) and (29) describes the SHG process for the pulse with super-short duration, which is about a few femtoseconds.

Further reduction of these equations is a consequence of the following assumptions.

The group velocities of pulses and their phase velocities differ insignificantly:

Difference between inverse values of the phase velocities of waves is many times less than unity:

As a result, (28) written with respect to the wave for basic frequency ω, reduces to and (29) takes the following form:

Second type of coordinate system transform

In a coordinate system moving with overage velocity the set of Eqs (21) and (22) could be re-written as follows:

If the assumptions (30)–(32) are valid, then (36), with respect to basic frequency ω, reduces to and (37), with respect to doubled frequency 2ω. takes the following form:

Below we use the First type of Coordinate system transform.

Dimensionless variables and problem statement

Let us introduce dimensionless variables:

Above the parameter τ is a duration of the incident pulse for basic frequency ω, I is its maximal intensity, the parameter a is a beam radius.

In new variables, the GNLSE set (33) and (34) takes the following form: which is considered in the following domain with the boundary conditions (BCs) for the functions A(ξ, x, t), j = 1, 2: or and with the initial conditions

Above γ is a parameter, which is inversely proportional to the doubled frequency 2ω and the pulse duration τ:

The parameter ν describes group-velocity mismatch (GVM) normalized by pulse duration τ:

Parameters D21, D22 are equal to a ratio between the diffraction length of the beam and a dispersion length of the pulse for each wave:

They characterize the second order dispersion of a wave packet. Parameter D is equal to unity for chosen normalization of spatial coordinate along which the laser pulse propagates. Parameter α describes the nonlinear coupling of waves and is defined in the following form which is written using the following assumption: which does not restrict our consideration. Let us note that if the equality (51) is not valid, then one can transform the equations with respect to new variables, in which we obtain only single coefficient (it characterizes the nonlinear coupling of waves for both equations, see [21]).

Usually, at theoretical analysis of the problem (42), the exponentially decaying functions A(ξ, x, t), j = 1, 2 are used as the initial laser pulse distributions: and instead of the BCs (43) and (44) the set of equations (42) is solved in the following unbounded (increasing of |x|, |t| and z) domain

Below, for definiteness, we use BCs (43). In addition, for convenience, below we use the previous notation for a longitudinal coordinate: z instead of ξ.

Equations transform

For transforming the equations set (42) to the form, which does not contain a time derivative of the nonlinear response of medium and mixed derivatives on time and spatial coordinate, we use the following equalities:

Also, due to a finiteness of the initial distributions of complex amplitudes and boundedness of the laser pulse propagation distance, let us state the following additional conditions:

Using the transforms (53) and condition (54), GNLSE set (42) can be written as follows (see detailed derivation in Appendix [A]):

It should be stressed that these equations will be also used for writing the conservation laws.

For further reduction of the equations, let us introduce the functions: or

Keeping in mind the relations (53) and (58), the set of Eqs (55) and (56) can be written as follows

Further, let us introduce the new functions E(z, x, t), j = 1, 2 in accordance with a rule: and the functions , j = 1, 2 as

Then, (59) and (60) can be written in the following manner:

These equations should be solved together with the relaxation equations with respect to the functions E(z, x, t), , A(z, x, t):

Taking into account the BCs (43) or (43′), and initial conditions (45) for the functions A, j = 1, 2, we obtain the following BCs with respect to the functions E, on time coordinate:

As consequence of BCs (43), we obtain the following conditions: if we solve the problem with BCs (43′).

Obviously, the BCs for these functions on x-coordinate are

The initial conditions for the functions E, , j = 1, 2 are

All conditions written above will be used at the construction of the conservation laws.

Problem invariants

Using new variables, some conservation laws can be derived.

The conservation law of energy

Theorem 1. The problem (42) with BCs (43) and initial conditions (45) possesses the conservation law of energy

Proof. Let us suppose that the functions A satisfy the conditions (43) and (54). Then, we multiply (55) and (56) by , and the equations, conjugated to (55) and (56),—by A. Then, we sum the obtained equations and integrate the resulted equation with respect to x, t coordinates in the domain Ω. Thus, we obtain the following equation in Appendix [B]:

For the validity of the invariant (73) one has to show that the last integrals in (74) are equal to zero:

In fact, from (55) and (56) at the time moment t = L and from the BCs (54), it follows that

This yields to the relations

It means that the Fourier harmonic amplitudes at the frequencies and vary with x-coordinate growing. This property occurs even if we consider the laser pulse linear propagation (in (42) the parameter α equals zero). Obviously, the physical reason of such harmonic amplitude evolution is absent. Therefore, the functions C1(z), C2(z) have to be equal zero. If C1(z) will be equal to nonzero constant then the energy of the laser beam will increase with a propagation distance. Obviously, a physical reason of this is absent for the case under consideration. Consequently, the amplitudes of Fourier harmonics at the frequencies and must be equal zero:

Thus, the invariant (73) is valid. In a case of the laser pulse nonlinear interaction (α ≠ 0), the same conditions take place.

Corollary 1. The Fourier harmonics at the frequencies and must be absent in the incident pulses spectra the energy conservation law preservation for the set of equations (42):

Let us note that if the parameter γ tends to zero, then the set of equations (42) reduces to the set of the NLSEs and the energy invariant (73) is valid. If the spectral harmonics or are absent in the incident pulse spectrum, then these spectral harmonics will be absent in the pulse spectrum computed in any section of a medium. This statement validity is stated by the spectral invariants, which are formulated below.

Spectral invariants

Spectral invariants describe an evolution of spectral harmonic amplitudes at the frequencies and along z-coordinate and show that the amplitudes should not increase during the laser pulses interaction. For writing of the Spectral invariants, let us introduce the following additional conditions:

Actually, in laser physics the propagation distance is boundedness and the domain under consideration can be chosen in such a way that the additional conditions are valid in this domain. We note that these conditions are important only for the invariant obtaining and they are not mandatory for the problem statement. The problem could be solved without additional conditions. Obviously, instead of BCs (80) one can use BCs (52). We stress that the laser sources generate the finite beam distribution in spatial coordinates. Therefore, the statements mentioned above are valid with respect to the domain choice.

Theorem 2. The problem (42) with BCs (43) and additional conditions (78) and (80), and initial conditions (45) possesses the spectral invariants

Proof. Let us consider (63) at a time moment t = L and integrate them with respect to x-coordinate:

Taking into account the relaxation Eqs (64) and (65) at time moment t = L and additional conditions (78) we obtain:

It is easy to see that the third terms in (83) and (84) are equal to zero:

The last integral (88) equals zero due to BCs (43). Thus, we obtain the following problems: with the initial conditions

Therefore, the invariants (81) and (82) can be obtained as the solutions of the problems (89)–(92).

Hamiltonian of the problem

Let us use the following substitution for the function A2(z, x, t):

Then, the functions E2(z, x, t), can be replaced as follows:

Using the substitution (93)– (95), the set of equations (42) can be written in the following form: which does not contain the terms with exp(±iΔkz). Consequently, with respect to the substitutions (94) and (95), the set of equations (96) takes the form:

Theorem 3. The problem (42) with BCs (43) and additional conditions (78), and initial conditions (45) possesses a Hamiltonian

Proof. Let us multiply (63) by and the equation, conjugated to it,—by . Then, we sum the obtained equations and integrate the resulting expression with respect to x, t coordinates in the domain Ω under consideration. Then, we write an equation in the form of the sum of integrals:

Analysis of the integrals is made in Appendix [C]. On its base, one can obtain the Hamiltonian of the problem under consideration.

Corollary 2. Using the inverse transforms for (93)–(95) the Hamiltonian takes the following form:

Discussion

A few words about applicable range of the laser pulse parameters and a medium at which it is necessary to use the GNLSEs for a description of the process under consideration. It should be noticed that the first experiment, which demonstrated self-steepening of one pulse, was made in paper [91] for the pulse with picosecond duration propagating in optical fiber about 5 km long. The pulse, propagated this distance, possessed non-symmetrical spectrum. Therefore, a key role for application of GNLSEs plays a relation between the dispersion length of the pulse and the laser pulse propagation distance. However, the pulse self-steepening appears at short length (about a few centimeters) of a medium if the pulse duration is about 10-50 fs.

Another case of using GNLSE for the laser pulse propagation description corresponds to falling of incident laser pulse with sharp intensity distribution on a medium because an influence of a time derivative of the nonlinear response of the medium enhances many times in comparison with the Gaussian pulse propagation. However, contrast temporal nonlinear response can be induced by the laser pulse if an optical bistability occurs. As is well-known, in this case the explosive changing of the characteristics of a medium occurs. As a result, the temporal structure with strong gradient appears. Therefore, a time derivative of the nonlinear response of a medium increases its influence on the phase modulation.

A detail analysis of SHG describing in the framework of the GNLSEs has to be made using computer simulation. Nevertheless, one can suppose some characteristics features of the frequency conversion in such conditions. First, the group velocity of each of wave packets will depend on complex amplitude of another wave packet. Second, obviously that the pulse spectra will distort and become non-symmetrical. Third, the nonlinear length, that defines the distance at which the energy of base wave transfers to the energy of the wave with doubled frequency, will differ for front of the pulse and its trailing part. Fourth, because of the presence of mixed derivatives, the optical beam diffraction will depend on time moment of propagating pulse.

Conclusions

In the framework of the SEWA approximation, we derived the GNLSEs, describing the SHG in a medium with a quadratic nonlinear response for the pulses, containing a few cycles. The main feature of the equations set concludes in a presence of the second order time derivative of the nonlinear response of a medium (dispersion of a medium nonlinear response) as well as the mixed derivatives on time and spatial coordinate.

We proposed the equations transform, which reduced these equations to the other ones, containing neither mixed derivatives of complex amplitudes nor time derivatives of the nonlinear response. These equations are more convenient for the computer simulation and for theoretical analysis of the frequency doubling process of the optical pulses, containing a few cycles.

Based on this transform we derived some conservation laws (invariants) for the SHG problem. We showed an existence of two specific frequencies (singularities in the Fourier space) inherent to the problem. They may cause an appearance of non-physical absolute instability of the problem solution if the spectral invariants are not taken into account. It should be mentioned that a presence of such singularities was also discussed in [32] at analysis of the modulation instability occurring in χ(2)-medium. We showed that the energy conservation law is valid at certain conditions on the incident pulse spectra: the spectra must not contain non-zero spectral amplitudes at two specific frequencies. Their zero-value amplitudes at the pulse propagation are provided by the spectral invariants. We derived also the Hamiltonian (the third invariant) of the problem.

All invariants mentioned above should be taken into account at least for developing of the conservative finite-difference schemes at computer simulation of the problem under consideration.

Appendix

A Equation transform derivation

Taking into account the transforms (53), one can write – for the first equation of the set (42), and – for the second equation of the set (42).

Then, multiplying (101) by , and (102) by , and integrating them with respect to t-coordinate, we obtain: correspondingly.

Due to a boundedness of the initial distribution and boundedness of the laser pulse propagation distance, let us state the following additional type of the conditions:

Then, (103) and (104) become the following ones: which contain the first derivative of the quadratic nonlinearity. Multiplying (106) by , and (107) by , we obtain the following equations:

Let us stress that these equations are used also for writing of conservation laws.

B Energy invariant derivation

Using (55) and (56) we obtain the equations for the functions A and , j = 1, 2, which can be re-written as a sum of integrals:

Below we discuss transforms of the integrals mentioned above. The first one contains the exact differential of functions A1, A2 intensities. Therefore, it is easy to see that the integral I1 transforms to the following expression:

Also, the integral I2 contains the time derivative of function A2 intensity: which equals zero, due to BCs (43). As a consequence of the integration by parts, the equalities take place under the conditions (43) validity.

Let us show, that the integral of the nonlinear terms I5 equals zero. The terms in integral I5 without the parameter γ cancel each other and integral I5 can be written as

For the integrals I41, I42, taking into account the functions P, j = 1, 2 definitions (57) and the BC (43), we obtain:

Thus, (110) transforms to the following equation:

C Hamiltonian derivation

As mentioned in (99), we obtain a sum of the integrals. Let us transform them. First, taking into account (64) and (78) and integrating of the I1 by parts we obtain:

The integrals I2 are transformed into the following expression:

It is easy to see that the integrals J21 and J22 are transformed into the following integrals:

Thus, the integrals I2 take the form:

The integrals I3 can be written in the following way: where and

Let us note that multiplying equation, conjugated to (63), by , and taking into account the BCs (43) and additional conditions (78), we obtain, that the second term in (123) equals zero. Thus,

The integrals J33 can be written as follows

Also, the integral I4, by using (64) and integrating by parts, reduces to the following one:

The integrals in I5 can be written as

Using the condition (78) and integrating by parts, one can transform the integrals I51, I52 to the form and

Obviously, the integral I6 can be written as follows:

Substituting the expressions (116), (120), (121), (126), (127) and (130) into (99), we obtain the invariant (98).