Optimization-based calculation of optical nonlinear processes in a micro-resonator

We present a new method of calculating the performance of nonlinear processes in a resonator. An optimization-based approach, conceptually similar to techniques used in nonlinear circuit analysis, is formulated and used to find the wave magnitudes that satisfy all of the boundary conditions and account for nonlinear optical effects. Unlike previous solution methods, this technique is applicable to any nonlinear process (second-order, third-order, etc.) and multiple coupled resonators, maintains the phase relations between the waves, and is exact. Examples are given for second-order nonlinear processes in a one-dimensional resonator. © 2006 Optical Society of America OCIS codes: (190.4360) Nonlinear Optics, devices; (190.2620) Frequency conversion; (230.3990) Microstructure devices References and links 1. J.A. Armstrong, N. Bloembergen, J. Ducuing and P. S. Pershan, “Interactions between light waves in a nonlinear dielectric,” Phys. Rev. 127, 1918–1939 (1962). 2. A. Ashkin, G. D. Boyd, and J. M. Dziedzic, “Resonant optical second harmonic generation and mixing,” IEEE J. Quantum Electron. QE-2, 109–124 (1966). 3. V. Berger, “Second-harmonic generation in monolithic cavities,” J. Opt. Soc. Am. B 14, 1351–1360 (1997). 4. C. Simonneau, et al, “Second-harmonic generation in a doubly resonant semiconductor microcavity,” Opt. Lett. 22, 1775–1777 (1997). 5. G. Klemens, C.-H. Chen, and Y. Fainman, “Design of optimized dispersive resonant cavities for nonlinear wave mixing,” Opt. Express 13, 9388–9397 (2005). http://www.opticsinfobase.org/abstract.cfm?URI=oe-13-23-9388 6. D.S. Bethune, “Optical harmonic generation and mixing in multilayer media: analysis using optical transfer matrix techniques,” J. Opt. Soc. Am. B 6, 910–916 (1989). 7. A. Calderone and J. Vigneron, “Computation of the electromagnetic harmonics generation by stratified systems containing nonlinear layers,” Int. J. Quantum Chem., 70, 763–770 (1988). 8. M.G. Martemyanov, T.V. Dolgova and A.A. Fedyanin,“Optical third-harmonic generation in one-dimensional photonic crystals and microcavities,” J. Exp. Theor. Phys., 98, 463–477 (2003). 9. S. El-Rabaie, V. F. Fusco, and C. Stewart,“Harmonic balance evaluation of nonlinear microwave circuits–a tutorial approach,” IEEE Trans. on Education, 31, 181–192 (1988). 10. I. Shoji, T. Kondo, A. Kitamoto, M. Shirane and R. Ito, “Absolute scale of second-order nonlinear-optical coeffients,” J. Opt. Soc. Am. B 14, 2268-2294 (1997). 11. W. J. Tropf, M. E. Thomas, and T. J. Harris, ”Properties of crystals and glasses,” in Handbook of Optics: Volume II, McGraw-Hill (1995). 12. K. E. Atkinson, An Introduction to Numerical Analysis, ch. 8, John Wiley (1989). 13. C. H. Chen, K. Tetz, W. Nakagawa, and Y. Fainman, “Wide-field-of-view GaAs/AlxOy one-dimensional photonic crystal filter,” Appl. Opt. 44, 1503–1511 (2005). #74053 $15.00 USD Received 14 August 2006; revised 25 September 2006; accepted 25 September 2006 (C) 2006 OSA 16 October 2006 / Vol. 14, No. 21 / OPTICS EXPRESS 9864 14. M.M. Fejer, G.A. Magel, D.H. Jundt, and R.L. Byer, “Quasi-Phase-Matched Second Harmonic Generation : Tuning and Tolerances,” IEEE J. Quantum Electron. 28, 2631–2654 (1992).


