Dynamics and stability of a new class of periodic solutions of the optical parametric oscillator

Optics Communications 240 (2004) 423–436 www.elsevier.com/locate/optcom

Dynamics and stability of a new class of periodic solutions of the optical parametric oscillator Sarah E. Hewitt *, Karen Intrachat, J. Nathan Kutz Department of Applied Mathematics, University of Washington, Seattle, WA 98195-2420, USA Received 12 April 2004; received in revised form 18 June 2004; accepted 21 June 2004

Abstract The stability and dynamics of a new class of periodic solutions is investigated when a degenerate optical parametric oscillator system is forced by an external pumping field with a periodic spatial profile modeled by Jacobi elliptic functions. Both sinusoidal behavior as well as localized hyperbolic (front and pulse) behavior can be considered in this model. The stability and bifurcation behaviors of these transverse electromagnetic structures are studied numerically. The periodic solutions are shown to be stabilized by the nonlinear parametric interaction between the pump and signal fields interacting with the cavity diffraction, attenuation, and periodic external pumping. Specifically, sinusoidal solutions result in robust and stable configurations while well-separated and more localized field structures often undergo bifurcation to new steady-state solutions having the same period as the external forcing. Extensive numerical simulations and studies of the solutions are provided. Ó 2004 Elsevier B.V. All rights reserved. PACS: 42.65.Y; 05.45; 07.05.T Keywords: Parametric ocsillators; Nonlinear dynamical systems; Modeling; Jacobi elliptic functions

1. Introduction and overview The demonstration of optical frequency conversion marked the early studies of nonlinear optics * Corresponding author. Tel.: +1-206-685-9395; fax: +1-206685-1440. E-mail address: [email protected] (S.E. Hewitt).

[1]. In the four decades since, major advances and breakthroughs have been made in a wide range of applications where a quadratic or cubic nonlinearity plays a critical role. Recently, the phenomena of frequency conversion by means of optical parametric oscillators (OPOs) has come to be regarded as a promising source of broadly tunable coherent radiation [2]. Consequently theoretical, computational, and experimental studies

0030-4018/$ - see front matter Ó 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.optcom.2004.06.043

424

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

have attempted to develop an understanding of the nontrivial spatial structures arising in the transverse electric field due to the parametric interaction of the pump and signal fields in the presence of diffraction, external pumping and attenuation. In this paper, the stability of periodic transverse structures are considered when a degenerate OPO system is forced by an external pumping field with a periodic spatial profile. The periodic spatial forcing is modeled by Jacobi elliptic functions [3,4] which can model sinusoidal behavior as well as hyperbolic (fronts and pulses) behavior. A new class of exact periodic solutions is found using this forcing and their stability and bifurcation structures are studied analytically and numerically. Few analytic solutions to the governing OPO equations are known. Aside from the well-known plane wave solutions [5], asymptotic reductions of the OPO equations are required to yield simplified models and their analytic results [6]. However, computational studies show that the parametric interaction in the OPO, when coupled with cavity diffraction, can lead to extended periodic patterns which are a result of a dissipative modulational instability [7,8]. It has also been predicted that diffractive parametric mixing in passive, dissipative cavities support localized (nonperiodic) fields which are often referred to as cavity solitons [9]. These cavity solitons can be bright spots or one-dimensional dark or ring waves and are often stable even under perturbation [10–17]. In some cases, the topological solitons can be proven to be stable analytically [18]. This motivates the consideration of a spatially periodic external pumping field which, remarkably, gives rise to three families of exact solutions modeled by Jacobi elliptic functions [3,4]. These periodic solutions are valid well beyond the onset of instability, i.e. in the strongly nonlinear regime. Further, linear stability can be applied to these periodic structures to understand the onset of secondary instabilities from these periodic patterns. The solutions allow for the consideration of sinusoidal patterns as well as wellspaced localized transverse structures (e.g. cavity solitons). The onset of periodic patterns in OPO systems are a well-known phenomena which have been widely observed in simulation [19,20]. The under-

