1. Introduction
Four-wave mixing (FWM) is a fundamental nonlinear optical process, in which four photons interact by means of the third order nonlinear susceptibility of the medium. FWM can be observed in a wide range of materials, including optical fibers [
1
G. P. Agrawal, Nonlinear Fiber Optics , (Academic Press, San Diego, 2001).
]. Modulational instability (MI), also known as Benjamin-Feir instability [
2
T. B. Benjamin, K. Hasselmann, T. B. Benjamin, and J. E. Feir, “The disintegration of wave trains on deep water. Part 1. Theory,” J. Fluid Mech.
27, 417 (1967). [CrossRef]
], is a special kind of FWM instability, and it is typical of nonlinear dispersive systems. It is known to exist in different branches of physics such as fluid-dynamics, plasma physics, solid-state physics and nonlinear optics [
2
T. B. Benjamin, K. Hasselmann, T. B. Benjamin, and J. E. Feir, “The disintegration of wave trains on deep water. Part 1. Theory,” J. Fluid Mech.
27, 417 (1967). [CrossRef]
,
3
T. Taniuti and H. Washimi, “Self-Trapping and Instability of Hydromagnetic Waves Along the Magnetic Field in a Cold Plasma,” Phys. Rev. Lett.
21, 209 (1968). [CrossRef]
,
4
A. Hasegawa, Plasma Instabilities and Nonlinear Effects (Springer-Verlag, Heidelberg, 1975). [CrossRef]
]. It manifests itself as an exponential growth of two spectrally symmetric sidebands when launching a continuous wave (CW), by virtue of the interplay between Kerr nonlinearity and group velocity dispersion (GVD). In optical fibers it leads to the breakup of a CW or quasi-CW beam into a train of ultrashort pulses [
5
G. Millot, E. Seve, S. Wabnitz, and M. Haelterman, “Observation of induced modulational polarization instabilities and pulse-train generation in the normal-dispersion regime of a birefringent optical fiber,” J. Opt. Soc. Am. B , 1266 (1998).
], and in the spatial case to the formation of patterns [
6
J. Klinger, H. Martin, and Z. Chen, “Experiments on induced modulational instability of an incoherent optical beam,” Opt. Lett.
26, 271 (2001). [CrossRef]
]. The long-term dynamics of MI can be very complicated, leading for instance to recurrence phenomena [
7
E. Fermi, J. Pasta, H. C. Ulam, M. Tsingou, G. Van Symaeys, Ph. Emplit, and M. Haelterman, “Experimental Demonstration of the Fermi-Pasta-Ulam Recurrence in a Modulationally Unstable Optical Wave,” Phys. Rev. Lett.
87, 033902 (2001). [CrossRef]
], and it is strictly connected to the existence of solitary waves [
8
N. J. Zabusky and M. D. Kruskal, “Interaction of ‘Solitons’ in a Collisionless Plasma and the Recurrence of Initial States,” Phys. Rev. Lett.
15, 240 (1965). [CrossRef]
].
In the context of nonlinear optics MI can be interpreted as
degenerate FWM: two pump photons of identical frequency
ω are converted into two new photons of different frequencies
ω –
δ and
ω +
δ, symmetric with respect to
ω, called respectively
Stokes and
anti-Stokes photons [
1
G. P. Agrawal, Nonlinear Fiber Optics , (Academic Press, San Diego, 2001).
]. The existence and efficiency of frequency conversion by FWM crucially depends on the so-called
phase-matching conditions [
1
G. P. Agrawal, Nonlinear Fiber Optics , (Academic Press, San Diego, 2001).
], which are associated to the momentum conservation of the four photons involved in the process. Such conditions are determined in large part by the linear characteristics of the system, i.e. its GVD, even though the entire process is nonlinear in its nature.
Traditionally, FWM instabilities are analyzed in nonlinear optics in the framework of the nonlinear Schrödinger equation (NLSE) [
1
G. P. Agrawal, Nonlinear Fiber Optics , (Academic Press, San Diego, 2001).
]. This equation (and its various extensions) provides the simplest, and in many cases even an accurate model for predicting the unstable frequency regions. In this work we show that, although NLSE (and its various generalizations, known as GNLSE, see [
1
G. P. Agrawal, Nonlinear Fiber Optics , (Academic Press, San Diego, 2001).
,
9
Y. S. Kivshar and G. P. Agrawal, Optical Solitons , (Academic Press, San Diego, 2003).
]) works extremely well in the context of glass optical fibers [
10
F. Biancalana, D. V. Skryabin, and P. St. J. Russell, “Four-wave mixing instabilities in photonic-crystal and tapered fibers,” Phys. Rev. E
68, 046603 (2003). [CrossRef]
], there are physical situations in which, in addition to the term containing the second order derivative in time (which models the temporal dispersion), it is necessary to include a spatial dispersion coefficient, proportional to the Laplacian operator. This ‘spatiotemporal’ dispersion may be due, for example, to materials that exhibit the presence of elementary matter excitations, such as
excitons [
11
F. Bassani and G. Pastori Parravicini, Electronic States and Optical Transitions in Solids, (Pergamon Press, Oxford, 1976).
]. The transverse components of the spatial dispersion term due to excitons contribute to the diffraction of the light beam propagating in the medium. This diffraction can be compensated for example by using common semiconductor waveguides [
14
T. Pertsch, T. Zentgraf, U. Peschel, A. Brauer, and F. Lederer, “Anomalous Refraction and Diffraction in Discrete Optical Systems,” Phys. Rev. Lett.
88, 093901 (2002). [CrossRef] [PubMed]
]. However, the longitudinal component of the spatial dispersion cannot be compensated by the waveguide, and it introduces a new element in the dynamics of CWs and pulses, that, as we show here, strongly affects the nonlinear optical instabilities.
The structure of the paper is as follows. In section 2 we propose a propagation equation, by extending the (1+1)-NLSE to include both spatial and temporal dispersion coefficients, and we analyze the linear and nonlinear behavior of carrier plane waves in detail. Our equation must be distinguished from the non-paraxial NLSE [
12
P. Chamorro-Posada, G. S. McDonald, and G. H. C. New, “Non-Paraxial Solitons,” J. Mod. Opt.
45, 1111 (1998). [CrossRef]
,
13
J. M. Christian, G. S. McDonald, P. Chamorro-Posada, J. M. Christian, G. S. McDonald, and P. Chamorro-Posada, “Bistable Helmholtz Solitons in Cubic-Quintic Materials,” Phys. Rev. A
76, 033833 (2007). [CrossRef]
], due to a crucial difference in the sign of the spatial dispersion coefficient, which provides a totally different physical interpretation of the equation. Appendix A gives a physical justification of
Eq. (1), derived in a system pumped in the proximity of an exciton resonance in a specially designed superlattice structure described in the literature. In section 3 we analyze the nonlinear parametric instabilities of the master equation, and we find some rather unconventional and unexpected properties. In section 4 we show the existence of an analytical family of bright and dark soliton solutions of the master equation, and we study their existence regions in the parameter space. Conclusions and final remarks are given in section 5.
2. Nonlinear evolution equation
In the present work we study the following evolution equation:
Eq. (1) is a NLSE with a spatially dispersive term proportional to the coefficient
D, which acts
along the longitudinal coordinate z. Parameter
s = ±1 depending on whether we have anomalous or normal temporal dispersion, respectively.
ψ is a field proportional to the slowly varying amplitude of the electric field. The last term models nonlinear cubic (Kerr) nonlinearity. We wish not to distract the reader from the main purpose of this work, thus the full physical justification of the model described by
Eq. (1) shall be given in the Appendix A, where we study a waveguide of semiconductor material populated by exciton-polaritons possessing a
negative effective mass, when optically excited slightly off-resonance.
Eq. (1) is formally equivalent to the equation for non-paraxial solitons [
12
P. Chamorro-Posada, G. S. McDonald, and G. H. C. New, “Non-Paraxial Solitons,” J. Mod. Opt.
45, 1111 (1998). [CrossRef]
,
13
J. M. Christian, G. S. McDonald, P. Chamorro-Posada, J. M. Christian, G. S. McDonald, and P. Chamorro-Posada, “Bistable Helmholtz Solitons in Cubic-Quintic Materials,” Phys. Rev. A
76, 033833 (2007). [CrossRef]
], with the fundamental difference that parameter
D can assume both signs, and that this parameter does not depend on the transverse profile of the field, see also discussions in Appendix A.
To explore the parametric instabilities of
Eq. (1), the knowledge of the stationary states of the system is necessary. For plane waves, this is accomplished in the Fourier domain by the expansion
ψ =
ψ
0
e
ikz
-
iωt
, where
ψ
0 is a real constant field. This gives a dispersion relation 𝓓(
k,
ω,
I) = 0, where
I≡
ψ
2
0 is the dimensionless electric field intensity. In scattering problems usually one wants to keep
ω real by solving the dispersion relation for complex
k.
For D = 0, 𝓓 reduces to the well-known NLSE dispersion: 𝓓(k,ω,I)=-k-sω
2/2+I = 0, which has only one solution in terms of k = k(ω). Thus, only one steady state is physically excitable for NLSE. In the case D ≠ 0, we have 𝓓(k,ω,I)=-k+Dk
2-sω
2/2+I = 0, which has the two solutions k
α
(ω,I)=[1+α(1-4D[I-sω
2/2])1/2]/(2D), where α = ±1 defines the upper and lower branches respectively.
Two cases must be distinguished:
D > 0 and
D < 0. The former represents the most interesting situation for this paper, due to the unconventional instability content of
Eq. (1). When
D < 0 the modulational instabilities of
Eq. (1) show no significant difference with the case of the non-paraxial NLSE [
12
P. Chamorro-Posada, G. S. McDonald, and G. H. C. New, “Non-Paraxial Solitons,” J. Mod. Opt.
45, 1111 (1998). [CrossRef]
,
13
J. M. Christian, G. S. McDonald, P. Chamorro-Posada, J. M. Christian, G. S. McDonald, and P. Chamorro-Posada, “Bistable Helmholtz Solitons in Cubic-Quintic Materials,” Phys. Rev. A
76, 033833 (2007). [CrossRef]
]. Therefore here we only treat the case
D > 0. A physical interpretation of these two cases in terms of the sign of effective excitonic mass is considered in Appendix A.
Let us first examine the linear excitation case (
I = 0). In anomalous dispersion conditions, i.e. for
s = +1, the linear wave number
k
L
≡
k(
ω,
I = 0) is always real and positive. Dispersion is made of two branches,
k
+ and
k
-, separated by a
k-bandgap, the extension of which is Δ
k = 1/
D, see
Fig. 1(a). For
s = -1, i.e. in normal dispersion conditions, we have an
inverted ω-bandgap (IBG), i.e. a range of frequencies in which
k
α
are purely real, see
Fig. 1(b). The extension of this gap is Δ
ω = [2/
D]
1/2. Outside this bandgap, i.e. for |
ω| > [2
D]
-1/2, wave numbers become purely imaginary, and stationary states are not possible, due to strong field absorption or gain (depending on whether Im{
k} > 0 or Im{
k}<0 respectively). Therefore, in contrast to conventional bandgaps, the IBG is the only region which allows indefinite propagation of linear plane waves.
Fig. 1. Dispersion k
± of the two branches of the carrier wave in the linear excitation case (I = 0) for (a) anomalous dispersion regime (s = +1) and (b) normal dispersion regime (s = -1). D = 1 in both cases. Solid (dashed) lines indicate the real (imaginary) parts of k
±. Blue (red) lines refer to k
+ (k
-). Gray shaded regions in (a) and (b) indicate the k-bandgap and the ω-inverted bandgap, respectively.
Let us now examine the dispersion relations in the nonlinear regime (
I ≠ 0). For
s = +1, dispersion does not change in an essential way for
I <
I
t
≡1/(4
D). At this threshold intensity, the two branches touch each other. For
I >
I
t
we have the situation shown in
Fig. 2(a). A true nonlinearly-induced bandgap in
ω appears in the dispersion (with extension Δ
ω = [2/
D-8I]
1/2), where eigenwaves possessing the full nonlinear wave number cannot propagate without being either suppressed or amplified. For
s = -1, dispersion also has a qualitative change for
I >
I
t
, see
Fig. 2(b), where the closure of the IBG is shown, and field absorption (or gain) is dominant at all frequencies.
It is thus possible to go from a dispersion which exhibits a
k-bandgap [
Fig. 1(a)] to a dispersion which shows an
ω-bandgap [
Fig. 2(a)], or from a dispersion which exhibits an inverted
ω-bandgap [
Fig. 1(b)] to a dispersion showing total absorption/gain [
Fig. 2(b)], by just increasing the incident intensity of the carrier plane wave. This effect is made possible by the simultaneous presence of both spatial and temporal dispersion terms in
Eq. (1).
3. Parametric instabilities
Analysis of parametric wave generation in the system described by
Eq. (1) follows the wellknown Benjamin-Feir scheme [
2
T. B. Benjamin, K. Hasselmann, T. B. Benjamin, and J. E. Feir, “The disintegration of wave trains on deep water. Part 1. Theory,” J. Fluid Mech.
27, 417 (1967). [CrossRef]
]. The stationary states found above are disturbed by a small perturbation field
f (
z,
t), and the eigenvalues of the operator describing the small fluctuations are monitored in the frequency domain. Thus, in
Eq. (1) we make the substitution
ψ = (
ψ
0+
ε
f)
e
ik
(
ω,
I)
z-
iωt
, where
ε is a small parameter. The leading term in the perturbation theory leads to the following equations:
where with the bar we always indicate complex conjugation. The plus (minus) sign of
α in
Eq. (2) refers to the perturbation of the steady state in the upper (lower) dispersion branch. Defining the new vector
v≡[
f,
f̄], we can write
Eq. (2) in the short form
where Â≡diag(a,ā), a≡[1-4D(I-sω
2/2)]1/2 is the complex discriminant, σ̂
z
≡diag(1,-1) is the third Pauli matrix (which commutes with Â), and matrix operator M。 is given by
In
Eq. (4), the dispersion operator induced by the perturbation is defined by
D̂(
i∂
t
)≡-
isω∂
t
+(1/2)
s∂
2
t
. We are now interested in the spatial evolution of the small perturbational plane waves (propagating on the carrier plane wave with frequency
ω and intensity
I) governed by
Eq. (3). In order to extract the stability content of
Eq. (3), we express
v as
, where
w
AS
and
w
S
are respectively the anti-Stokes (high frequency) and the Stokes (low frequency) modes. By directly substituting into
Eq. (3) we immediately obtain the secular equation for the complex wave number detuning of the perturbational wave
κ, as a function of its (real) frequency detuning
δ:
In the following we shall always assume that we work in regions of
ω that make
a real and positive, so that
 =
