# Patent

## Title
Efficient terahertz sources based on difference-frequency generation in triply-resonant photonic resonators

## Abstract

A system for efficient generation of THz radiation is provided that includes a triply-resonant nonlinear photonic resonator coupled to at least one near-infrared (NIR) or optical waveguide and to at least one THz waveguide. The energy traveling through the at least one near-infrared (NIR) or optical waveguide is converted to THz radiation inside the triply-resonant photonic resonator via a nonlinear difference frequency generation (DFG) process.

## Background
PRIORITY INFORMATION
This application claims priority from provisional application Ser. No. 61/224,566 filed Jul. 10, 2010, which is incorporated herein by reference in its entirety.


GOVERNMENT SPONSORSHIP
This invention was made with government support under grant numbers DMR0819762 and PHY0646094 awarded by the National Science Foundation, and grant number W911NF-07-D-0004 awarded by the Army Research Office. The government has certain rights in this invention.


BACKGROUND OF THE INVENTION
The invention is related to terahertz sources, and in particular to efficient terahertz sources based on difference-frequency generation in triply-resonant photonic resonators.
Nonlinear optical frequency conversion is widely used for the generation of light in parts of the spectrum for which there are no convenient sources. In particular, nonlinear processes are regarded as a promising route to generation of coherent radiation in the THz frequency range. Recently there has been a renewed interest in using cavity-enhanced nonlinear frequency conversion to enable greater flexibility and THz powers, as improved designs and fabrication techniques have paved the way for the realization of wavelength-scale cavities, thus allowing efficient conversion at increasingly low powers.
SUMMARY OF THE INVENTION
According to one aspect of the invention, there is provided a system for efficient generation of THz radiation. The system includes a triply-resonant nonlinear photonic resonator coupled to at least one near-infrared (NIR) or optical waveguide and to at least one THz waveguide. The energy traveling through the at least one near-infrared (NIR) or optical waveguide is converted to THz radiation inside the triply-resonant photonic resonator via a nonlinear difference frequency generation (DFG) process.
According to another aspect of the invention, there is provided a method for generating efficient THz radiation. The method includes providing a triply-resonant nonlinear photonic resonator coupled to at least one near-infrared (NIR) or optical waveguide and to at least one THz waveguide. Also, the method includes converting the THz radiation inside the triply-resonant photonic resonator via a nonlinear difference frequency generation (DFG) process.



BRIEF DESCRIPTION OF THE DRAWINGS
FIGS. 1A-1B are schematic diagrams of two alternative configurations of the invention;
FIGS. 2A-2C are schematic diagrams illustrating the first physical implementation for efficient THz generation;
FIG. 3A-3F are schematic diagrams illustrating the second physical implementation for efficient THz generation;
FIG. 4 is a schematic diagram illustrating the CMT framework of an all-resonant DFG used in accordance with the invention;
FIGS. 5A-5C are graphs illustrating maximum and minimum normalized efficiency parameter (Eff) for stable a CW DFG;
FIGS. 6A-6B are graphs illustrating the effect of geometry on the conversion efficiency and stability of steady state solutions;
FIG. 7 is a graph illustrating the suggested trajectory in the space for THz generation in a single device;
FIG. 8 is a graph illustrating the conversion efficiencies for three frequencies involved in the nonlinear conversion process as a function of the total input power; and
FIG. 9A-9B are graphs illustrating dependence on the product of the cavity quality factors of the input power yielding optimal efficiency and the corresponding THz output power.