lying mechanism giving rise to such spatial structures is an unstable growth mode (eigenvalue) of a steady-state solution which has a positive real part and an oscillatory (imaginary) component. These and other field structures, or patterns, arise in one-, two- and three-dimensional geometries and are part of a broader class of spatially extended nonlinear systems which are driven beyond equilibrium [21]. In the study of such systems, the stability of nontrivial solutions and their bifurcations are central issues [22]. Linear stability analysis is the typical analytic method used to characterize the onset of growth modes and their resulting patterns. This analysis is capable of capturing the dominant spatial frequencies and patterns observed near the onset of instability. Strictly speaking, this analysis applies only near the bifurcation point and typically only in the weakly-nonlinear regime, i.e. when the generated patternÕs amplitude fluctuations are small. In this work, we construct exact periodic solutions to the governing OPO equations which are valid well beyond the onset of instability, thus greatly extending the range of validity of the periodic solutions. Further, the Jacobi elliptic solutions allow for the consideration of the onset and growth of localized structures and secondary instabilities. This paper aims to extend the consideration of periodic patterns by providing a new set of periodic solutions which can describe very general sinusoidal and localized behavior. These solutions augment the existing literature concerning the well-known observation of periodic solutions which can only be treated analytically near the onset of instability. Here a full treatment can be given well beyond the onset of instability. Because exact solutions are known, secondary instabilities from periodic profiles can be explicitly and analytically considered via a linear stability analysis. This shows that periodic patterns can lead to secondary bifurcations which show a variety of interesting, and hitherto unobserved, spatial periodic profiles. The paper is outlined as follows: In Section 2 the governing equations along with the form of the external pumping field in terms of Jacobi elliptic functions are presented. Three solution types are also derived in this section. A linear stability analysis of the exact solutions is performed in Sec-

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

tion 3 describing the onset of instabilityÕs dependence on system parameters. Specifically, the elliptic modulus of the solutions is found to be the key bifurcation parameter. In Section 4 extensive numerical simulations are performed which help determine the stability and bifurcation structure of the derived periodic solutions. Section 5 concludes this work with a brief summary of the problem and the results.

2. Governing equations and solutions Frequency conversion in the OPO occurs when a pump field at frequency xp is directed into a birefringent quadratic crystal. The second order polarization of the v(2) material creates a nonlinear quadratic field interaction which produces two new waves at frequencies xs (signal field) and xi (idler field) [23]. The relationship between frequencies of the signal, idler and pump waves is described by the simple relation: xp ¼ xs þ xi :

ð1Þ

The frequencies of the signal and idler waves can be tuned by adjusting the angle that the pump wave is directed into the v(2) material [23]. In this paper, the degenerate regime is considered for which the signal and idler frequencies are equivalent. A variety of physical effects influence the propagation of the signal and pump fields in an OPO crystal. Of primary interest are effects arising from the interaction of diffraction, parametric coupling between the signal, idler and pump fields, attenuation, and the external pumping (driving) of the second harmonic field. For the degenerate case where the signal and idler fields coalesce, the equations governing the leading order behavior for a continuous wave can be derived from MaxwellÕs equations via a high frequency expansion [24,25]. It should be noted that the nonlinear Schro¨dinger equation for fiber optics can be derived in a similar fashion [26]. The three key asymptotic reductions for this long wavelength expansion come from applying the slowly varying envelope approximation, the paraxial approximation and rotating wave approximation in succession. The normalized system

425

of coupled, nonlinear partial differential equations in one dimension is [5] oU i o2 U ¼ þ VU
 ð1 þ iD1 ÞU ; ot 2 ox2

ð2aÞ

oV i o2 V ¼ q 2  U 2  ða þ iD2 ÞV þ Sðx; tÞ; ot 2 ox

ð2bÞ

where U(x,t) is the envelope of the signal field at frequency x and V(x,t) is the envelope of the pump field at frequency 2x. The system is pumped (driven) externally by the field S(x,t) which will be considered to be spatially modulated in x for this work. The parameter q measures the diffraction ratio between signal and pump fields while a determines the pump-to-signal loss ratio. Detuning between the signal and pump fields is measured by the parameters D1 and D2. In general it is difficult to find solutions to the governing equations (2). However, two solution forms are known. When the OPO crystal is externally pumped by a plane wave of a prescribed amplitude, plane wave solutions for the first and second harmonics can be supported. Below the instability threshold, the solution for the signal field is the trivial field (U = 0), while above threshold, there is a constant amplitude plane wave solution [5]. These solutions do not account for the wide range of spatial patterns observed in computational simulations. The two known analytic solutions for the pump and signal fields take the general form of the external pumping S(x,t). This motivates the consideration of a spatially modulated (periodic) external pumping field. Although trigonometric functions are obvious candidates for modeling the modulation, we consider the more general periodic modulation provided by Jacobi elliptic functions [3,4] (see Appendix A for a further explanation of elliptic functions). When S(x,t) is constructed using a combination of Jacobi elliptic functions, exact periodic solutions for the signal and pump fields can be generated. Specifically, the one parameter family of solutions is given by Sðx; tÞ ¼ ð3qk 4 sn4 ðx; kÞ þ C 1 sn2 ðx; kÞ þ C 2 Þe2ixt ; ð3aÞ

426

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

8 ixt > < snðx; kÞe ; U ðx; tÞ ¼ A cnðx; kÞeixt ; > : dnðx; kÞeixt ;

ð3bÞ

V ðx; tÞ ¼ ð1  ik 2 sn2 ðx; kÞÞe2ixt ;

ð3cÞ