a1̂, where 1̂ is the identity matrix. For given
D,
s and
ω, this gives a constraint on the maximum intensity,
I < 1/(4
D)+s
ω
2/2. This constraint is due to the fact that we are only interested in carrier waves which allow stationary states, and therefore cannot have any gain or loss in the dispersion. This condition will be clarified in section 4 when solitary wave solutions will be discussed.
Equation (5) is a set of two quartic equations in
κ (one equation for the upper branch,
α = +1 and one for the lower branch,
α = -1). However,
Eq. (5) is unchanged under the discrete symmetry
α→-
α,
a→
ā, δ→-
δ, and due to the fact that instability spectra are always symmetric with respect to the
δ = 0 axis, it is sufficient to consider the case
α = +1 without any loss of generality. The imaginary parts of the four solutions
κ
j
(
δ) (
j = {1,…,4}) of
Eq. (5) determine the linear stability of the system. More specifically, if Im{
κ
j
} is negative, then frequency
ω +
δ will be unstable.
Fig. 2. Dispersion k
± of the two branches of the carrier wave in the nonlinear excitation case (I = 0.3 > I
t
) for (a) anomalous dispersion regime (s = +1) and (b) normal dispersion regime (s = -1). D = 1 in both cases. Solid (dashed) lines indicate the real (imaginary) parts of k
±. Blue (red) lines refer to k
+ (k-).
Let us now review all possible unstable scenarios shown by the solutions of
Eq. (5) for
D > 0. First of all, let us explore the anomalous dispersion case,
s = +1. In the following, we always pose for simplicity
ω = 0. As we have seen, for light intensities
I <
I
t
, the two dispersion branches of the carrier wave are always real, see
Fig. 3(a,c). Absolute values of the imaginary parts of all unstable roots (i.e. the ones for which Im {
κ
j
(
δ)}<0) as functions of
δ are shown in
Fig. 3(b,d), for two representative values of the spatial dispersion parameters
D. For
D <
D
t
≡1/(6
I), one recovers the usual two lobes of modulational instability which is known for optical fibers [
1
G. P. Agrawal, Nonlinear Fiber Optics , (Academic Press, San Diego, 2001).
], see
Fig. 3(c). For increasing values of
D, such that
D >
D
t
, novel effects due to spatial dispersion term come into play, and the zero-detuning part of the roots becomes unstable, and the instability plot resembles the shape shown in
Fig. 3(d).
Fig. 3. Parametric instabilities for anomalous dispersion regime, s = +1. (a,c) Nonlinear dispersion branches of the carrier wave, for D = 0.7<D
t
< D
c
[(a)], and for D = 0.85 > D
t
< D
c
[(c)]. Solid (dashed) lines are real (imaginary) parts of k
±. Blue (red) lines refer to k
+ (k-). (b,d) Nonlinear instability gain Im{κ
j
(δ)}<0 for the parameters described in (a) and (c) respectively. The carrier wave intensity is I = 0.22 in both cases. The values of D
t
and D
c
are 0.7576 and 1.1364 respectively.
The unstable behavior in normal dispersion conditions (
s = -1) is somewhat more interesting. As we have seen before, the full nonlinear carrier wave shows an IBG of extension Δ
ω = [2/
D-8
I]
1/2, centered around
ω = 0, see
Fig. 4(a,b). In the limit
D→0, this extension goes to infinity, together with the critical intensity
I
t
. For a fixed value of
I, and for
D∈[0,
D
t
], there are two unstable frequency regions near the boundaries of the IBG window. Outside this window, the imaginary parts of the carrier wave number is complex, and as a consequence there cannot be any steady state, due to photon absorption/gain. This is shown in
Fig. 4(a,b). The two unstable regions tend to merge at
D =
D
t
. For
D∈[
D
t
,
D
c
], where
D
c
≡1/(4
I) >
D
t
, a qualitative change occurs, and we have a situation similar to that shown in
Fig. 4(d), where it is interesting to observe a multi-valued function for the gain plot. In the range of frequency detunings where the instability gain is multi-valued, the most negative imaginary part will dominate the exponential growth of perturbations.
For D > D
c
the IBG window closes completely, and steady states over which unstable noise can grow cannot exist.
Fig. 4. Parametric instabilities for normal dispersion regime,
s = -1. (a,c) Nonlinear dispersion branches of the carrier wave, for
D = 0.7<
D
t
<
D
c
[(a)], and for
D = 0.85 >
D
t
<
D
c
[(c)]. Solid (dashed) lines are real (imaginary) parts of
k
±. Blue (red) lines refer to
k
+ (
k-). (b,d) Nonlinear instability gain Im{
κ
j
(
δ)κ<0 for the parameters described in (a) and (c) respectively. The carrier wave intensity is
I = 0.22 in both cases.
D
t
and
D
c
are the same as in
Fig. 3.
4. Bright and dark solitons
The linear stability analysis performed in the previous section shows a non-trivial unstable content for plane waves travelling according to
Eq. (1). It is well-known that parametric instabilities are connected to the existence of solitary waves in nonlinear systems. In this section we clarify this connection for
Eq. (1), and we give explicit analytical solutions for bright solitons, which vanish at the space-time boundaries. Of course an essential requirement for bright solitary wave solutions of
Eq. (1) to exist, is that they must be localized in both space (
z) and time (
t), because Kerr nonlinearity must balance both spatial and temporal dispersion terms in the equation.
In order to find bright soliton solutions for
Eq. (1), we have used the following ansatz:
where
A
0 is the soliton amplitude,
ν is the soliton velocity and
z
0 is the soliton width. After plugging expression (6) into
Eq. (1), and imposing it to vanish for any value of
z and
t, one obtains
Again,
α = ±1 in
Eqs. (7)–
(9) refers to upper and lower dispersion branches respectively. It is possible to see that solution (7)–(9) represents a family of 2-parameter solitons, parameterized by (
ν,
ω). Moreover, coefficients
A
0 and
z
0 in
Eqs. (7)–
(9) must be real, and this imposes restrictions on the values that
v ad
ω can assume. In addition, it is possible to prove with simple algebra that the following extra condition must hold:
Fig. 5. (a) Region of existence on the parameter space (
ν,
ω) for bright solitons (
s = +1,
D = 1), calculated by imposing reality conditions on
Eqs. (7)–
(9), plus the extra condition
Eq. (10). (b) Regions of existence on the parameter space (
ν,
ω) for dark (black + gray) solitons (
s = -1,
D = 1), calculated by imposing reality conditions on
Eqs. (12–
14), plus the extra condition
sνω = ±
α|
ν||
ω|, valid for the black and the gray solitons respectively.
Constraints given above are quite restrictive. Importantly, it can be demonstrated with a bit of elementary algebra that, analogously to the case of nonlinear fiber optics [
1
G. P. Agrawal, Nonlinear Fiber Optics , (Academic Press, San Diego, 2001).
,
9
Y. S. Kivshar and G. P. Agrawal, Optical Solitons , (Academic Press, San Diego, 2003).
], bright solitons can never exist for
s = -1, irrespective of the values of the parameters (
ν,
ω), in any branch of the dispersion relation. The regions of existence for the family of bright solitons found above is displayed in
Fig. 5 for the upper branch of the dispersion
k
+. This region of the parameter space is, however, profoundly different from the analogous region found for fiber solitons [
1
G. P. Agrawal, Nonlinear Fiber Optics , (Academic Press, San Diego, 2001).
,
9
Y. S. Kivshar and G. P. Agrawal, Optical Solitons , (Academic Press, San Diego, 2003).
]. Black regions in
Fig. 5 are those regions in the parameter space where reality of (7)–(9) and constraint (10) are simultaneously satisfied. For solitons living on the lower branch
k-, the region of existence (not shown here) is obtained by taking the mirror image of
Fig. 5 with respect of the
ν = 0 axis. A noteworthy feature is the onset of a velocity gap, which forbids the existence of bright solitons inside the velocity gap -[2
D]
1/2<
ν < [2
D]
1/2. This gap is filled with dark solitons, as we shall see shortly. Solitons possessing a smaller velocity
v (but such that |
v| > [2
D]
1/2), exist for a wide range of frequencies (positive or negative, depending on the sign of
v and on the specific dispersion branch
k
α
), while large velocity solutions can only exist near a narrow region around
ω = 0.
We have also found dark soliton solutions for
Eq. (1). We chose the following ansatz:
where
θ is an angular parameter which regulates the depth of the dark soliton. For
θ = 0 one has a so-called ‘black soliton’, which has a vanishing minimum of the amplitude, to be distinguished from the ‘gray soliton’, which has
θ ≠ 0. This distinction is made not only for classification purposes, but because in general black and gray solitons may have different dynamical and stability properties [
9
Y. S. Kivshar and G. P. Agrawal, Optical Solitons , (Academic Press, San Diego, 2003).
]. We have found by direct substitution of
Eq. (11) into
Eq. (1) that parameters of solution
Eq. (11) are given by the following expressions:
Note that, given
ν and
ω, the value of
θ is not arbitrary, and is determined by
Eq. (15). Also, it is clear from
Eq. (15) that the black soliton solutions (
θ = 0) and gray soliton solutions (
θ ≠ 0) have different existence regions in the parameter space (
ν,
ω). In fact the black soliton exists under the condition
sνω =
α|
ν||
ω|, while the gray soliton exists if
sνω = -
α|
ν||
ω|, where
α = ±1 correspond to the upper and lower branch of the dispersion, respectively. This can be also seen in
Fig. 5(b), where the existence region for the ‘gray’ soliton (in the upper dispersion branch) is shown by the light gray area, while the existence region for the ‘black’ soliton is shown with a black area. Note that the above dark solitons fill the velocity gap which appeared in the parameter space of the bright solitons family, see
Fig. 5(a).
A. Derivation of Eq. (1) for an off-resonance excitonic transition
In this Appendix we demonstrate that the dynamics of intense light propagation near an excitonic resonance in a 1D semiconductor waveguide in presence of Kerr nonlinearity, is regulated by
Eq. (1). The starting point of our derivation is the dielectric function of an isotropic nonmagnetic material under excitation of monochromatic light of frequency ω, which, in the proximity of an excitonic resonance reads:
In the above expression,
ε
b
is the background dielectric constant,
ω̃
0≡[
ω
2
0-
γ
2]
1/2, with
ω
0 and
γ the exciton resonance and damping, respectively, Δ≡
Ω
2
c
/(2
ω̃
0), Ω
c being a frequency proportional to the oscillator strength for the coherent exciton-photon interaction. The parameter Γ≡
h̄/(2
M*
x
), where
M*
x
is the effective exciton mass, takes into account the excitonic spatial dispersion through the wavevector dependence of the excitonic transition,
ω =
ω
0+Γ|
k|
2. Note that in expression (16) only the resonant part of the dielectric function is taken into account, while the non-resonant part is generally totally negligible [
15
H. Haug and S. W. Koch, Quantum Theory of the Optical and Electronic Properties of Semiconductors , (World Scientific, Singapore, 2004).
]. Solving Maxwell equations with the constitutive relation
D=
ε
k,
ω
E, with
ε
k,
ω given by
Eq. (16), yields the following equation for transverse modes:
which represents the dispersion of exciton-polariton waves, i.e., the true mixed modes propagating in the medium, resulting from the interaction between excitons and the retarded electromagnetic field. As first pointed out by Pekar [
16
S. I. Pekar, “Supplementary light waves in crystals and exciton absorption,” Sov. Phys. Uspekhi
5, 515 (1962). [CrossRef]
], the inclusion of the spatial dispersion in the excitonic dielectric function
Eq. (16), leads in general to a twofold solution of
Eq. (17) at a given
ω: an upper transverse polariton branch
ω
UT=
ω
UT(|
k|) and a lower transverse one
ω
LT=
ω
LT(|
k|). Polaritons in bulk semiconductors (as well as in confined structures) are well-established concepts, and the existence of the above branches has been demonstrated for various crystals in the proximity of the Wannier exciton resonances (see e.g. Ref. [
17
V. M. Agranovich and V. L. Ginzburg, Crystal Optics with Spatial Dispersion, and Excitons , (Springer, Berlin, 1984).
]). In semiconductor crystals, a Wannier exciton has typically a positive effective mass (
M*
x
> 0), and the
ω
LT(|
k|) branch is thus characterized by a positive group velocity
ν
gr≡
∂ω
LT/
∂|
k| > 0. In some special cases, however, the lower branch
ω
LT(|
k|) is characterized by a negative group velocity,
ν
gr<0, or, in other words, a negative effective mass,
M*
x
<0.
Let us now assume that we pump an intense laser beam into a semiconductor waveguide which provides an effective transverse confinement, therefore making the problem effectively one-dimensional. k = [0,0,k] is the wave vector for the effective 1D problem, the only non-vanishing component of which points along the longitudinal propagation direction Z. Then the linear constitutive relation between the Fourier components of polarization P
k
,
ω
and electric field E
k
,
ω
can be written as
where
χ
X
,
W
are the exciton and waveguide contributions to the susceptibility, respectively. Note that
χ
W
does not depend on
k, due to the absence of spatial non-locality in conventional semiconductor single-mode waveguides [
9
Y. S. Kivshar and G. P. Agrawal, Optical Solitons , (Academic Press, San Diego, 2003).
].
The linear propagation of Fourier components of the electric field (due to the exciton part of the susceptibility only) is regulated by the following equation:
where
E
(±)
k
,
ω
are respectively the forward and backward components of the electric field (which is supposed to be linearly polarized along one of the waveguide’s principal axis), and
c is the speed of light in vacuum. Taking equations only for the forward component corresponds to performing the well-known slowly-varying envelope approximation (SVEA), which is one of the cornerstones of nonlinear optics [
1
G. P. Agrawal, Nonlinear Fiber Optics , (Academic Press, San Diego, 2001).
].
Substitution of expression (16) into
Eq. (19) yields the following equation for the
forward propagation:
where
δω≡
ω-
ω̃
0. We now assume that the central frequency of the input laser is tuned slightly off-resonance, but such that |
δω|≫|-Γ
k
2+
i
γ| and |δ
ω|≫Δ, i.e. on the ‘tail’ of lower exciton-polariton branch. In this physical situation we do not need to take into account the coupled nature of the electric and polarization fields, and the knowledge of the dynamics of the electric field is sufficient to accurately model light propagation. Also, the oscillator damping constant will be negligible in this regime, and in the following we set it equal to zero for simplicity. In this limit
Eq. (20) can be approximated by
By taking the inverse Fourier transform in
k and
ω of
Eq. (21), and including the waveguide dispersion contribution as given by
Eq. (18), we obtain
where Z is the longitudinal propagation direction, T is time, A≡A(Z,T) is the slow-varying envelope of the linearly polarized electric field, ω
0 is the central pulse frequency, and γNL is the Kerr nonlinear coefficient exhibited by the semiconductor waveguide. The operator
models the waveguide temporal dispersion, and contains terms such as the group velocity [proportional to the
β
1 coefficient in
Eq. (23)] and the group velocity dispersion (GVD, proportional to the
β
2 coefficient), plus higher-order dispersion terms (
β
j
≥3) if necessary [
18
S. V. Chernikov and P. V. Mamyshev, “Femtosecond soliton propagation in fibers with slowly decreasing dispersion,” J. Opt. Soc. Am. B
8, 1633 (1991). [CrossRef]
]. For our purposes, it is sufficient to make the expansion
D̂
W
(
i∂
T
)≃
iβ
1
∂
T
-(1/2)
β
2
∂
2
T
, where
β
1,
2 are the inverse of the group velocity and the GVD parameters respectively.
β
1 can always be eliminated by a proper choice of moving coordinates without loss of generality (mixed derivatives that can arise are typically very small due to the fact that the exciton transverse dispersion is much smaller than the waveguide dispersion), and we therefore consider only
β
2 in the calculation.
Eq. (22) is the conventional generalized nonlinear Schroedinger equation (GNLSE, [
1
G. P. Agrawal, Nonlinear Fiber Optics , (Academic Press, San Diego, 2001).
,
9
Y. S. Kivshar and G. P. Agrawal, Optical Solitons , (Academic Press, San Diego, 2003).
,
18
S. V. Chernikov and P. V. Mamyshev, “Femtosecond soliton propagation in fibers with slowly decreasing dispersion,” J. Opt. Soc. Am. B
8, 1633 (1991). [CrossRef]
]), but with the additional contribution of the longitudinal spatial dispersion of the exciton, proportional to the second order derivative in space.
By introducing dimensionless variables
z≡
Z/
z
0,
t≡
T/
t
0, and
ψ≡
A/
A
0, by imposing the scalings
z
0=
t
2
0/|
β
2|,
A
0=[
γ
NL
z
0]
-1/2,
t
0 being the input pulse duration, and defining
and
, we finally obtain our master
Eq. (1).
Two physically different situations can be realized for
Eq. (1). The first regime occurs for
D < 0, i.e. when the effective exciton mass
M*
x
is positive. We name this case
normal spatial dispersion regime. This is a standard case which can be realized, for instance, for the very common semiconductor materials exhibiting conventional bulk Wannier exciton resonances. This case does not present particular interest for us, due to its equivalence to the non-paraxial NLSE case which can be found in the literature [
12
P. Chamorro-Posada, G. S. McDonald, and G. H. C. New, “Non-Paraxial Solitons,” J. Mod. Opt.
45, 1111 (1998). [CrossRef]
,
13
J. M. Christian, G. S. McDonald, P. Chamorro-Posada, J. M. Christian, G. S. McDonald, and P. Chamorro-Posada, “Bistable Helmholtz Solitons in Cubic-Quintic Materials,” Phys. Rev. A
76, 033833 (2007). [CrossRef]
], and it will be not analyzed in the present paper. The second regime, when
D > 0 (which we name
anomalous spatial dispersion regime), is characterized by a negative effective mass
M*
x
< 0 and is the most interesting one for this paper, due to the new and unconventional unstable content which can result from
Eq. (1). This case is physically realizable, for instance, in some organic molecular crystals [
19
A. S. Davydov, Theory of Molecular Excitons , (Plenum Press, New York-London, 1971).
,
20
V. M. Agranovich, Y. R. Shen, R. H. Baughman, and A. A. Zakhidov, “Optical bulk and surface waves with negative refraction,” J. Lumin.
110, 167 (2004). [CrossRef]
,
21
V. M. Agranovich and Y. N. Gartstein, “Spatial dispersion and negative refraction of light,” Phys. Usp.
49, 1029 (2006). [CrossRef]
]. However, here we consider as our representative example a structure made of a strained ZnCdSe/ZnSe superlattice, investigated in a series of recent works, where the authors clearly show that the experimental data (taken at room temperature) are fully consistent with a negative exciton-polariton dispersion, induced by high compressive strain [
22
A. Je. Semjonow, U. W. Pohl, A. Je. Semjonow, U. W. Pohl, and R. Engelhardt, “Negative exciton mass in ZnCdSe/SnSSe superlattices observed by excitonic polariton Raman scattering,” J. Phys. Condens. Matter
11, 1735 (1999). [CrossRef]
]. The excitonic feature for the type-II superlattice of [
22
A. Je. Semjonow, U. W. Pohl, A. Je. Semjonow, U. W. Pohl, and R. Engelhardt, “Negative exciton mass in ZnCdSe/SnSSe superlattices observed by excitonic polariton Raman scattering,” J. Phys. Condens. Matter
11, 1735 (1999). [CrossRef]
] is located in the visible range with a central value of
ω
0≃2.51 eV, i.e.
λ
0≃493 nm. The background dielectric constant for this material is
ε
b≃8.66, the coupling frequency is Ω
c
≃85.4 meV, and the exciton mass assumes the value
M*
x
= -0.65
m
0, where
m
0 is the free electron mass. With these values of parameters, one can calculate Γ≃-8.9 × 10
-5 m
2/sec, and Δ≃1.5 meV. By choosing
t
0≃200 fsec, and |
β
2| = 0.1 psec
2/m, and pumping the input pulse with a detuning equal to
δω≃
ω0/250, we obtain a value of the dimensionless spatial dispersion parameter equal to
D≃1.18 × 10
-2, which is relatively large and should be observable by means of the nonlinear propagation regimes described in the previous sections. The Kerr nonlinear coefficient for ZnSe can be estimated to assume the large value
γNL~0.8 m
-1W
-1 for a pulse transverse effective area of about
A
eff
~100
µm
2 [
23
A. A. Said, M. Sheik-Bahae, D. J. Hagan, T. H. Wei, J. Wang, J. Young, and E. W. Van Stryland, “Determination of bound-electronic and free-carrier nonlinearities in ZnSe, GaAs, CdTe and ZnTe,” J. Opt. Soc. Am. B
9, 405 (1992). [CrossRef]
].
One final thing to be noticed (which we have already mentioned) is the formal similarity of
Eq. (1) with the governing equation of non-paraxial solitons intensively investigated by Chamorro-Posada, see
Eq. (3) of [
12
P. Chamorro-Posada, G. S. McDonald, and G. H. C. New, “Non-Paraxial Solitons,” J. Mod. Opt.
45, 1111 (1998). [CrossRef]
]. However, one must notice the profound physical difference between the non-paraxial term, proportional to the second order derivative in the longitudinal coordinate, and our analogous term in
Eq. (1). In fact, the non-paraxial spatial dispersion coefficient is given in our notation by
D
nonpar
≡-(
λ
0/
w
0)
2/(4
π
2
ε
b
), where
w
0 is the transverse pulse width (
A
eff
≡
w
2
0).
D
nonpar
is always negative and becomes appreciable only for
w
0→
λ
0. This is in striking contrast with the physical origin of coefficient
D in
Eq. (1), attributed to exciton-polaritons, which can have both signs and is independent from the transverse width of the pulse. This means that the non-paraxial contribution is negligible when compared with the exciton-polariton contribution (i.e. |
D
nonpar
|≪|D|) whenever
w
2
0≫
λ
2
0/[4
π
2
εb
D]
1/2, a condition that in the example described above translates to the relation
w
2
0≫0.06
µm
2, which is easily satisfied for our choice of
A
eff
. It is interesting to note that non-paraxiality can be provided also by other physical situations, other than narrow beams, like for instance oblique propagation of the beam [
13
J. M. Christian, G. S. McDonald, P. Chamorro-Posada, J. M. Christian, G. S. McDonald, and P. Chamorro-Posada, “Bistable Helmholtz Solitons in Cubic-Quintic Materials,” Phys. Rev. A
76, 033833 (2007). [CrossRef]
]. However, in the present work these interpretations do not apply, due to the fact that our model has one spatial coordinate only. Moreover, solutions
Eqs. (6) and
(11) are qualitatively different from the solutions found in [
12
P. Chamorro-Posada, G. S. McDonald, and G. H. C. New, “Non-Paraxial Solitons,” J. Mod. Opt.
45, 1111 (1998). [CrossRef]
,
13
J. M. Christian, G. S. McDonald, P. Chamorro-Posada, J. M. Christian, G. S. McDonald, and P. Chamorro-Posada, “Bistable Helmholtz Solitons in Cubic-Quintic Materials,” Phys. Rev. A
76, 033833 (2007). [CrossRef]
]. In these latter works, for
D
nonpar
→0, solutions collapse to the expressions for the NLSE soliton, while this is not the case for
D→0 in
Eqs. (6) and
(11).
Our model is also suitable of further extensions, such as the inclusion of dispersive terms of order higher than the second in
Eq. (23) (that become important for short pulse durations), and of more realistic nonlinear terms, such as Raman and saturable nonlinearities.