DETAILED DESCRIPTION OF THE INVENTION
The invention involves a system for efficient generation of THz via difference frequency generation (DFG) caused by a second order (X(2)) nonlinear polarization in a triply-resonant photonic structure. In contrast with most of previous works in this context, the case where both the input pump (highest frequency) and idler waves are resonant, as well as the generated signal wave at the difference frequency. The invention shows that continuous-wave (CW) conversion from pump to signal with quantum-limited efficiency (complete depletion of the pump) is possible in such a system for any power at the pump frequency, provided there is a critical input power at the idler frequency. The steady-state solutions to the coupled wave equations and their corresponding conversion efficiencies are shown to only depend on universal dimensionless parameters, indicating that their structure is qualitatively the same in all geometries. While the stability of these solutions is the only feature that is geometry dependent, one can demonstrate that in a certain region of the pump-idler power-space, there is only one steady-state solution to the coupled wave equations that is stable in all geometries. Within this range, quantum-limited conversion corresponds to the only stable solution for all triply resonant geometries, at a critical idler power. Outside of this region, it is shown that bistability exists. The effect of linear losses (e.g. scattering and absorption) results in a rescaling of the conversion efficiency and input powers, with the dynamics remaining qualitatively the same. The maximum efficiency is attained in the over-coupled limit, when the Q-factor is limited by coupling to the input/output port.
FIGS. 1A-1B show a schematic layout of two alternative configurations 2, 4 of the invention. In the configuration 2 displayed in FIG. 1A, the power carried by two optical or near-infrared beams 12 is transferred to a first photonic microresonator 6 by means of an optical or near-infrared waveguide 10 coupled evanescently to the microresonator 6. This first photonic microresonator 6 includes a material with second-order (X(2)) optical nonlinearity. Inside the microresonator 6, the optical or near-infrared energy is converted into THz energy by means by way of a nonlinear second-order (X(2)) nonlinear difference-frequency generation.
The efficiency of this process is enhanced dramatically when the whole system is designed so that the frequencies of the two optical or near-infrared beams 12 coincides with two resonant frequencies of the first photonic resonator 6 and, at the same time, the frequency of the emitted THz radiation coincides with the frequency of a nonlinear photonic resonator 8. Finally, the THz radiation is coupled to a THz waveguide that permits extracting efficiently the THz energy generated inside the first microresonator 6. The description of the configuration 4 shown in FIG. 1B is the same as FIG. 1A, but now there are two optical or near-infrared waveguides 14, 16 coupled to the microresonator 6; one of them acting as the input port 16 (right waveguide in the figure), while the other acting as a output port 14 for the optical or near-infrared beams 12.
The triply-resonant photonic resonator can include a nonlinear photonic resonator 8 for NIR or optical frequencies embedded in or placed in the proximity of the photonic resonator for THz frequencies, in such a way that if the NIR cavity is removed a THz cavity mode would be significantly modified. The nonlinear photonic resonator 8 for NIR or optical frequencies can include a doubly-resonant photonic resonator. The doubly-resonant photonic resonator comprises two resonant wavelengths λ1 and λ2 in air that support two resonant modes characterized by quality factors Q1 and Q2, respectively, and by modal volumes V1 and V2, respectively; with 400 nm<λ1<20 μm and 400 nm<λ2<20 μm. Also, the doubly-resonant photonic resonator can include a photonic crystal cavity or a longitudinally perturbed ridge waveguide or an undercut ridge waveguide. The photonic resonator for THz frequencies can include a resonant wavelength in air λT, a quality factor Q, and a modal volume VT; with 30 μm<λT<3000 μm. The photonic resonator can have a longitudinally perturbed ridge waveguide or undercut ridge waveguide or photonic crystal cavity. The THz waveguide can have a dielectric waveguide or a photonic crystal waveguide.
There are also various physical implementations of the invention. Several different physical implementations of the invention are discussed, as shown in FIGS. 2A-2F. These implementations include two nested cavities. The first cavity 24, designed to operate at telecom wavelengths and realized in nonlinear material (e.g. III-V), produces nonlinear-polarization at THz, as shown in FIG. 2A. The cavity 24 is placed in the proximity of a second 26, much larger, single-mode cavity with the resonance in THz. This implementation takes advantage of the vast difference in the lengthscales of the telecom/optical and THz cavities to allow the designs to be completely decoupled.
The system can be viewed as a classical version of a two-level atom strongly coupled to the optical cavity, with the dual mode telecom/visible cavity making a negligible perturbation on the THz mode. A proof-of-concept example of this type of design is explored here in detail and shown in FIG. 2A. Here, both the telecom and THz cavities are photonic crystal “nanobeam” cavities 22, 24. The telecom cavity has both TE-like and TM-like high-Q fundamental modes. One can suspend the telecom band cavity slightly above the THz cavity near the field maximum, sufficiently far that the telecom modes are not affected by the presence of the THz cavity, but sufficiently close that the long-wavelength THz mode extends sufficiently into this cavity. In this way one can achieve a high overlap between the three confined modes by mixing three mutually orthogonal polarizations through the strong X(2)ijk(i≠j≠k) terms in III/V semiconductors (e.g. GaAs, GaP). The invention takes advantage of recently demonstrated doubly-resonant ultra-high Q photonic crystal nanobeam resonators 22, 24 with mutually orthogonally polarized modes (i.e. one is TE-like one is TM-like). Both cavities are preferentially coupled to a corresponding waveguide extending from one end, with the length of the Bragg mirror of holes at one end used to tune the strength of the coupling. The Bragg mirror at the other end is made sufficiently long that leakage in that direction is much smaller than out of plane losses. An alternative design would be based on a telecom PhCNC embedded at the center of the THz cavity composed of different material (e.g. polymers, ceramics).
FIG. 2B shows nested 2-D photonic crystal cavities 28, 30. Telecom cavity supports 30 two high-Q closely spaced modes and is embedded inside of (or suspended just above) a much larger 2D photonic crystal cavity 28 supporting THz modes. The THz cavity 28 can also be a 1-D nanobeam cavity in this geometry. FIG. 2C shows another variant of this concept having a dual-mode nanobeam cavity 34 placed inside the gap of a metamaterial split ring resonator cavity 32.
The second physical implementation are disclosed FIGS. 3A-3F. In this type of design, the telecom cavity is designed as a component of the THz cavity, its presence having a considerable effect on the THz modes. The specific geometry 40 in FIG. 3A is examined in detail. The power carried by two near-infrared beams of wavelengths λ1=1542 nm and λ1=1550 nm (their corresponding power is denoted in FIG. 3A by P1i, and P2i, respectively) is coupled to a nonlinear ring 44 resonator by means of an index-guided waveguide 46. The width of the waveguide 46 is w=450 nm, while the external and internal radius of the ring 44 are Rext=40.1 μm and Rint=30.1 μm, respectively.
These values for and Rext and Rint were optimized so the ring resonator 44 can support two high-order whispering gallery modes (WGM) at ω1 and ω2, such that ωT=ω1−ω2. Specifically, the proposed system is created from GaAs, at λ1=1542 nm the ring resonator supports a TE-polarized WGM characterized by an angular moment m1=575, whereas for λ2=1550 nm there exits a TM-polarized WGM with m2=572. The inset 48 of FIG. 3B displays a contour plot of the product E1x(x, y)E*2y(x, y), where it is apparent the azimuthal dependence exp[i(m2−m1)θ]. As can be also deduced from the inset 48 of FIG. 3B, the high-order value of angular moment of the two considered WGM yields in both cases an electric field profile that is confined in a small region close to edge of the ring resonator. This is illustrated in bottom index of FIG. 3B, which renders the radial dependence of |(x, y)|, where the vertical red line corresponds to Rext.
The described ring resonator is in turn embedded as a dipole-like defect in an otherwise perfectly periodic THz-wavelength scale photonic crystal formed by a square lattice of GaAs rods of radius r=20.4 μm, the lattice constant being α=102 μm. The above described structure has been optimized for efficient emission at ωT/2π=1 THz. In order to extract efficiently the THz power generated inside the ring resonator, a coupled-cavity waveguide (CCW) is added to the THz-scale PhC. In particular, this CCW THz waveguide has been created by reducing to the radius of a row rods to rd1=10.2 μm. Furthermore, to avoid the coupling of the emitted THz radiation to the degenerate dipole-mode orthogonal to the one displayed in FIG. 3B (which can not be coupled to the THz CCW modes) the radius of the two of the nearest neighbors rods of the disk resonator to rd2=13.2 μm. FIGS. 3C-3F display alternative implementations of this scheme.
In particular, FIG. 3C shows a monopole defect mode 50 in a THz-scale nonlinear photonic crystal 52 is phase-matched to two whispering gallery modes (WGM) circulating inside the defect rod by introducing a periodic poling in the defect rod. FIG. 3D shows another implementation 54 having a one-dimensional (1D) periodic array of nonlinear rods 56, in which the central rod has been replaced by a ring resonator 58 supporting two WGMs at the idler and pump frequencies. The nonlinear interaction that takes place between both WGMs produces a current distribution that radiates THz energy. The emitted THz energy is collected by a THz waveguide 60 placed in the proximity of the PhC configuration 56 and underneath the PhC configuration 56 is a NIR or optical waveguide 62.
FIG. 3E show a schematics diagram of another implementation 66 based on a triply-resonant photonic crystal configuration. The power carried by two optical or NIR beams traveling through an optical or NIR waveguide 68 (which, in turn, is embedded in a THz waveguide 74) gets trapped inside a 1D doubly-resonant asymmetric PhC cavity 72. This cavity 72 is embedded in a THz-scale cavity 70 created by reducing the radius of the central rod of two 1D arrays of dielectric rods placed at each side of the THz waveguide 74.
FIG. 3F shows another implementation 80 for efficient THz generation. A circular array of nonlinear dielectric cylinders 84 is placed close to a THz waveguide 84. The THz radiation is generated inside each cylinder 88 by nonlinear frequency mixing of the radiation emitted by two optical or NIR emitters 86 embedded in each one of the cylinders. The evanescent coupling between the THz-scale cylinders 88 gives rise to a waveguide confined mode that circulates inside the cylinder array 84 and decays into the THz waveguide 82.
The general framework 94 of the temporal CMT model is shown in FIG. 4. The cavity has three resonant modes at frequencies, ω1 (pump), ω2 (idler), ωT (signal) which satisfy ωT=ω1−ω2. The electric field in each cavity mode, Ek(r,t), is decomposed into Ek(r,t)=Ek,0(r)αk(t)exp(iωkt), where Ek,0 are the normalized mode profiles (∫d3rε|Ek,0|2=1) and αk(t) are the slowly varying wave amplitudes, normalized so that |αk|2 is the energy stored in each resonant mode. Energy is coupled in/out of each cavity mode through a designated input/output port (e.g. waveguide). The propagating modes in the input/output ports are similarly represented. The relevant incoming and outgoing wave amplitudes are represented by the variables sk+(t) and sk−(t) respectively, where |sk±(t)|2 is the power of the propagating mode. Leakage from a cavity mode is described by the loss rate, γk=ωk/2Qk, where Qk is the Quality factor of that mode. One can decompose this cavity loss rate into γk=γk,s+γk,l, where, γk,1 represents intrinsic losses (absorption/scattering), while γk,s represents coupling to the input/output port. The incoming and outgoing waves are related by

sk−=−sk++√{square root over (2γk,s)}αk   (1)