purely sinusoidal case, allowing for the consideration of a wider variety of configurations. Note that all solutions assume periodic boundary conditions in the transverse direction. Thus the model assumes that the crystal is wide enough to accommodate many periods of the spatial modulation. And although the external forcing S(x,t) takes a very specific form, in practice the spatial structure only needs to be close to the specified form to observe the given behaviors. Thus the solutions are valid under even moderately strong perturbation of S(x,t).

3 2.5 2 1.5 1 0.5 0

2.1. CN solution For the signal field solution U(x,t) taking the form of the cn(x;k) function, the coefficients of the pump field V(x,t), external pumping S(x,t), and rotation frequency x in Eq. (3) are given by C 1 ¼ A2  k 2 ð2D1 þ 1  2qð1 þ k 2 Þ  D2 þ iaÞ; ð4aÞ Signal Field, U Pump FIeld, V External Pump Field, S

k=0.1

−2π

−π

0

π

2π

−π

0

π

2π

−π

0

π

2π

2 1.5 k=0.9

1 0.5 0 Field Magnitude

Field Magnitude

Field Magnitude

where A is a free parameter. The coefficients C1 and C2 as well as the frequency rotation x are determined by the physical parameters of the system, the value of the elliptic modulus 0 6 k 6 1, and the type of elliptic function solution of signal field U(x,t). In the limit k ! 0, the external pumping field S(x,t) in Eq. (3a) models a purely sinusoidal external pumping. For intermediate values (e.g. k < 0.9) the external potential still closely resembles the sinusoidal behavior and can be thought of as a slight perturbation away from the sinusoidal case. However, as k ! 1, the external pumping models a series of well-separated hyperbolic secant forcings. The freedom in choosing k lets us consider more general spatial pumpings than the

−2π

4 3

k=0.999

2 1 0

−2π

Fig. 1. Magnitudes of the signal U(x,t), pump V(x,t) and external pump S(x,t) fields for the exact cn(x;k) solution to the OPO governing equations for various values of the elliptic modulus k. The amplitude of the signal field is set at A = 1 and the minimum magnitude of the external pump field S(x,t) is lowered to match the minimum of the other two fields illustrated. The solutions have been scaled to have a period length of 2p and four periods are illustrated.

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

427

C 2 ¼ A2 þ a  qk 2  ið1 þ 2D1  D2 Þ;

ð4bÞ

2.3. SN solution

x ¼ D1 þ 12:

ð4cÞ

For the signal field solution U(x,t) taking the form of the sn(x;k) function, the coefficients of the pump field V(x,t), external pumping S(x,t), and rotation frequency x in Eq. (3) are given by

These solutions are depicted in Fig. 1 for various values of the elliptic modulus k.

C 1 ¼ A2  ik 2 ða  iðð1  2qÞð1 þ k 2 Þ þ 2D1  D2 ÞÞ;

2.2. DN solution For the signal field solution U(x,t) taking the form of the dn(x;k) function, the coefficients of the pump field V(x,t), external pumping S(x,t), and rotation frequency x in Eq. (3) are given by C 1 ¼ k 2 ðA2 þ 2D1 þ k 2  2qð1 þ k 2 Þ  D2 þ iaÞ; ð5aÞ C 2 ¼ A2 þ a  qk 2  iðk 2 þ 2D1  D2 Þ;

ð5bÞ

x ¼ D1 þ 12k 2 :

ð5cÞ

C 2 ¼ a  qk 2  ið1 þ k 2 þ 2D1  D2 Þ;

ð6bÞ

x ¼ D1 þ 12ð1 þ k 2 Þ:

ð6cÞ

These solutions are depicted in Fig. 3 for various values of the elliptic modulus k.

3. Linear stability The stability of the exact solutions of Section 2 can be determined from a linear stability analysis. Signal Field, U Pump FIeld, V External Pump Field, S

1.01 k=0.1 1.0 0.99

Field Magnitude

Field Magnitude

Field Magnitude

These solutions are depicted in Fig. 2 for various values of the elliptic modulus k.

ð6aÞ

−2 π

−π

0

π

2π

−π

0

π

2π

−π

0

π

2π

2.0 1.5

k=0.9

1.0 0.5 0.0

−2 π

4.0 3.0

k=0.999

2.0 1.0 0.0

−2 π

Fig. 2. Magnitudes of the signal U(x,t), pump V(x,t) and external pump S(x,t) fields for the exact dn(x;k) solution to the OPO governing equations for various values of the elliptic modulus k. The amplitude of the signal field is set at A = 1 and the minimum magnitude of the external pump field S(x,t) is lowered to match the minimum of the other two fields illustrated. The solutions have been scaled to have a period length of 2p and two periods are illustrated.

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

Field Magnitude

Field Magnitude

Field Magnitude

428

Signal Field, U Pump Field, V External Signal Field, S

3 2

k=0.1

1 0

−2π