Introduction
Resonant structures have an important role in the enhancement of nonlinear optical processes, although techniques for the analysis of such structures have not been well developed.Secondharmonic generation (SHG) in a cavity was described in 1962 [1] and measured experimentally soon after [2].For SHG, achieving resonance at the fundamental frequency would raise the input field magnitudes, and designing the cavity to resonate at the harmonic frequency would increase output power by multiple passes.The concept was later extended to include doublyresonant cavities, which achieve resonance at the fundamental and the harmonic [3,4], and triply-resonant cavities for three-wave parametric processes [5].
The analyses of these structures has followed an ad hoc approach, rather than a general technique applicable to all cavities and nonlinear processes.One approximate approach applicable to the case of SHG in a singly-resonant cavity is to introduce a multiplier based on the number of passes a beam makes in the cavity (the cavity finesse) [2].Another method is to first use a linear solution method to calculate the fundamental wave magnitude, then calculate the nonlinear polarization in the cavity from the fundamental wave.The harmonic wave magnitude is then solved using the nonlinear polarization as a distributed source [6,7].The case of third-harmonic generation can be similarly solved [8].This method has also been applied to second-order parametric processes, in which first the pump and idler wave magnitudes are found and then used to calculate the signal wave magnitude [5].These methods are inadequate for a general analysis method because they do not account for the full interference effects within the cavity.
We introduce a new and general method for the calculation of fields in a resonator, based on optimization.In this new method, values of the field magnitudes entering and exiting the cavity at steady state are assumed, and those fields are then propagated through the cavity using the transfer matrix method for the dielectric interfaces and using coupled differential equations for the nonlinear wave mixing.If the assumed values are correct, the calculations will produce results that are self-consistent, as defined below.If the assumed values are incorrect, then a residual is produced, and an optimization routine can repeat the process until the residual is brought close to zero.Note that this procedure is conceptually similar to the harmonic balance method which is used to analyze resonant nonlinear circuits and is superior to transient analysis for these cases [9].In that method, the branch currents between the linear and nonlinear circuit portions are assumed, and then the resulting voltages calculated, producing either a selfconsistent result or a residual to be minimized.This method converges much faster (by orders of magnitude) than a transient analysis for large nonlinear circuits forming resonators.Similarly, when applied to a nonlinear optical cavity, this method bypasses the slow and complicated transient analysis, and calculates the steady-state directly.
Below we demonstrate the use of this method for the case of second-order nonlinear processes in a single cavity.We assume normal incidence for the incoming beams, although the method can be extended to TE or TM incident waves at any angle of incidence.In the next section we offer a cursory review of the solution for a linear resonator using transfer matrices.We then show how the method can be modified to account for nonlinearity in the cavity, and present a solution procedure.In Section 3 we present a case of second harmonic generation in a microcavity as a numerical example, followed by concluding remarks in Section 4.
Layer 2 (Cavity) Layer 1 (Exterior) Layer 3 (Exterior) Fig. 1.Field analysis of a single cavity.The operator A represents a phase-propagation matrix in the linear case, or coupled differential equations for a nonlinear dielectric.Analysis from left-to-right or right-to-left will produce self-consistent results for a correct set of vectors.