The invention requires the assumption that the nonlinear polarizations are sufficiently weak to be treated in the perturbative limit (δPk,NL<<Pk). By extending the formalism to non-degenerate three-wave mixing, where one can assume there is no input at the signal frequency (sT+=0), the coupled wave equations can be written as





  
    
      
        
          
            d
            
              a
              1
            
          
          
            d
            t
          
        
        =
        
          
            
              -
              
                γ
                1
              
            
            ⁢
            
              a
              1
            
          
          -
          
            i
            ⁢
            
              
            
            ⁢
            
              ω
              1
            
            ⁢
            
              β
              1
            
            ⁢
            
              a
              2
            
            ⁢
            
              a
              T
            
          
          +
          
            
              
                2
                ⁢
                
                  γ
                  1
                
                ⁢
                
                  Γ
                  1
                
              
            
            ⁢
            
              s
              1
              +
            
          
        
      
    
    
      
        (
        2
        )
      
    
  
  
    
      
        
          
            d
            
              a
              2
            
          
          
            d
            t
          
        
        =
        
          
            
              -
              
                γ
                2
              
            
            ⁢
            
              a
              2
            
          
          -
          
            
              iω
              2
            
            ⁢
            
              β
              2
            
            ⁢
            
              a
              1
            
            ⁢
            
              a
              T
              *
            
          
          +
          
            
              
                2
                ⁢
                
                  γ
                  2
                
                ⁢
                
                  Γ
                  2
                
              
            
            ⁢
            
              s
              2
              +
            
          
        
      
    
    
      
        (
        3
        )
      
    
  
  
    
      
        
          
            d
            
              a
              T
            
          
          
            d
            t
          
        
        =
        
          
            
              -
              
                γ
                T
              
            
            ⁢
            
              a
              T
            
          
          -
          
            i
            ⁢
            
              
            
            ⁢
            
              ω
              T
            
            ⁢
            
              β
              T
            
            ⁢
            
              a
              1
            
            ⁢
            
              a
              2
              *
            
          
        
      
    
    
      
        (
        4
        )
      
    
  




where Γk≡γk,s/γk is the fraction of power escaping from the cavity that is collected into the input/output waveguide. The coupling constants (βk) are determined by first order perturbation theory, where one can also normalize so that |αk|2 is the energy stored in each mode, remaining consistent with Eqs. (1-4). This yields the following expression for the coupling coefficients:





  
    
      
        
          β
          ≡
          
            β
            1
          
        
        =
        
          
            β
            2
            *
          
          =
          
            
              β
              T
              *
            
            =
            
              
                1
                4
              
              ⁢
              
                
                  ∫
                  
                    
                      
                        d
                        3
                      
                      ⁢
                      
                        x
                        ρ
                      
                    
                    ⁢
                    
                      
                        ∑
                        
                          i
                          ,
                          j
                          ,
                          k
                        
                      
                      ⁢
                      
                        
                      
                      ⁢
                      
                        
                          ɛ
                          0
                        
                        ⁢
                        
                          χ
                          ijk
                          
                            (
                            2
                            )
                          
                        
                        ⁢
                        
                          
                            E
                            
                              1
                              ,
                              i
                            
                            *
                          
                          ⁡
                          
                            (
                            
                              
                                
                                  E
                                  
                                    2
                                    ,
                                    j
                                  
                                
                                ⁢
                                
                                  E
                                  
                                    T
                                    ,
                                    k
                                  
                                
                              
                              +
                              
                                
                                  E
                                  
                                    T
                                    ,
                                    j
                                  
                                
                                ⁢
                                
                                  E
                                  
                                    2
                                    ,
                                    k
                                  
                                
                              
                            
                            )
                          
                        
                      
                    
                  
                
                
                  
                    
                      
                        ∫
                        
                          
                            
                              d
                              3
                            
                            ⁢
                            
                              x
                              ρ
                            
                          
                          ⁢
                          ɛ
                        
                      
                      ❘
                      
                        
                          
                            E
                            ρ
                          
                          1
                        
                        ⁢
                        
                          ❘
                          2
                        
                      
                    
                  
                  ⁢
                  
                    
                      
                        ∫
                        
                          
                            
                              d
                              3
                            
                            ⁢
                            
                              x
                              ρ
                            
                          
                          ⁢
                          ɛ
                        
                      
                      ❘
                      
                        
                          
                            E
                            ρ
                          
                          2
                        
                        ⁢
                        
                          ❘
                          2
                        
                      
                    
                  
                  ⁢
                  
                    
                      
                        ∫
                        
                          
                            
                              d
                              3
                            
                            ⁢
                            
                              x
                              ρ
                            
                          
                          ⁢
                          ɛ
                        
                      
                      ❘
                      
                        
                          
                            E
                            ρ
                          
                          T
                        
                        ⁢
                        
                          ❘
                          2
                        
                      
                    
                  
                
              
            
          
        
      
    
    
      
        (
        5
        )
      
    
  




where ε=εrε0 is the dielectric constant of the considered material at the corresponding resonance frequency. The above equations (subject to Kleimann symmetry in X(2)) satisfy the energy conservation condition, β1=β2*=βT*. Hereafter the parameters, βk, will be replaced with the single parameter, β, which is defined in Eq. (5), with β represented in units of J−1/2. Optimization of the mode overlap, β, in Eq. (5) is analogous to satisfying the phase-matching condition in the nonlinear mixing of propagating modes.

DFG efficiency can be assessed by analyzing the steady state solutions (dαk(t)/dt=dsk±(t)/dt=0) of the coupled wave equations (2-4). In contrast with degenerate nonlinear frequency conversion (e.g. second harmonic generation) [21], the non-degenerate case (e.g. DFG) does not allow for 100% conversion of power from ω1 to ωT. This arises because the destruction of a photon at ω1 and the creation of a photon at ωT is always accompanied by the creation of a second photon at ω2. This fundamental constraint is described by the Manley-Rowe equations [1]. Considering also losses (Γk), the quantum limit of conversion efficiency is given by:




  
    
      
        
          
            
              
              
                s
                T
                -
              
              
            
            2
          
          
            
              
              
                s
                1
                +
              
              
            
            2
          
        
        ≤
        
          
            
              ω
              T
            
            
              ω
              1
            
          
          ⁢
          
            Γ
            1
          
          ⁢
          
            Γ
            T
          
        
      
    
    
      
        (
        6
        )
      
    
  




This condition can be derived directly from the steady-state coupled wave equations, and will be shown once Eqs. (2-4) are re-expressed in dimensionless form. The conversion efficiency is maximized by maximizing the ratios, Γ1 and ΓT. This means that overcoupling the cavity modes to the input/output waveguide is essential to achieving high conversion efficiency. Namely, it is more important that the quality factors be limited by leakage to the desired output port than that they be high. However, notice that the loss ratio in the idler, Γ2, does not appear in Eq. (6), i.e. does not affect the overall efficiency. This is because power at ω2 is being input but is not down-converted; therefore losses can be compensated for by pumping in more power at this frequency. One can find a critical relationship between the input powers, |s1+|2 and |s2+|2 that allows for maximum conversion efficiency (complete depletion of the pump) by imposing either of the following equivalent constraints: s1−=0 or equality in Eq. (6). This relationship is given by:





  
    
      
        
          
            
            
              s
              2
              +
            
            
          
          2
        
        =
        
          
            
              ω
              2
            
            ⁢
            
              
                
                
                  1
                  -
                  
                    
                      4
                      
                        ω
                        1
                      
                    
                    ⁢
                    
                      
                        
                        β
                        
                      
                      2
                    
                    ⁢
                    
                      Q
                      1
                    
                    ⁢
                    
                      Q
                      2
                    
                    ⁢
                    
                      Q
                      T
                    
                    ⁢
                    
                      Γ
                      1
                    
                    ⁢
                    
                      
                        
                        
                          s
                          1
                          +
                        
                        
                      
                      2
                    
                  
                
                
              
              2
            
          
          
            16
            ⁢
            
              
                
                β
                
              
              2
            
            ⁢
            
              Q
              1
            
            ⁢
            
              Q
              2
            
            ⁢
            
              Q
              T
            
            ⁢
            
              Γ
              2
            
          
        
      
    
    
      
        (
        7
        )
      
    
  