−π

0

π

2π

−π

0

π

2π

−π

0 X

π

2π

5 4 k=0.9

3 2 1 0

−2π

8 6 k=0.999 4 2 0

−2π

Fig. 3. Magnitudes of the signal U(x,t), pump V(x,t) and external pump S(x,t) fields for the exact sn(x;k) solution to the OPO governing equations for various values of the elliptic modulus k. The amplitude of the signal field is set at A = 1 and the minimum magnitude of the external pump field S(x,t) is lowered to match the minimum of the other two fields illustrated. The solutions have been scaled to have a period length of 2p and two periods are illustrated.

The analysis is critical in determining which solutions are potentially observable in an experimental setting. Further, the onset of secondary instabilities can be determined. Linearization occurs by assuming a solution of the form U ðx; tÞ ¼ ðU 0 ðxÞ þ ~ uðx; tÞÞeixt ;

ð7aÞ

V ðx; tÞ ¼ ðV 0 ðxÞ þ ~vðx; tÞÞe2ixt ;

ð7bÞ

W t ¼ MW ;

ð8Þ

where the matrix M is given by 0 2 L U 0 ðxÞ B L 0 0 þ B M¼B @ 2U 0 ðxÞ 0 a 0

where U0 and V0 are the transverse portions of the steady-state solutions given in Eqs. (3). Specifically, U0 = A sn(x;k), A cn(x;k) or A dn(x;k) with corresponding V0 from Eq. (3c). Here ~ u; ~v  1 are small perturbations to the exact solutions. The growth or decay of these perturbations is the central concern of the linear stability analysis. Inserting Eqs. (7) into (2) and keeping only terms of size Oð~ uÞ and Oð~vÞ results in the linearized system

2U 0 ðxÞ

L0

1 0 U 0 ðxÞ C C C L0 A a ð9Þ

and T

W ¼ ðRf~ug If~ug Rf~vg If~vgÞ :

ð10Þ

Thus the growth of W, i.e. the growth of the real and imaginary parts of the perturbation, determines the stability of a given solution. The self-adjoint linear operators L± and L0 are given by L0 ¼ 

q o2  2x þ D2 ; 2 ox2

ð11aÞ

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

L ¼

1 o2 2  k 2 ðU 0 =AÞ þ ðx  D1 Þ: 2 ox2

ð11bÞ

The linerization for the three solution types, sn(x;k), cn(x;k) and dn(x;k), are easily constructed by replacing the value of U0 in Eqs. (9) and (11b) with its appropriate solution type.The stability of the solution can then be determined by letting W ¼ V ðxÞ expðktÞ;

ð12Þ

which yields the eigenvalue problem MV ¼ kV :

ð13Þ

The onset of instability occurs when the real part of any eigenvalue is greater than zero, i.e. Rfkg P 0. Analytically, it is difficult to construct the eigenvalue spectrum of the matrix (9). A great deal of this difficulty lies in the linearized operators L± which are Hill type operators that produce bang-gap spectral structures. However, the eigenvalue structure can be computed numerically. Approximating the second derivative in the linearized operators (11) with a fourth-order central differencing scheme produces a matrix whose eigenvalues can be numerically computed using standard eigenvalue solvers. Two important observations are made concerning the eigenvalue structure. First, the operator L± is independent of D1, D2 and q. Indeed, it is only a function of the elliptic modulus k. This can be easily verified by noting that the difference (x  D1) in the definition of L± is only k dependent. Thus the stability will depend critically on the elliptic modulus due to its role in L±. Second, the operator L0 depends on q, D1, D2 and k. Thus all the stability dependence on the diffraction ratio and detuning parameters is manifested through the operator L0. To illustrate the stability properties of the various solutions and their dependence upon the system parameters, the eigenvalues of (9) are numerically computed. As a specific example, the linear stability of the cn(x;k) solution is considered as a function of the elliptic modulus k with the parameters D1 = 1, D2 = 1, a = 1, and q = 1/2. Fig. 4 depicts the location of the eigenvalues in the complex plane. Note that the bulk of the spectra lies at Rfkg ¼ 1 due to the choice of linear damping a = 1. However, there are eigenvalues which move

429