Review of linear cavity calculation method
The transfer matrix method can be used to find the field distribution in a set of dielectric slabs with incident field of frequency ω.Besides an overall reflection and transmission coefficient, this method provides the forward and backward traveling wave in each layer.The method is diagrammed in Fig. 1 for the case of a single dielectric layer.Following convention, we represent the fields in each layer by a vector consisting of the complex magnitude of the forward-traveling wave (E + ) and that of the backward-traveling wave (E − ).Assuming the spaces outside the cavity are linear, then the incoming and outgoing waves are completely represented by vectors where E 1 is valid immediately to the left of the first interface, and E 3 applies to the immediate right of the last interface.The vector E 2 , representing the fields internal to the cavity, varies along the cavity length.
From electromagnetic boundary conditions, the wave vectors as defined above on each side of an interface are related via a matrix M, defined by with ρ and τ the boundary reflection and transmission coefficients, respectively.Propagation across the cavity is characterized by an operator, A(•), which in the linear case is a phaseshifting matrix given by where d is the cavity length and k is the propagation constant for a wave in the cavity at frequency ω (for cavity index n and vacuum speed of light c, the relation is k = ωn/c).
Using the transfer matrices to propagate through the cavity from left to right produces the characteristic equation where we define M to be the transfer matrix from inside the cavity going out, and M −1 is the transfer matrix going into the cavity.For the case of a wave entering from the left, and no and solve for r and t, the overall reflection and transmission coefficients.A valid set of field magnitudes will also be self-consistent.Propagating vector E 1 from left-to-right will produce vector E 3 , and propagating vector E 3 from right-to-left will produce vector E 1 .This suggests that an alternative method of analysis is to try different values of field magnitudes until a selfconsistent set is found.For the linear case, however, this method is relatively unwieldy since a simple matrix algebra based solution is available.

Generalization to nonlinear cavities
The transfer matrix method of the previous section can be generalized to account for nonlinear processes by the definition of the operator A(•).Unlike the linear case, this operator is not, in general, a matrix; therefore, the straightforward solution using linear algebra is not applicable.Specifically, for the case of second-harmonic generation, propagation through the cavity is described by a set of coupled differential equations [1]: where A m , n m and α m are the field magnitudes, the refractive indices and the absorption coefficients, with m = 0, 1 for pump and harmonic, respectively.The phase mismatch term, Δ k is defined as Δ k = k 1 − 2k 0 , with k 0 and k 1 the propagation constants of the fundamental and harmonic, respectively.For parametric processes involving waves of three frequencies, the operator A(•) becomes a set of three coupled differential equations [1]: with subscripts p, s, and i, denoting pump, signal and idler, respectively, and phase mismatch term Δ k = k p − k s − k i , with k p , k s and k i the propagation constants of the pump, signal and idler waves, respectively.Similar equations can be written for higher-order nonlinear processes.Note that the slowly-varying envelope equation has been used to remove second-derivative terms.For cases where this approximation does not apply, the second-derivative terms can be included in the definition of operator A(•), at the expense of computation time.
Although the calculation method has changed, the characteristic relations still apply to any valid set of fields, where E 1 and E 3 now represent vectors with more than two elements, for rightward and leftward traveling waves of all the frequencies involved.The updated analysis method is diagrammed in Fig. 1.If we define a residual, R, then the characteristic relations are satisfied when R is at its minimum value of zero.Alternative definitions of the residual, such as the root-mean-square of the difference values, are also possible.Squaring the difference terms would increase the relative weight of the large field values, such as the pump, in the residual.The matrices M and M −1 , and, under the assumption of lossless dielectrics, the operators A(•) and A −1 (•) conserve energy between their inputs and outputs, so the field vectors E 1 and E 3 will contain the same energy when the self-consistency condition is met.A term comparing the energy at the two sides of the cavity can still be included in the residual to further enforce conservation of energy.The solution process is now formulated as an optimization problem, in which the residual is evaluated for different assumed field magnitudes, and those magnitudes are then adjusted to lower the residual.The convergence criteria is application-specific, but a useful general condition is that the residual should be low enough to resolve the angle of the weakest field.If an adequate optimization method is used, the fields values will converge to the correct values.

Second-harmonic generation, singly-resonant cavity
As examples of using this technique we consider second-order nonlinear wave mixing in a microcavity of gallium arsenide (GaAs).The nonlinear coefficient of GaAs is relatively high at 100 pm/V [10], and Sellmeier coefficients characterizing the dispersion are readily available [11].Once the initial assumption of field magnitudes had been made, we used the conjugate gradient method [12] to find new field values that lowered the residual.The unknowns to be found are the reflected and transmitted waves for the pump and harmonic waves.Considering either real and imaginary or magnitude and phase, there are eight unknowns.The implementation of this technique that we used is summarized by these steps: 1.The external fields at the left, E 0 left , E 1 left , and the right, E 0 right , E 1 right , of the cavity are transformed using Equation 2 into the fields within the cavity, E 0 * left , E 1 * left , E 0 * right , E 1 * right .These vectors are made up of the forward-traveling and backward-traveling waves as in Equation 1.