Notice that the power input into the cavity at the idler frequency, |s2+|2, depletes the signal (ωT) when nonlinearly converted (via sum-frequency generation), and produces power at the pump (ω1). Therefore, in the case in which the pump is completely depleted (s1−=0), no net power from the idler frequency is converted and |s2−|2=|s2+|2+|s1+|2Γ1Γ2ω2/ω1. The conversion of the total input power, |s1+|2+|s2+|2, is thus maximized in the limit |s2+|2→0, |s1+|2=ω1/(4|β|2Q1Q2QTΓ1). Note, however, that Eqs (2-4) require a non-zero value (this can be arbitrarily small) of |s2+|2 for non-zero conversion efficiency. Eq. (7) also allows for quantum-limited conversion for an arbitrarily small value of |s1+|2 as |s2+|2→ω2/(16|β|2Q1Q2QTΓ2). This is a unique property of the non-degeneracy of DFG, and could prove useful for low power applications. A similar result was found in the quantum limit for doubly resonant DFG, specifically that a single photon could be converted at high efficiency from ω1 to ωT, given an appropriate pump power at ω2. Note that in both limits the optimal input power is proportional to 1/Q1Q2QT. A similar result was found in the degenerate case of second harmonic generation, where the critical power scaled as 1/Q1Q22. This optimal power also depends on the cavity mode-volumes implicitly through the nonlinear overlap, β. In the specific case of THz generation through DFG of telecom-band modes (ω1≈ω2), assuming the two telecom-band modes have the same mode volume (V1≈V2), then under the conditions of best possible overlap, the overlap would scale as β ̃1/VT1/2. This means that in triply resonant THz generation, the mode volumes of the telecom-band cavities do not affect the strength of the nonlinear coupling, provided they match each other.
The structure of Eqs. (2-4) is analyzed in detail. The dynamics that are qualitatively universal to all triply resonant DFG processes is to be identified, i.e. that do not depend explicitly on {Qk,ωk,β}. To do this, it is helpful to redefine the variables {αk, sk±} in a dimensionless form that removes as many of the parameters, {Qk,ωk,β}, as possible. This is best accomplished by resealing the variables, {t, αk, sk±}, as follows:




  
    
      
        
          T
          ≡
          
            
              γ
              1
            
            ⁢
            t
          
        
        =
        
          
            
              ω
              1
            
            ⁢
            t
          
          
            2
            ⁢
            
              Q
              1
            
          
        
      
    
    
      
        (
        8
        )
      
    
  
  
    
      
        
          A
          k
        
        ≡
        
          2
          ⁢
          
            
              
                
                  
                    Q
                    1
                  
                  ⁢
                  
                    Q
                    2
                  
                  ⁢
                  
                    Q
                    T
                  
                
                
                  Q
                  k
                
              
            
            ⁡
            
              [
              
                
                  β
                  *
                
                +
                
                  
                    δ
                    
                      k
                      ,
                      T
                    
                  
                  ⁡
                  
                    (
                    
                      
                        
                          -
                          i
                        
                        ⁢
                        
                          
                        
                        ⁢
                        β
                      
                      -
                      
                        β
                        *
                      
                    
                    )
                  
                
              
              ]
            
          
          ⁢
          
            a
            k
          
        
      
    
    
      
        (
        9
        )
      
    
  
  
    
      
        
          S
          k
        
        ≡
        
          4
          ⁢
          
            
              
                
                  Q
                  1
                
                ⁢
                
                  Q
                  2
                
                ⁢
                
                  Q
                  T
                
                ⁢
                
                  Γ
                  k
                
              
              
                ω
                k
              
            
          
          ⁢
          
            β
            *
          
          ⁢
          
            s
            k
            +
          
        
      
    
    
      
        (
        10
        )
      
    
  




where δk,T is the Kronecker delta. One can also define a dimensionless parameter that quantifies the DFG efficiency:





  
    
      
        
          
            E
            ff
          
          ≡
          
            
              
                ω
                1
              
              
                
                  Γ
                  1
                
                ⁢
                
                  Γ
                  T
                
                ⁢
                
                  ω
                  T
                
              
            
            ⁢
            
              
                
                  
                  
                    s
                    T
                    -
                  
                  
                
                2
              
              
                
                  
                  
                    s
                    1
                    +
                  
                  
                
                2
              
            
          
        
        =
        
          4
          ⁢
          
            
              U
              T
            
            
              P
              1
            
          
        
      
    
    
      
        (
        11
        )
      
    
  




where the notation, |Ak|2≡Uk, |Sk|2≡Pk is introduced for the normalized energy and power, respectively. In this notation, the quantum limit corresponds to Eff=1. The dimensionless coupled mode equations are given by





  
    
      
        
          
            d
            
              A
              1
            
          
          
            d
            T
          
        
        =
        
          
            -
            
              A
              1
            
          
          +
          
            
              A
              2
            
            ⁢
            
              A
              T
            
          
          +
          
            S
            1
          
        
      
    
    
      
        (
        12
        )
      
    
  
  
    
      
        
          
            d
            
              A
              2
            
          
          
            d
            T
          
        
        =
        
          
            r
            2
          
          ⁡
          
            [
            
              
                -
                
                  A
                  2
                
              
              -
              
                
                  A
                  1
                
                ⁢
                
                  A
                  T
                  *
                
              
              +
              
                S
                2
              
            
            ]
          
        
      
    
    
      
        (
        13
        )
      
    
  
  
    
      
        
          
            d
            
              A
              T
            
          
          
            d
            T
          
        
        =
        
          
            r
            T
          
          ⁡
          
            [
            
              
                -
                
                  A
                  T
                
              
              -
              
                
                  A
                  1
                
                ⁢
                
                  A
                  2
                  *
                
              
            
            ]
          
        
      
    
    
      
        (
        14
        )
      
    
  




where rk=γk/γ1. Since S1 and S2 denote amplitudes of waves with different frequencies, the time origin can be redefined so that both S1 and S2 are positive real numbers (i.e. the waves are in phase at T=0), assuming that the waves are not phase modulated. It follows that the normalized mode amplitudes, A1, A2 and AT are purely real in any steady state solution. In the steady state, Eqs. (12-14) simplify to,





  
    
      
        
          A
          1
        
        =
        
          
            S
            1
          
          
            (
            
              1
              +
              
                U
                2
              
            
            )
          
        
      
    
    
      
        (
        15
        )
      
    
  
  
    
      
        
          A
          2
        
        =
        
          
            S
            2
          
          
            (
            
              1
              -
              
                U
                1
              
            
            )
          
        
      
    
    
      
        (
        16
        )
      
    
  




AT=−A1A2   (17)

Note that, the steady state solutions (dAk/dT=0) do not depend on the geometry-specific parameters, rk, and thus can be treated in a general form for an arbitrary system. Specifically the conditions for efficient conversion will depend only on the generalized input parameters, Sk, whose magnitude is related to the input power by Eq. (10). However the stability of the steady-state solutions may in general depend on the specific geometry (rk), as shown below. An analogous simplification was also found in third harmonic generation (X(3)) in doubly resonant systems, however in this case the competition between frequency conversion and self-action processes (self-phase modulation and cross-phase modulation), prevented the structure of the steady-state solutions from being completely geometry-independent. Using Eqs. (15-17), the dimensionless efficiency parameter can be expressed in the steady state as purely a function of the energy in the idler cavity, U2:





  
    
      
        
          E
          ff
        
        =
        
          
            4
            ⁢
            
              U
              2
            
          
          
            
              (
              
                1
                +
                
                  U
                  2
                
              
              )
            
            2
          
        
      
    
    
      
        (
        18
        )
      
    
  




One can now easily verify that the quantum limit is enforced by the coupled wave equations (2-4), that is Eff≦1 in the steady state. Further, the quantum-limited conversion is achieved when





  
    
      
        
          
            U
            2
          
          =
          1
        
        ,
        
          
            U
            1
          
          =
          
            
              U
              T
            
            =
            
              
                P
                1
              
              4
            
          
        
        ,
        
          
            P
            2
          
          =
          
            
              (
              
                1
                -
                
                  
                    P
                    1
                  
                  4
                
              
              )
            
            2
          
        
      
    
    
      
        (
        19
        )
      
    
  