along the real axis according to the value of the elliptic modulus k. In particular, as k increases from k = 0, pairs of eigenvalues are ejected from Rfkg ¼ 1 in opposite directions. For the case of k = 0.7 (circles), the eigenvalues which are ejected remain in region Rfkg 6 0 so that no instability results. However, for the case of k = 0.9 (stars), two eigenvalues have crossed over into the right half plane where Rfkg P 0 and instability occurs. In this case, the onset of instability is expected. Moreover, the form of the instability is expected to take the shape of the eigenfunctions associated with the unstable modes. Inset in Fig. 4 is the shape of the two unstable eigenmodes for the case of k = 0.9. It will shown in the full numerical computations of Eq. (2) that follow that these eigenfunction shapes are indeed characteristic of the unstable growth. To more accurately characterize the transition between stable and unstable solutions, the eigenvalue with the largest real part can be tracked as a function of the elliptic modulus. Fig. 5 gives the real part of this eigenvalue as a function of the elliptic modulus for the cn(x;k) solution with D1 = 1, D2 = 1, a = 1, and q = 1/2. This shows that the onset of instability occurs when k  0.77. Thus bifurcation to a secondary solution will occur at this point. The onset of instability is a function not only of k, but also of the other physical system parameters. Typically, q = 1/2 from physical considerations and we can scale a = 1 without loss of generality. This then leaves three free parameters: D1, D2 and k. By keeping track of the real part of the largest eigenvalue as a function of these three parameters, a three-dimensional surface plot can be made which characterizes the separation between the regions of stability and instability. Fig. 6 illustrates the three dimensional surface separating the stable (below) and unstable (above) regions in the parameter space D1, D2 and k for the cn(x;k) solution. Thus the linear stability considerations give the critical information necessary to ascertain the practical feasibility of observing the generated solutions (3) in experiment. The stability surfaces for the sn(x;k) and dn(x;k) solutions can be determined in a similar fashion. Note that although the external forcing S(x,t) takes a very specific form, in practice the spatial structure only

430

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

Eigenfunction

|u| Spectrum 10000 k=0.7, stable k=0.9, unstable

|v|

Imaginary{λ}

5000

0

x

5000

-10000 -2.5

-2

-1.5

-1 Real{λ}

-0.5

0

0.5

Fig. 4. The numerically computed spectrum of the matrix M, given in Eq. (9), linearized about the cn solution with system parameters D1 = D2 = 1, q = 1/2, and a = 1. For elliptic modulus k = 0.7 (circles) the eigenvalues lie in the left half plane, thus ensuring stability. For k = 0.9 (stars) two eigenvalues have positive real parts which leads to instability. The inset depicts the modulus of the unstable eigenfunctions for one period length in the transverse dimension.

0.35

Eigenvalue with most positive real part

Real Part of Eigenvalue

0.3 0.25 0.2 0.15 0.1 0.05 0

0.7

0.8

0.9

1

Elliptic Modulus, k

Fig. 5. The value of the real part of the eigenvalue with largest growth rate is plotted as a function of elliptic modulus k. The eigenvalues of the matrix M are computed for the cn(x;k) solution with system parameters D1 = D2 = 1, q = 1/2, and a = 1.

Fig. 6. Stability surface which separates the region of stable behavior (below) from the unstable region (above). The minimum elliptic modulus value that induces bifurcation for the cn(x;k) solution is computed as a function of the detuning parameters D1 and D2. All cn(x;k) solutions with elliptic modulus above the surface for a detuning parameter will bifurcate. Note that a = 1 and q = 1/2.

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

needs to be close to the specified form to observe the given behaviors. Thus the solutions are valid under even moderately strong perturbation of S(x,t).

4. Computational results New solutions (3) and their linear stability properties have been studied in Sections 2 and 3. This section provides computational evidence for the onset of instability which agrees well with the linear stability analysis of the previous section. Further, the computations allow us to observe the development of secondary instabilities and the formation of new, hitherto unobserved, steady-state solutions. The stability and dynamics of the sn(x;k), cn(x;k) and dn(x;k) solutions are all considered under a wide range of parameter space. The linear stability analysis suggests this is appropriate given the dependence of stability on the physical parameters of the OPO system. To be plotted are the moduli of the signal and second harmonic fields along with the spatially modulated pumping field. In-depth study is provided for particularly interesting parameter regimes. In order to limit the size of the parameter space, the physical parameters a = 1 and q = 1/2 are fixed [5]. The periodic nature of the solutions considered here suggest the implementation of Fourier spectral methods [27]. Thus, a filtered pseudo-spectral method in space is coupled with a fourth-order Runge–Kutta algorithm in time t. This numerical procedure combines the advantages of split-step and explicit Runge–Kutta methods. All simulations were completed for large spatial domains containing at least four periods of the solution. The figures presented crop the results to two periods in order to more effectively demonstrate the dynamics. For certain parameter ranges, the solutions given by Eq. (3) bifurcate [22], i.e. the onset of instability is observed. The bifurcation point is defined numerically as the point for which a ten percent discrepancy exists in the norm measuring the difference in the absolute value be-

431