The forward-traveling components of E 0 *
left and E 1 * left are propagated through the cavity using Equations 6.This yields two complex numbers, E 0 * + le f t and E 1 * + le f t , which are the forward-traveling waves at the right interface.Also, the backward-traveling components of E 0 * right and E 1 * right are propagated through the cavity, yielding E 0 * − right and E 1 * − right .
3. The residual is assembled by comparing E 0 * + le f t and E 1 * + le f t with the forward-traveling components of E 0 * right and E 1 * right , and E 0 * − right and E 1 * − right with the backward-traveling components of E 0 * left and E 1 * left .The effect of cavity resonance is demonstrated in Fig. 2, where the optimization method was used to make the calculations.The generated harmonic energy as a function of wavelength for several values of mirror reflection coefficients shows the enhancement provided by the cavity resonance.In this example the mirror reflectivity is purely real and applied only to the pump, producing a singly-resonant cavity.It is seen that the conversion efficiency is proportional to the cavity finesse, which is equivalent to the number of passes the pump beam makes in the cavity.

Second-harmonic generation, doubly-resonant cavity
A further example of varying the cavity length demonstrates the importance of using an exact solution method that preserves phase relations over an approximate method.In this case we  assume the cavity mirrors are designed so that their complex reflection coefficients cancel the phase accumulation for each pass through the cavity for both the pump and the harmonic.For a wave with propagation relation exp(− jkd), we design the mirror reflection coefficient to have phase term exp(+ jkd).By resetting the phases of the waves at each mirror, the Fabry-Perot resonance condition can be satisfied for any length cavity.Furthermore, the nonlinear phase mismatch is compensated.The mirrors can be created, for example, with stacks of dielectric layers [5].The field magnitudes within the dielectric stack of the mirrors are assumed to be much smaller than in the cavity, allowing us to neglect the nonlinear conversion happening within the mirror layers.For a more detailed analysis, the mirror layers can be included in the modeling as coupled cavities.Any number of cavities can be included in the analysis, at the expense of solution time.This method of phase-mismatch compensation can be understood as a generalization of periodic poling (quasi-phase matching).In that method, the non-centrosymmetric lattice is reversed periodically to create a phase shift of π radians in order to reset the phase mismatch [14].Using designed dielectric stacks for cavity mirrors allows for the creation of phase shifts other than π radians.
The generated harmonic field is plotted as a function of cavity length for different values of mirror reflectivities in Fig. 3, using standard dispersion parameters for GaAs [11].The nonlinear coherence length for conversion of optical radiation from a wavelength of 2 micrometer (μm) to 1 μm in GaAs is approximately 6 μm.In the SHG process, energy transfers from the pump wave to the harmonic wave during propagation until half of a coherence length is reached.At that point the direction of energy flow is from the harmonic to the pump wave, due to the phase mismatch.The optimal cavity length is therefore half of a coherence length.For low mirror reflectivities the variation of harmonic conversion efficiency with cavity length is straightforward.For high finesse cavities, however, the reflected waves inside the cavity are large enough in magnitude to produce interference effects.This is seen by the added minima  in Fig. 3 for high mirror reflectivity.This effect has not been described in the previously-used calculation methods or periodic-poling experiments [14].The analogy of a cavity to a folded periodically-poled structure does not apply to a high-finesse cavity with strongly interfering waves.
The high-finesse case can be analyzed by considering the phase of the harmonic field photons being reflected at the cavity mirrors.The angle of the mirror reflection coefficients was set to kd for the pump and the harmonic in order to reset the phase of any beams that travel across the entire cavity.Photons of harmonic field generated within the cavity, however, have traveled a distance less than d when they reach the mirrors.A phase-compensation scheme that accounts for this is presented in the Appendix.A more effective optimization technique, as opposed to the simple phase-compensation scheme used above, is to search the space of mirror phase reflectivities for both frequencies.In Fig. 4, the phase of the reflection coefficients of the mirrors at the fundamental and the harmonic frequencies were both varied from 0 to 2π, and the maximum conversion achieved is plotted.The plots show that while a simple phase compensation method is adequate for cavity lengths of half of a coherence length, other cavity lengths can only be optimized by the analysis methods presented here.Such mirrors, with designed complex reflection coefficients at the fundamental and the harmonic wavelengths, could require over 20 layers, as described in [5,13].It should also be noted that as the finesse increases, at some point the slowly-varying envelope approximation will no longer be valid.For such cases the second derivative terms can be included in the solution.These added terms will affect the phase of the solution, so the same conversion can be achieved, but the required mirror reflection coefficients will be shifted in angle from the SVEA case.