Stability of steady state solutions in a nonlinear dynamical system can be analyzed qualitatively based on the eigenvalues of the points' Jacobian matrices. A steady state is stable if all the eigenvalues of its Jacobian have a negative real part. The Jacobian can then be written as:





  
    
      
        J
        =
        
          [
          
            
              
                
                  -
                  1
                
              
              
                
                  A
                  T
                
              
              
                
                  A
                  2
                
              
            
            
              
                
                  
                    -
                    
                      r
                      2
                    
                  
                  ⁢
                  
                    A
                    T
                  
                
              
              
                
                  -
                  
                    r
                    2
                  
                
              
              
                
                  
                    -
                    
                      r
                      2
                    
                  
                  ⁢
                  
                    A
                    1
                  
                
              
            
            
              
                
                  
                    -
                    
                      r
                      T
                    
                  
                  ⁢
                  
                    A
                    2
                  
                
              
              
                
                  
                    -
                    
                      r
                      T
                    
                  
                  ⁢
                  
                    A
                    1
                  
                
              
              
                
                  -
                  
                    r
                    T
                  
                
              
            
          
          ]
        
      
    
    
      
        (
        20
        )
      
    
  




The stable solutions with maximum- and minimum-efficiency are shown in FIGS. 5A-5B respectively where steady state solutions for r2=rT=0.4. The dotted line 97 in FIGS. 5A-5B indicates a saddle-node bifurcation. In the region to the left of this line 97 (lower values of P1) there is only one steady state solution. To the right of this line 97, there are three possible steady states for every combination of CW input powers, two of which are stable.
Given that there exists a large region of the {P1, P2} parameter space with only one solution to Eqs. (15-17), one would expect there to exist a subset of this region where this solution is stable for any geometry {r2, rT}. In fact, it can be shown that in the entire region this solution is stable in all geometries (in the region where there is only one solution). This can be proven by analyzing the Jacobian's characteristic polynomial using the Routh-Hurwitz algorithm. The characteristic polynomial of Eq. (20) can be written in the form λ3+B λ2+Cλ+D=0. Using Eqs. (15-17) to simplify, one can find that:

B=1+r2+rT   (21)





  
    
      
        C
        =
        
          
            
              r
              2
            
            ⁢
            
              
                r
                T
              
              ⁡
              
                (
                
                  1
                  -
                  
                    U
                    1
                  
                
                )
              
            
          
          +
          
            
              r
              2
            
            ⁡
            
              (
              
                1
                +
                
                  
                    1
                    4
                  
                  ⁢
                  
                    E
                    ff
                  
                  ⁢
                  
                    P
                    1
                  
                
              
              )
            
          
          +
          
            
              r
              T
            
            (
            
              
            
            ⁢
            
              1
              +
              
                U
                2
              
            
            )
          
        
      
    
    
      
        (
        22
        )
      
    
  
  
    
      
        D
        =
        
          
            r
            2
          
          ⁢
          
            
              r
              T
            
            ⁡
            
              (
              
                1
                -
                
                  U
                  1
                
                +
                
                  U
                  2
                
                +
                
                  
                    3
                    4
                  
                  ⁢
                  
                    E
                    ff
                  
                  ⁢
                  
                    P
                    1
                  
                
              
              )
            
          
        
      
    
    
      
        (
        23
        )
      
    
  




Note that all variables in Eqs. (21-23) can only have real positive values. The signs of the real parts of the Jacobian's eigenvalues can be determined by analyzing the first column of the Routh-Hurwitz matrix. If all entries in the column have the same sign then there are no eigenvalues with positive real parts. For the Jacobian this column is given by,





  
    
      
        
          {
          
            H
            
              i
              ,
              1
            
          
          }
        
        =
        
          [
          
            
              
                1
              
            
            
              
                B
              
            
            
              
                
                  
                    1
                    B
                  
                  ⁢
                  
                    (
                    
                      BC
                      -
                      D
                    
                    )
                  
                
              
            
            
              
                D
              
            
          
          ]
        
      
    
    
      
        (
        24
        )
      
    
  




The additional constraint that no eigenvalues have a real part equal to zero imposes the further constraints, (BC−D), D≠0. If one requires U1≦1, all the coefficients (B-D) are strictly positive. In this case, both constraints are reduced to a single condition, (BC−D>0. The term, BC−D, can be expressed as





  
    
      
        
          BC
          -
          D
        
        =
        
          
            
              r
              2
            
            ⁢
            
              
                r
                T
              
              ⁡
              
                [
                
                  
                    
                      (
                      
                        
                          r
                          2
                        
                        +
                        
                          r
                          T
                        
                      
                      )
                    
                    ⁢
                    
                      (
                      
                        1
                        -
                        
                          U
                          1
                        
                      
                      )
                    
                  
                  +
                  
                    2
                    ⁢
                    
                      (
                      
                        1
                        -
                        
                          
                            
                              E
                              ff
                            
                            ⁢
                            
                              P
                              1
                            
                          
                          4
                        
                      
                      )
                    
                  
                
                ]
              
            
          
          +
          
            
              (
              
                1
                +
                
                  
                    
                      E
                      ff
                    
                    ⁢
                    
                      P
                      1
                    
                  
                  4
                
              
              )
            
            ⁢
            
              (
              
                1
                +
                
                  r
                  2
                
              
              )
            
            ⁢
            
              r
              2
            
          
          +
          
            
              (
              
                1
                +
                
                  U
                  2
                
              
              )
            
            ⁢
            
              (
              
                1
                +
                
                  r
                  T
                
              
              )
            
            ⁢
            
              r
              T
            
          
        
      
    
    
      
        (
        25
        )
      
    
  




Therefore, a steady-state solution is stable in all geometries when U1<1, P1≦4/E. It follows from Eq. (15) that all steady states are stable when P1≦1. In fact numerically it is found that these conditions are satisfied for the steady state solution everywhere in the region of mono-stability (FIGS. 5A-5B, to the left of the dotted line 97). Thus, one can expect universal dynamics in this region of {P1, P2} space, with the conversion efficiency stabilizing to a value determined purely by the input powers, P1, P2, shown in FIGS. 5A-5B. Furthermore, for all P1<4, the solution with quantum-limited efficiency given in Eq. (19), is universally stable (for all r2 and rT}). For Eq. (19) with P1>4, stability is guaranteed by (21-25) when rT<0.414 (21/2−1).
FIGS. 6A-6B demonstrate how the stability of the high efficiency solutions varies with the geometric parameters, r2 and rT, for P1>4 when P1 and P2 are related by Eq. (19). When r2 and rT are both large, there remains only one stable solution (shown by the black curve in FIGS. 6A-6B) for all P1, which has quantum-limited efficiency for P1<4.
Stability is also independent of the parameters, {r2, rT}, in the limit P2<<1. It is expected that there will be multiple steady states for P1>1. One can analytically determine the efficiencies and stabilities of all steady state solutions in this limit. For any value of P1≠1, there is one solution where U1→P1, U2→0, and Eff→0 (see Eqs. (15-17)). This will be the only solution for P1<1. In particular, it follows directly from Eqs. (21-25) that this solution is universally stable for P1<1 and universally unstable for P1>1. For P1>1, there is an additional pair of solutions in which U1→1+, and U1→1−. These solutions differ only by the sign of A2 (phase), both having U2=P11/2−1, and Eff=4(P11/2−1)/P1. It can also be shown directly from Eqs. (21-25) that both of these solutions are stable for all geometries. This means that to a good approximation, all geometries display only one stable conversion efficiency as a function of P1 for P2<<1.
In addition to verifying the stability of high conversion-efficiency solutions, it is important to assess how one can excite these solutions. In the monostable region, the efficiency will stabilize at the value shown in FIGS. 5A-5B for all excitations, however different excitation patterns can lead to different long-time stabilities in the bistable region. This principle was demonstrated for multi-stable conversion in doubly resonant third harmonic generation. FIG. 5C shows the asymptotic solution reached after a simple step excitation (Uk(T=0)=0, Pk(T≧0)=constant, Pk(T<0)=0, γ2=γT=0.4). For P1<4, quantum limited conversion is achieved from a step excitation when P2 is given by Eq. (19), however this is not true for P1>4. A more complicated excitation scheme would be required to excite the high-efficiency solution in this case. The particular class of “easily accessible” steady-state solutions shown in FIG. 5C has a slightly increased domain of universal stability, U1<1, P1≦4/Eff, indicated by the region to the left of the dashed line 98 in FIG. 5C, which extends beyond the region of monostability. In particular, universal stability is guaranteed in this set of solutions for all P1<4.
It is important for many applications that a single cavity system (fixed Qk, ωk, β) is capable of operating efficiently over a broad range of power levels. Operation of a single system is optimized when the pump and idler powers follow a trajectory that optimizes the total power conversion (|sT−|2/[|s1+|2+|s2+|2]). In general this trajectory can be transcribed onto the plot in FIG. 7 in normalized co-ordinates, however this curve is not unique in the space of normalized powers {P1, P2} as they will depend on the ratio of the pump and idler frequencies ω1/ω2. For THz generation based on the difference frequency of telecom or optical modes it is safe to make the assumption ω1/ω2≈1. FIG. 7 is an expanded view of FIG. 5A showing on the dashed line 100 the optimized trajectory for single cavity operation in the case of THz generation. In the monostable region, this trajectory optimizes total power conversion. In the multistable region ease of excitation and stability are considered. While for P1>4, it is possible to excite slightly higher efficiency solutions ( ̃10%), these solutions are bistable and a composite excitation is required to excite these solutions (lower efficiency solutions are more easily excited). These high efficiency solutions are also only stable for a limited range of Q-factor ratios (r2, rT). Beyond P1=4, keeping the idler power low (P2<<1) ensures high efficiency monostable conversion.
Now the regions of the dimensionless parameter space are accessible for a given geometry (Qk, ωk, β). Sufficient resonance enhancement of the density of optical modes at each of the three frequencies (Qk>>1) is required in order to effectively suppress coupling to other frequencies (e.g. the sum frequency, ω1+ω2). CMT relies on the slowly-varying amplitude approximation (|dαk/dt|<<|ωkωk|). This condition can be mapped into the dimensionless parameter space, imposing the following conditions on the cavity quality factors:




  
    
      
        
          
            
              U
              1
            
            
              2
              ⁢
              
                Q
                2
              
            
          
          ⪡
          1
        
        ,
        
          
            
              U
              2
            
            
              Q
              1
            
          
          ⪡
          1
        
      
    
    
      
        (
        26
        )
      
    
  