tween the magnitude of the numerical solution and its corresponding analytical solution. As noted in the last section, a key parameter responsible for inducing this bifurcation is the elliptic modulus k. This is due to the dependence of the linearized operator L± on k. Thus when the solutions are in the hyperbolic limit, thus exhibiting pulse and front-type behavior, they are more likely to be unstable and bifurcate to a new steady-state solution not described by the analytic solutions found in this work, i.e. they are above the surface shown in Fig. 6 which are unstable. Of course, the detuning parameters D1 and D2 modify this to some extent as illustrated in Fig. 6. Stable evolution is found in the sinusoidal limit for the solutions (U) given by Eqs. (3) and (4), (5), and (6), respectively along with the corresponding pump field (V) and external forcing (S) are stable. Specifically, with the physical parameters of the system D1 = D2 = 1, and elliptic modulus k = 0.5, no bifurcation occurs and the solutions are stable. These solutions resemble the small k limit of Figs. 1–3. This

Fig. 7. Evolution of the cn(x;k) signal field solution (U) given by Eqs. (3) and (4) along with the corresponding pump field (V) and external forcing (S(x,t)). The solution bifurcates at t = 80 with the original minimums becoming maximums and vice versa. Although the fields develop an oscillatory envelope, they retain the same period as the external pump field S. The physical parameters of the system are D1 = D2 = 1. The elliptic modulus is k = 0.999.

432

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

agrees with the linear stability analysis of the last section. The stability of these sinusoidal signal and pump fields holds for the elliptic modulus 0 6 k  0.79 (see Fig. 13). As the elliptic modulus is increased towards unity, the elliptic functions approach the hyperbolic limit of fronts and pulses given by hyperbolic secants and hyperbolic tangents [3,4]. In this case, a wide variety of phenomena can occur as the number of unstable modes increases. The linear stability analysis shows that in this case one or more eigenvalues cross into the right half plane and result in growth modes and the onset of instability. Figs. 7–10 show the different bifurcations which result in four different parameter regimes for the cn(x;k) solution. These parameter regimes consider different values of the cavity detuning parameters D1 and D2. All the solutions are found to

Fig. 8. Evolution of the cn(x;k) signal field solution (U) given by Eqs. (3) and (4) along with the corresponding pump field (V) and external forcing (S(x,t)). The solution bifurcates at t = 80 with the original minimums becoming maximums and vice versa as in Fig. 7. In this case, the transverse width containing a single peak now has a local maximum and minimum after bifurcation. Although the fields develop an oscillatory envelope, they retain the same period as the external pump field S. The physical parameters of the system are D1 = 1, D2 = 0. The elliptic modulus is k = 0.999.

bifurcate at t  80. Hallmark features of the new steady-state solutions are a small amplitude modulation across the spatial extent of the periodic solutions. Note that each new steady-state solution has a distinct transverse spatial structure which is a result of the interaction with the solution and the unstable growth mode eigenfunction (see, for instance, Fig. 4). Further, for the cn(x;k) solutions the maximum and minimum switch locations between the original solution and the new bifurcated steady-state. The dn(x;k) and sn(x;k) solutions also bifurcate for sufficiently high elliptic modulus k. Figs. 11 and 12 depict the bifurcation of these two solutions to their new steady-state solution. Unlike the cn(x;k) solution the maxima and minima do not exchange places. Nevertheless, the new steady-state solutions exhibit some of the same features such as those associated with the cn(x;k) solution bifurcations. The numerical simulations suggest the critical role played by the elliptic modulus on the stability and bifurcation structure of the solutions. This simply confirms the linear stability analysis for which the L± operators play a

Fig. 9. Evolution of the cn(x;k) signal field solution (U) given by Eqs. (3) and (4) along with the corresponding pump field (V) and external forcing (S(x,t)). The solution bifurcates at t = 80 with the original minimums becoming maximums and vice versa as in Figs. 7 and 8. Although the fields develop an oscillatory envelope, they retain the same period as the external pump field S. The physical parameters of the system are D1 = D2 = 0. The elliptic modulus is k = 0.999.

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

Fig. 10. Evolution of the cn(x;k) signal field solution (U) given by Eqs. (3) and (4) along with the corresponding pump field (V) and external forcing (S(x,t)). The solution bifurcates at t = 80 with the original minimums becoming maximums and vice versa as in Figs. 7–9. Although the fields develop an oscillatory envelope, they retain the same period as the external pump field S. The physical parameters of the system are D1 = 0, D2 = 1. The elliptic modulus is k = 0.999.

433

Fig. 12. Evolution of the sn(x;k) signal field solution (U) given by Eqs. (3) and (6) along with the corresponding pump field (V) and external forcing (S(x,t)). While the sn(x;k) solution does bifurcate at t = 170, it retains more of its original structure than the cn(x;k) and dn(x;k) cases shown in Figs. 7 and 11. A small amplitude oscillatory envelope develops on the pump field while the previously flat regions of the signal field become more rounded with an oscillatory envelope. The physical parameters of the system are D1 = D2 = 1. The elliptic modulus is k = 0.999.

Bifurcation Position (t)

1000 Cn Analytic Limit Dn Analytic Limit Cn Bifurcation Position Dn Bifurcation Position