Conclusion
We have presented a new method for nonlinear cavity calculations based on optimization methods.Similar to the harmonic balance method in nonlinear circuit analysis, a residual is calculated from the assumed field values, and the values are then adjusted to lower the residual.The calculations for each set of field values are made in the steady-state condition, where two counter-propagating waves exist in the cavity.The nonlinearity is represented by coupled differential equations acting on the steady-state waves.This method is general since any number of waves and nonlinear effects that can be characterized by differential equations can be modeled.Multiple cavities can be modeled with the addition of more unknowns.We have presented the example of this method applied to second-harmonic generation in a resonant cavity.Since this method finds the full phases of the fields exactly in the resonator, new effects due to interference can be described.

Appendix
We can find an example of an alternate phase-compensation scheme by considering the phase of the harmonic field photons being reflected at the cavity mirrors.As explained earlier, the angle of the mirror reflection coefficients was set to kd for the pump and the harmonic in order to reset the phase of any beams that travel across the entire cavity.Photons of harmonic field generated within the cavity, however, have traveled a distance less than d when they reach the mirrors.A photon of harmonic field generated at distance z within the cavity has phase term exp(− j2k 0 z), where k 0 is the propagation constant of the pump.When the mirror is reached, the phase term of that photon is exp[− j2k 0 z − jk 1 (d − z)], where k 1 is the harmonic propagation constant.This term can be integrated from z = 0 to z = d to find the the complex phase term of the total harmonic field at the mirror: Adjusting the mirror reflection coefficient of the mirrors to compensate for this phase term produces the harmonic conversion curves of Fig. 5.The minima that appear in Fig. 3 are removed by this compensation scheme.The conversion at most other cavity lengths is less than before, due to the effect of the new compensation scheme on the Fabry-Perot resonance at the harmonic frequency.Some optimization is therefore useful for each cavity length, as shown in Fig. 4.

Fig. 2 .
Fig. 2. Calculated harmonic intensity for a singly-resonant 1.5 micrometer cavity of GaAs.The pump intensity is 10 MW/m 2 .Shown are results for different values of cavity finesse at the pump frequency.

Fig. 3 .
Fig. 3. Calculated harmonic intensity for a cavity of GaAs with phase-compensating mirrors.The pump intensity is 10 MW/m 2 .Shown are results for different values of cavity finesse.

Fig. 4 .
Fig.4.Calculated harmonic intensity for a cavity of GaAs with phase-compensating mirrors.Parameters are the same as in the previous case, but unlike the simple phasecompensation technique, the mirror phases at the fundamental and the harmonic frequencies are searched for the optimum values at each cavity length.

Fig. 5 .
Fig.5.Calculated harmonic intensity for a cavity of GaAs with phase-compensating mirrors.The compensation is according to Equation13.One curve with the same phase compensation of Fig.3is included for comparison.The pump intensity is 10 MW/m 2 .Shown are results for different values of cavity finesse.