In the case that Eff=1, Eq. (26) reduces to,

P1<<8Q2   (27)

When Eqs. (26-27) are not satisfied, it implies that the nonlinear terms in Eqs. (2-4) are sufficiently large that the perturbative treatment of the nonlinearity and the slowly varying amplitude approximation are no longer valid. A final consideration that must be made is that the energy in the cavity is sufficiently below the threshold for the appearance of higher order regimes or material damage. While Eqs. (26-27) ensure that effects owing to X(2) remain in the perturbative limit, the nonlinear coupling strength (quantified by β) must also be sufficiently large for high efficiency conversion to be achievable below the damage threshold. The average energy density in the kth mode can be related to the corresponding dimensionless energy parameter (Uk) by,





  
    
      
        
          ξ
          k
        
        =
        
          
            
              U
              k
            
            ⁢
            
              Q
              k
            
          
          
            4
            ⁢
            
              Q
              1
            
            ⁢
            
              Q
              2
            
            ⁢
            
              Q
              T
            
            ⁢
            
              
                
                β
                
              
              2
            
            ⁢
            
              V
              k
            
          
        
      
    
    
      
        (
        28
        )
      
    
  




where Vk is the mode volume at ωk. Recall that |β|2 has units of inverse energy. This energy density can be compared to the material and frequency dependent threshold for the appearance of higher order nonlinear effects or material damage (breakdown) for a specific geometry.

The above described theory is applied to analyze efficient THz generation by the first physical implementation disclosed here. This implementation makes use of the orders-of-magnitude difference in wavelengths of THz and visible/NIR radiation to allow independent design of the NIR and THz cavities. The key requirement of such a design is that the NIR cavity be sufficiently small in comparison to the THz cavity that the THz mode is not significantly altered by placing the NIR cavity nearby. One can examine in detail a specific implementation of this concept, shown in FIG. 2A, where the NIR cavity is a dual-mode photonic crystal nanobeam cavity placed in the proximity of a perpendicularly oriented THz nanobeam cavity.
In order to estimate the coupling constant, β, a geometry similar to that shown in FIG. 2A is considered: a photonic crystal “nanobeam” (PhCNB) cavity 24, with closely-spaced TE and TM modes in the telecom spectral region (ω1/2π, ω2/2π ̃200 THz) positioned on top of larger THz cavity 26. These cavities, both telecom and THz, are based on a free-standing ridge waveguide patterned with a one-dimensional (1D) lattice of holes, utilizing a tapered photonic crystal cavity approach that has been recently demonstrated. By engineering the impedance matching between the cavity mode and evanescent Bloch-mode that exists in the photonic crystal (Bragg) mirror demonstrated, both experimentally and theoretically, Quality factors on the order of 106.
Moreover, the inventive monolithic cavity design supports orthogonal TE and TM modes, which is necessary to achieve efficient nonlinear conversion as dictated by the nonlinear tensor of a III-V semiconductor material. A key design point is that the TE and TM cavity modes can be tuned independently by varying the cross-sectional aspect ratio of the ridge. As a practical implementation, a PhCNB cavity is designed with fundamental TE mode at 193.9 THz (V=0.6(λ/n)3) and TM mode at 195 THz (V=1.34 (λ/n)3). A finite-difference time-domain code (FDTD) was used to independently design telecom and THz cavities, and a CMT model was used to describe the nonlinear frequency conversion. The Q factors in the cavities are designed to be limited by out-of plane scattering. Then, on one side, the mirror (number of holes) is shortened in order to over-couple the cavity to the input/output waveguide, in order to strongly couple the cavity to desired channel which in turn reduces the maximum values of Q by an order of magnitude. This over-coupling minimizes the reduction in efficiency caused by losses, and enables operation close to quantum-limited conversion.
The case where both cavities are composed of GaAs (n≈3.4) was considered, and the telecom cavity is suspended slightly above the THz cavity at height, d<<λT (again ωT/2π ̃1.5 THz). In this geometry one can find a coupling constant, β≈1.5 W−1/2. The cavity Q-factors of all three modes can be optimized so that maximum conversion efficiency (normalized S12=4) occurs at the desired power level. Moreover, one can find that small adjustments in the geometry (e.g. over-coupling) allows for large tuning of the Q-factors with minimal effect on the coupling constant, β. Higher Q-factors allow for optimal conversion at lower powers.
The decoupled nature of the length-scales used in this design allows complete freedom in the design of the THz cavity—any linearly polarized cavity mode in the THz range can be coupled to a nearby dual-mode nanobeam cavity with the appropriate orientation. For example, the THz cavity could be a ring resonator, micro-toroid, 2D photonic crystal, or a metamaterial-based cavity. The dual mode cavity can also alternatively be a 2D photonic crystal cavity, which has also been shown to support multiple high-Q modes with different polarization.
The above described theoretical framework can be applied to analyze in detail the generation of THz that takes place in the structure displayed in FIG. 3A. To gain insight into the optimization of the nonlinear coupling coefficient for this design, one can start by rewriting in cylindrical coordinates 49 the expression for the relevant component of the nonlinear polarization vector. Using the definition for the coordinates (r, θ) 49 shown in top inset of FIG. 3B, and assuming that the structure is created from GaAs, one finds,