800

600

400

200

Fig. 11. Evolution of the dn(x;k) signal field solution (U) given by Eqs. (3) and (5) along with the corresponding pump field (V) and external forcing (S(x,t)). The solution bifurcates at t = 90. Unlike the cn(x;k) solution, the original minima and maxima locations are preserved. The physical parameters of the system are D1 = D2 = 1. The elliptic modulus is k = 0.999.

0 .77

.8

.85 .9 Elliptic Modulus (k)

.95

1.0

Fig. 13. This figure depicts the dependence the bifurcation position has on the value of the elliptic modulus k. The data was complied using the system parameters described previously.

434

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

critical role. Sinusoidal solutions for which 0 6 k  0.79 do not undergo bifurcations (provided D1 and D2 remain small) whereas well-separated pulses and fronts, which are characterized by k P 0.9, interact in such a way to bifurcate to new and interesting steady-state solutions not described previously. The importance of the elliptic modulus on the bifurcation structure is illustrated in Fig. 13. It portrays the dependence of the bifurcation point upon the value of the elliptic modulus k. These stability results are generated from computations of the governing OPO equations. Included in the figure is the linear stability analysis prediction of when transition from stable to unstable occurs. Note the excellent agreement of the linear stability analysis with the OPO simulations. It should also be noted that the free parameter A has an effect on the bifurcation position. Specifically, the larger the value of A, and corresponding solution amplitudes, the longer the solution propagates in t before bifurcating. 5. Conclusions The stability and dynamics of a new class of periodic, nontrivial field structures are examined in a degenerate OPO system forced by an external pumping field with a periodic spatial profile. The periodic spatial forcing is modeled by Jacobi elliptic functions [3,4] which are general enough to model both sinusoidal behavior as well as hyperbolic (fronts and pulses) behavior. Exact solutions are found using this forcing and their stability and bifurcation structures are studied. Analytically, a linear stability analysis is developed which characterizes the stability of the periodic solutions and the onset of secondary instabilities as a function of the physical parameters. These stability results are verified computationally with direct numerical

integration of the governing OPO equations. Thus a new class of periodic solutions for the OPO are found which are stabilized by the parametric interaction between the pump and signal fields interacting with the cavity diffraction, attenuation, and periodic external pumping. Many of these solutions resemble well-known localized cavity solitons in the limit as the elliptic modulus k ! 1, while some are simple sinusoidal behaviors. From a physical point of view, the new solutions are generated by modulating the external pump field in the transverse direction. The modulated pump field can also be generalized to twoand three-dimensions. This work will be presented elsewhere. Along with the signal and pump detuning, the elliptic modulus appears to control the bifurcation structure of the solutions. In particular, sinusoidal solutions result in robust and stable configurations while well-separated and more localized field structures often undergo bifurcation to new steady-state solutions having the same period as the external forcing. The field structures of the bifurcated solutions are complicated and not easily represented by simple elliptic, sinusoidal or hyperbolic functions. However, the form of the new steady-state solutions appear to be a combination of the periodic solution with the eigenfunctions of the unstable eigenvalues. The analytic stability of these new solutions, which is beyond the scope of the present work, is difficult to describe due to the periodic nature of the linearized operators and their band-gap spectra. However, one can imagine the onset of new instabilities along these secondary branches of solutions. Appendix A. Jacobi elliptic functions There are three different Jacobi elliptic functions that are considered in this paper: sn(x;k),

Table 1 The variables u, k, and / are those described in Eq. (A.1) Function

Trigonometric limit k ! 0

Hyperbolic limit k ! 1

sn(u;k) = sin(/) cn(u;k) = cos(/) qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi dnðu; kÞ ¼ 1  k 2 sin2 ð/Þ

sin(u) cos(u)

tanh(u) sech(u)

1

sech(u)

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

cn(x;k), and dn(x;k). They arise from the inverse integral of the first kind [3,4] Z / dh pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ; u¼ ðA:1Þ 0 1  k 2 sin2 h where each individual function is generated from this integral (see Table 1). There are two independent variables that are associated with this integral: the spatial profile u and the elliptic modulus k. The elliptic modulus, which ranges from 0 6 k 6 1, determines the general shape of the function. When k is close to unity, the function is said to be near its hyperbolic limit, where it will exhibit front (hyperbolic tangent) or pulse (hyperbolic secant) like behavior, depending on the specific elliptic function. For all other values of the elliptic modulus the function is said to be in the trigonometric limit, where the functions exhibit behavior similar to standard sinusoidal functions or constants. The general, local behavior in these limits is described in Table 1. Elliptic functions share some convenient properties of trigonometric functions. This, along with the fact that many types of behavior may be modeled with a single function, makes working with elliptic functions appealing. Here are some helpful identities [3,4] for generating solutions to Eq. (2). sn2 ðu; kÞ þ cn2 ðu; kÞ ¼ 1;