Pz(r, θ)=2 ε0 d14(r)[E1rE*2θcos2 θ−E1θE*2r sin2 θ+(E1rE*2r−E1θE*2θ) sin θ cos θ]  (29)

where Ei,r and Ei,θ stand for the radial and azimuthal electric field components of the two TE-polarized WGM modes considered in the scheme. Here the subindex i labels the corresponding WGM mode. Also note that for the structure displayed in FIG. 3B, d14(r) is non zero only in the region r<R.

Now, inserting in Eq. (29) the explicit dependence of Ei,r and Ei,θ on the azimuthal coordinate θ, given by Ei,r(r,θ)=φi,r(r)exp(−i mi θ) and Ei,θ(r,θ)=φi,θ(r)exp(i mi θ) (where mi is an integer), Pz(r, θ) can be written as

Pz(r,θ)=ε0 d14(r)exp[(m2−m1)θ] [h1(r) sin 2θ−4 h2(r) sin2 θ+h3(r)]  (30)

where one can define

h1(r)=φ1,r(r)φ*2,r(r)−φ1,θ(r)φ*2,θ(r)   (31)

h2(r)=φ1,r(r)φ*2,θ(r)−φ1,θ(r)φ*2,r(r)   (32)

h3(r)=φ1,r(r)φ*2,θ(r)   (33)

Thus, using Eqs. (30)-(33) the nonlinear coupling coefficient βT can be written as





  
    
      
        
          β
          T
        
        =
        
          
            
              
                
                  ∫
                  
                    
                      d
                      z
                    
                    ⁢
                    
                      ∫
                      
                        
                          d
                          θ
                        
                        ⁢
                        
                          ∫
                          
                            
                              d
                              rr
                            
                            ⁢
                            
                              
                            
                            ⁢
                            
                              ɛ
                              0
                            
                            ⁢
                            
                              
                                d
                                14
                              
                              ⁡
                              
                                (
                                r
                                )
                              
                            
                            ⁢
                            
                              
                                E
                                THz
                                *
                              
                              ⁡
                              
                                (
                                
                                  r
                                  ,
                                  θ
                                
                                )
                              
                            
                            ⁢
                            
                              exp
                              ⁡
                              
                                (
                                
                                  
                                    -
                                    i
                                  
                                  ⁢
                                  
                                    
                                  
                                  ⁢
                                  Δ
                                  ⁢
                                  
                                    
                                  
                                  ⁢
                                  m
                                  ⁢
                                  
                                    
                                  
                                  ⁢
                                  θ
                                
                                )
                              
                            
                          
                        
                      
                    
                  
                
              
            
            
              
                
                  [
                  
                    
                      
                        exp
                        ⁡
                        
                          (
                          
                            2
                            ⁢
                            i
                            ⁢
                            
                              
                            
                            ⁢
                            θ
                          
                          )
                        
                      
                      ⁢
                      
                        
                          g
                          1
                        
                        ⁡
                        
                          (
                          r
                          )
                        
                      
                    
                    +
                    
                      
                        exp
                        ⁡
                        
                          (
                          
                            
                              -
                              2
                            
                            ⁢
                            i
                            ⁢
                            
                              
                            
                            ⁢
                            θ
                          
                          )
                        
                      
                      ⁢
                      
                        
                          g
                          1
                          *
                        
                        ⁡
                        
                          (
                          r
                          )
                        
                      
                    
                    +
                    
                      
                        g
                        3
                      
                      ⁡
                      
                        (
                        r
                        )
                      
                    
                  
                  ]
                
              
            
          
          
            
              
                
                  
                    
                      ∫
                      
                        
                          d
                          
                            r
                            r
                          
                        
                        ⁢
                        
                          ɛ
                          0
                        
                        ⁢
                        
                          
                            n
                            1
                            2
                          
                          ⁡
                          
                            (
                            
                              r
                              r
                            
                            )
                          
                        
                        ⁢
                        
                          
                            
                            
                              
                                
                                  E
                                  r
                                
                                1
                              
                              ⁡
                              
                                (
                                
                                  r
                                  r
                                
                                )
                              
                            
                            
                          
                          2
                        
                      
                    
                  
                  ⁢
                  
                    
                      ∫
                      
                        
                          d
                          
                            r
                            r
                          
                        
                        ⁢
                        
                          ɛ
                          0
                        
                        ⁢
                        
                          
                            n
                            2
                            2
                          
                          ⁡
                          
                            (
                            
                              r
                              r
                            
                            )
                          
                        
                        ⁢
                        
                          
                            
                            
                              
                                
                                  E
                                  r
                                
                                2
                              
                              ⁡
                              
                                (
                                
                                  r
                                  r
                                
                                )
                              
                            
                            
                          
                          2
                        
                      
                    
                  
                
              
            
            
              
                
                  
                    ∫
                    
                      
                        d
                        
                          r
                          r
                        
                      
                      ⁢
                      
                        ɛ
                        0
                      
                      ⁢
                      
                        
                          n
                          T
                          2
                        
                        ⁡
                        
                          (
                          
                            r
                            r
                          
                          )
                        
                      
                      ⁢
                      
                        
                          
                          
                            
                              
                                E
                                r
                              
                              T
                            
                            ⁡
                            
                              (
                              
                                r
                                r
                              
                              )
                            
                          
                          
                        
                        2
                      
                    
                  
                
              
            
          
        
      
    
    
      
        (
        34
        )
      
    
  
  
    
      with
    
    
      
        
      
    
  
  
    
      
        
          
            g
            1
          
          ⁡
          
            (
            r
            )
          
        
        =
        
          
            
              h
              2
            
            ⁡
            
              (
              r
              )
            
          
          -
          
            
              i
              2
            
            ⁢
            
              
                h
                1
              
              ⁡
              
                (
                r
                )
              
            
          
        
      
    
    
      
        (
        35
        )
      
    
  




g2(r)=h3(r)−2h2(r)   (36)

Thus, from Eq. (34) it can be deduced a large value of the nonlinear coupling coefficient can be obtained when the dependence of E*THz(r,θ) on the azimuthal coordinate θ cancels the modulation introduced in the integrand of Eq. (34) by exp[i(Δm±2)]. This condition can be viewed as a generalization of the canonical phase-matching condition obtained for efficient nonlinear frequency mixing in standard dielectric waveguides. For GaAs, assuming that the refractive indexes at λ2=1.550 μm and λT=300 μm are n2=3.45 and nT=3.61, respectively, an optimal configuration is found that corresponds to two WGM's with m1=575 and m2=572 embedded in a dipole-like defect mode of a photonic crystal. Note that in this case m1−m2=3, he dipole THz-scale (characterized by m=1) phase matches the dependence exp[i(Δm±2)]. Thus, using the corresponding electric field components computed by means of a FDTD method, one can obtain a value for nonlinear coupling coefficient of β=1.2 J−1/2.
The THz output power in this system can be computed using from Eqs. (2)-(5), but now replacing the reflected field amplitudes corresponding to the pump and idler frequencies defined by Eq. (1) by the transmitted amplitudes, as shown in inset of FIG. 8 for the definition of the different magnitudes employed in the CMT description of this system. FIG. 8 displays the conversion efficiency, obtained by assuming that the Q-factors for the three resonant frequencies are the same and equal to 5×104. In addition, in this calculation it is assumed P1in=P2in. The dependence of the conversion efficiency on the ratio S22/S12 can be deduced from FIG. 7. As can be seen in FIG. 8, for P1in=P2in, the maximum conversion efficiency when no losses are included in the calculation is given by ηmax=ωT/(2ω1). As discussed for the quantum-limited efficiency, the effect of considering the absorption losses in THz will simply rescale the conversion efficiencies displayed in FIG. 8 by a factor Γ1ΓT.
It has been demonstrated that overlap constants (β) of 1.2 J−1/2 or greater using both of the above implementations for triply resonant THz generation. As mentioned in the herein, the power levels at which THz radiation can be generated most efficiently in a given implementation depends only on the three frequencies, the overlap constant and the product of the three mode Q-factors. Here one can assume that the cavities are overcoupled to their respective waveguides and thus neglect losses. THz generation with quantum limited efficiency requires that the product 16Q1Q2QT|β|2Pin/ω12 is of order 1 or greater. Total power conversion is optimized when this product is equal to 4. FIG. 9B shows the optimal input power as a function of the Q-factor product for β=1.2 J−1/2. For some applications where a single cavity must operate across a range of different power levels, the best performance is achieved by varying the two input powers (|s1+|2+|s1+|2) to fit the trajectory outlined by the dashed line 100 in FIG. 7.
This optimal output power is shown in FIG. 9A for the for Q1Q2QT ̃1×1015, a value that is attainable since for example Q-factors larger than 106 have been demonstrated for all three modes in the PhCNB cavity geometry.
The invention provides feasible platforms for high efficiency THz generation in triply resonant photonic resonators. First, by placing a double mode NIR/visible cavity with closely spaced resonant modes near the field maximum of a much larger cavity with a THz mode, the invention allows a high nonlinear mode overlap to be attained through independent design of the NIR and THz cavities, exploiting their vastly different length-scales. In addition, the invention provides an implementation where a dual-mode NIR cavity acts as a component of a THz-scale cavity, with the NIR cavity having a considerable effect on the THz modes. In both systems, complete pump depletion can occur at all powers provided the correct pump:idler power relationship. Total power conversion approaches the quantum limit at a critical power which can be tuned by varying the product of the three mode Q factors.
In addition, the invention uses a general analysis of nonlinear second order difference frequency generation (DFG) where all three modes are resonant. The equations governing the system can be scaled into a set of universal dimensionless parameters and that the steady state solutions (corresponding to CW conversion) only depend on the generalized parameters, (labeled as P1 and P2), which are directly proportional to the input powers of the pump and idler respectively. As P1 increases for any fixed P2, a saddle-node bifurcation is reached, below which there only one steady-state solution whose stability is universal, and above which, there exist three solutions whose stabilities can vary with geometry and in many cases bistability.
The CW conversion from pump to signal with quantum-limited efficiency (the pump beam is perfectly depleted and the pump:signal conversion ratio is given by the frequency ratio) is possible in such a system for any power of the pump wave, provided a specific input power at the idler frequency. Quantum-limited conversion can occur with arbitrarily small idler power when the pump power approaches the critical value |s1+|2=ω1/(4|β|2Q1Q2QTΓ1) (P1=4); where β quantifies the mode overlap, which implicitly depends on the cavity mode-volumes and the phase matching of the modes, Γ1 denotes losses at the pump frequency, and {Qk} are the cavity mode quality factors. This critical power also acts as a threshold, below-which quantum-limited conversion is always stable and above-which the stability depends on geometry.
Linear losses do not qualitatively change the dynamics, but do cause a reduction in conversion efficiency. The efficiency is maximized when theoretical quality factors are kept sufficiently low that they are limited by coupling to the input/output port. A similar universal framework is expected to exist for second-order sum-frequency generation as well.
Although the present invention has been shown and described with respect to several preferred embodiments thereof, various changes, omissions and additions to the form and detail thereof, may be made therein, without departing from the spirit and scope of the invention.

## Claims

1. A system for efficient generation of THz radiation comprising:
a triply-resonant nonlinear photonic resonator coupled to at least one near-infrared (NIR) or optical waveguide and to at least one THz waveguide, wherein energy traveling through the at least one near-infrared (NIR) or optical waveguide is converted to THz radiation inside the triply-resonant photonic resonator via a nonlinear difference frequency generation (DFG) process.

2. The system of claim 1, wherein the triply-resonant photonic resonator comprises a nonlinear photonic resonator for NIR or optical frequencies embedded in or placed in the proximity of a photonic resonator for THz frequencies, in such a way that if the NIR cavity is removed a THz cavity mode would be significantly modified.

3. The system of claim 2, wherein the nonlinear photonic resonator for NIR or optical frequencies comprises a doubly-resonant photonic resonator, the doubly-resonant photonic resonator comprises two resonant wavelengths λ1 and λ2 in air that support two resonant modes characterized by quality factors Q1 and Q2, respectively, and by modal volumes V1 and V2, respectively; with 400 nm<λ1<20 μm and 400 nm<λ2<20 μm.

4. The system of claim 2, wherein the photonic resonator for THz frequencies comprises a resonant wavelength in air λT, a quality factor Q, and a modal volume VT; with 30 μm<λT<3000 μm.

5. The system of claim 4, wherein Q1>100, V1<1000 λ13, Q2>100, V2<1000 λ23, QT>100, and, VT<1000 λ13.

6. The system of claim 3, wherein the doubly-resonant photonic resonator comprises a photonic crystal cavity or a longitudinally perturbed ridge waveguide or undercut ridge waveguide.

7. The system of claim 2, wherein the photonic resonator comprises a longitudinally perturbed ridge waveguide or undercut ridge waveguide or photonic crystal cavity.

8. The system of claim 3, wherein the resonant modes comprise whispering gallery modes.

9. The system of claim 1, wherein the at least one near-infrared (NIR) or optical waveguide comprises a dielectric waveguide or a photonic crystal waveguide.

10. The system of claim 1, wherein the at least one THz waveguide comprises a dielectric waveguide or a photonic crystal waveguide.

11. The system of claim 3, wherein the nonlinear photonic resonator for NIR or optical frequencies comprises a III/V semiconductor.

12. A method for generating efficient THz radiation comprising:
providing a triply-resonant nonlinear photonic resonator coupled to at least one near-infrared (NIR) or optical waveguide and to at least one THz waveguide; and
converting the THz radiation inside the triply-resonant photonic resonator via a nonlinear difference frequency generation (DFG) process.

13. The method of claim 12, wherein the triply-resonant photonic resonator comprises a nonlinear photonic resonator for NIR or optical frequencies embedded in or placed in the proximity of a photonic resonator for THz frequencies, in such a way that if the NIR cavity is removed a THz cavity mode would be significantly modified.

14. The method of claim 12, wherein the nonlinear photonic resonator for NIR or optical frequencies comprises a doubly-resonant photonic resonator, the doubly-resonant photonic resonator comprises two resonant wavelengths λ1 and λ2 in air that support two resonant modes characterized by quality factors Q1 and Q2, respectively, and by modal volumes V1 and V2, respectively; with 400 nm<λ1<20 μm and 400 nm<λ2<20 μm.

15. The method of claim 13, wherein the photonic resonator for THz frequencies comprises a resonant wavelength in air λT, a quality factor Q, and a modal volume VT; with 30 μm<λT<3000 μm.

16. The method of claim 15, wherein Q1>100 , V1<1000 λ13, Q2>100, V2<1000 λ23, QT>100, and, VT<1000 λ13.

17. The method of claim 14, wherein the doubly-resonant photonic resonator comprises a photonic crystal cavity or a longitudinally perturbed ridge waveguide or undercut ridge waveguide.

18. The method of claim 13, wherein the photonic resonator comprises a longitudinally perturbed ridge waveguide or undercut ridge waveguide or photonic crystal cavity.

19. The method of claim 14, wherein the resonant modes comprise whispering gallery modes.

20. The method of claim 12, wherein the at least one near-infrared (NIR) or optical waveguide comprises a dielectric waveguide or a photonic crystal waveguide.

21. The method of claim 12, wherein the at least one THz waveguide comprises a dielectric waveguide or a photonic crystal waveguide.

22. The system of claim 14, wherein the nonlinear photonic resonator for NIR or optical frequencies comprises a III/V semiconductor.