ðA:2aÞ

k 2 sn2 ðu; kÞ þ dn2 ðu; kÞ ¼ 1;

ðA:2bÞ

snð0; kÞ ¼ 0;

cnð0; kÞ ¼ 1;

dnð0; kÞ ¼ 1; ðA:2cÞ

d2 snðu; kÞ ¼ 2k 2 sn3 ðu; kÞ  ð1 þ k 2 Þsnðu; kÞ; du2 ðA:2dÞ d2 cnðu; kÞ ¼ ð2k 2  1Þcnðu; kÞ  2k 2 cn3 ðu; kÞ; du2 ðA:2eÞ d2 dnðu; kÞ ¼ ð2  k 2 Þdnðu; kÞ  2dn3 ðu; kÞ: du2 ðA:2fÞ

435

Of greatest importance is the fact that differentiating twice reproduces the function and its cube, thus allowing it to be used as a solution for the OPO and other nonlinear Schro¨dinger type models [28] which are dominated by a cubic nonlinearity.

References [1] N. Bloemberger, Nonlinear Optics, Benjamin, New York, NY, 1965. [2] See, for instance, the special issue on Optical Parametric Oscillation and Amplification J. Opt. Soc. Am. B 10 (1993) 1655; and Optical Parametric Devices. J. Opt. Soc. Am. B 12 (1995) 2083. [3] M. Abramowitz, I. Stegun, Handbook of Mathematical Functions, Dover, New York, 1972. [4] H.T. Davis, Introduction to Nonlinear Differential and Integral Equations, Dover, New York, 1962. [5] S. Trillo, M. Haelterman, A. Sheppard, Opt. Lett. 22 (1997) 970. [6] S. Trillo, M. Haelterman, Opt. Lett. 23 (1998) 1514. [7] G.L. Oppo, M. Brambilla, L.A. Lugiato, Phys. Rev. A 49 (1994) 2028. [8] R.A. Fuerst, D.M. Baboiu, B. Lawrence, W.E. Torruellas, G.I. Stegeman, S. Trillo, S. Wabnitz, Phys. Rev. Lett. 78 (1997) 2756. [9] See, e.g. N.N. Rozanov, in: E. Wolf (Ed.), Progress in Optics, vol. 35, 1996, p. 1. [10] S. Longhi, Phys. Scripta 56 (1997) 611. [11] K. Staliunas, V.J. Sa´nchez-Morcillo, Opt. Commun. 139 (1997) 306. [12] K. Staliunas, V.J. Sa´nchez-Morcillo, Phys. Rev. A 57 (1998) 1454. [13] S. Longhi, Opt. Commun. 149 (1998) 335. [14] M. Tidli, P. Mandel, R. Lefever, Phys. Rev. Lett. 81 (1998) 979. [15] C. Etrich, U. Peschel, F. Lederer, Phys. Rev. Lett. 79 (1997) 2454. [16] M. Santagiustina, E. Hernandez-Garcia, M. San-Miguel, A.J. Scroggie, G.-L. Oppo, Phys. Rev. E 65 (2002) 036610. [17] G.-L. Oppo, A.J. Scroggie, W. Firth, Phys. Rev. E 63 (2001) 066209. [18] J.N. Kutz, T. Erneux, S. Trillo, M. Haelterman, J. Opt. Soc. Am. B 16 (1999) 1936. [19] M. Tildi, P. Mandel, M. Halterman, Phys. Rev. E 56 (1997) 6524. [20] G.J. de Valca´rcel, K. Staliunas, E. Rolda´n, V.J. Sa´nchezMorcillo, Phys. Rev. A 54 (1996) 1609. [21] M.C. Cross, P.C. Hohenberg, Rev. Mod. Phys. 65 (1993) 851. [22] P.G. Drazin, Nonlinear Systems, Cambridge, UK, 1992.

436

S.E. Hewitt et al. / Optics Communications 240 (2004) 423–436

[23] Boyd, Nonlinear Optics, Academic Press, San Diego, CA, 2003. [24] G.A. Bukauskas, V.I. Kabelka, A.S. Piskarskas, A.Yu. Stabinis, Sov. J. Quant. Elec. 4 (1974) 290. [25] R. Danelyus, G. Dikchyus, V. Kabelka, A. Piskarskas, A. Stabinis, Sov. J. Quant. Elec. 7 (1977) 1360.

[26] G.P. Agrawal, Nonlinear Fiber Optics, Academic Press, San Diego, CA, 2001. [27] L.N. Trefethen, Spectral Methods in MATLAB, SIAM, Philadelphia, PA, 2000. [28] J.C. Bronski, L.D. Carr, B. Deconinck, J.N. Kutz, Phys. Rev. Lett. 86 (2001) 1402.