The design and modeling of microstructured optical fiber Steven G. Johnson MIT / Harvard University Outline What are these fibers (and why should I care)? The guiding mechanisms: index-guiding and band gaps Finding the guided modes Small corrections (with big impacts) Outline What are these fibers (and why should I care)? The guiding mechanisms: index-guiding and band gaps Finding the guided modes Small corrections (with big impacts) Optical Fibers Today (not to scale) losses ~ 0.2 dB/km more complex profiles to tune dispersion

high index doped-silica core n ~ 1.46 silica cladding n ~ 1.45 at =1.55m (amplifiers every 50100km) LP01 confined mode field diameter ~ 8m protective polymer sheath [ R. Ramaswami & K. N. Sivarajan, Optical Networks: A Practical Perspective ] but this is ~ as good as it gets

The Glass Ceiling: Limits of Silica Loss: amplifiers every 50100km limited by Rayleigh scattering (molecular entropy) cannot use exotic wavelengths like 10.6m Nonlinearities: after ~100km, cause dispersion, crosstalk, power limits (limited by mode area ~ single-mode, bending loss) also cannot be made (very) large for compact nonlinear devices Radical modifications to dispersion, polarization effects? tunability is limited by low index contrast Long Distances High Bit-Rates Compact Devices Dense Wavelength Multiplexing (DWDM) Breaking the Glass Ceiling: Hollow-core Bandgap Fibers Bragg fiber

1000x better loss/nonlinear limits [ Yeh et al., 1978 ] (from density) 1d crystal + omnidirectional = OmniGuides 2d crystal Photonic Crystal (You can also put stuff in here ) PCF [ Knight et al., 1998 ] Breaking the Glass Ceiling:

Hollow-core Bandgap Fibers Bragg fiber [ Yeh et al., 1978 ] [ figs courtesy Y. Fink et al., MIT ] white/grey = chalco/polymer + omnidirectional = OmniGuides silica [ R. F. Cregan et al., Science 285, 1537 (1999) ] 5m PCF [ Knight et al., 1998 ]

Breaking the Glass Ceiling: Hollow-core Bandgap Fibers Guiding @ 10.6m [ figs courtesy Y. Fink et al., MIT ] (high-power CO2 lasers) loss < 1 dB/m white/grey = chalco/polymer (material loss ~ 104 dB/m) [ Temelkuran et al., Nature 420, 650 (2002) ] silica [ R. F. Cregan et al., Science 285, 1537 (1999) ]

5m Guiding @ 1.55m loss ~ 13dB/km [ Smith, et al., Nature 424, 657 (2003) ] OFC 2004: 1.7dB/km BlazePhotonics Breaking the Glass Ceiling II: Solid-core Holey Fibers solid core holey cladding forms effective low-index material Can have much higher contrast than doped silica strong confinement = enhanced nonlinearities, birefringence,

[ J. C. Knight et al., Opt. Lett. 21, 1547 (1996) ] Breaking the Glass Ceiling II: Solid-core Holey Fibers endlessly single-mode [ T. A. Birks et al., Opt. Lett. 22, 961 (1997) ] polarization -maintaining [ K. Suzuki, Opt. Express 9, 676 (2001) ] nonlinear fibers [ Wadsworth et al., JOSA B 19, 2148 (2002) ] low-contrast linear fiber

(large area) [ J. C. Knight et al., Elec. Lett. 34, 1347 (1998) ] Outline What are these fibers (and why should I care)? The guiding mechanisms: index-guiding and band gaps Finding the guided modes Small corrections (with big impacts) Universal Truths: Conservation Laws an arbitrary-shaped fiber (1) Linear, time-invariant system: (nonlinearities are small correction) z frequency is conserved cladding (2) z-invariant system:

(bends etc. are small correction) wavenumber is conserved core electric (E) and magnetic (H) fields can be chosen: E(x,y) ei(z t), H(x,y) ei(z t) Sequence of Computation 1 Plot all solutions of infinite cladding as vs. light cone empty spaces (gaps): guiding possibilities 2

Core introduces new states in empty spaces plot () dispersion relation 3 Compute other stuff Conventional Fiber: Uniform Cladding c 2 2 = +kt n c n uniform cladding, index n kt (transverse wavevector)

light cone light line: =c/n Conventional Fiber: Uniform Cladding c 2 2 = +kt n c n uniform cladding, index n light cone higher-order

core with higher index n pulls down index-guided mode(s) fundamental = c / n' PCF: Periodic Cladding periodic cladding (x,y) a Blochs Theorem for periodic systems: fields can be written: E(x,y) ei(z+kt xt t), periodic functions on primitive cell H(x,y) ei(z+kt xt t) transverse (xy) Bloch wavevector kt

primitive cell 1 2 kt , kt , H = 2 H satisfies c eigenproblem (Hermitian constraint: k , H =0 t if lossless) where: kt , = +ikt +i z PCF: Cladding Eigensolution Finite cell discrete eigenvalues n Want to solve for n(kt, ),

& plot vs. for all n, kt 1 n2 kt , kt , Hn = 2 Hn c constraint: kt , H =0 z where: kt , = +ikt +i H(x,y) ei(z+kt xt t) 1 Limit range of kt: irreducible Brillouin zone 2

Limit degrees of freedom: expand H in finite basis 3 Efficiently solve eigenproblem: iterative methods PCF: Cladding Eigensolution 1 Limit range of kt: irreducible Brillouin zone Blochs theorem: solutions are periodic in kt K first Brillouin zone = minimum |kt| primitive cell M 4 a 3 ky

kx irreducible Brillouin zone: reduced by symmetry 2 Limit degrees of freedom: expand H in finite basis 3 Efficiently solve eigenproblem: iterative methods PCF: Cladding Eigensolution 1 Limit range of kt: irreducible Brillouin zone 2 Limit degrees of freedom: expand H in finite basis must satisfy constraint: kt , H =0 Planewave (FFT) basis Finite-element basis

constraint, boundary conditions: H(xt ) = HGe iGxt Ndlec elements G constraint: [ Ndlec, Numerische Math. 35, 315 (1980) ] HG (G +k+ z) =0 uniform grid, periodic boundaries, simple code, O(N log N) 3 [ figure: Peyrilloux et al.,

J. Lightwave Tech. 21, 536 (2003) ] nonuniform mesh, more arbitrary boundaries, complex code & mesh, O(N) Efficiently solve eigenproblem: iterative methods PCF: Cladding Eigensolution 1 Limit range of kt: irreducible Brillouin zone 2 Limit degrees of freedom: expand H in finite basis (N) N H =H(xt ) = hmbm(xt ) 2 solve: A H = H

m=1 finite matrix problem: f g =f* g 3 2 Ah= Bh Aml = bm A bl Bml = bm bl Efficiently solve eigenproblem: iterative methods PCF: Cladding Eigensolution 1 Limit range of kt: irreducible Brillouin zone 2

Limit degrees of freedom: expand H in finite basis 3 Efficiently solve eigenproblem: iterative methods 2 Ah= Bh Slow way: compute A & B, ask LAPACK for eigenvalues requires O(N2) storage, O(N3) time Faster way: start with initial guess eigenvector h0 iteratively improve O(Np) storage, ~ O(Np2) time for p eigenvectors (p smallest eigenvalues) PCF: Cladding Eigensolution 1 Limit range of kt: irreducible Brillouin zone 2

Limit degrees of freedom: expand H in finite basis 3 Efficiently solve eigenproblem: iterative methods 2 Ah= Bh Many iterative methods: Arnoldi, Lanczos, Davidson, Jacobi-Davidson, , Rayleigh-quotient minimization PCF: Cladding Eigensolution 1 Limit range of kt: irreducible Brillouin zone 2 Limit degrees of freedom: expand H in finite basis 3 Efficiently solve eigenproblem: iterative methods

2 Ah= Bh Many iterative methods: Arnoldi, Lanczos, Davidson, Jacobi-Davidson, , Rayleigh-quotient minimization for Hermitian matrices, smallest eigenvalue 0 minimizes: variational theorem h' Ah =min h h' Bh 2 0 minimize by conjugate-gradient, (or multigrid, etc.) PCF: Holey Silica Cladding r = 0.1a 2r

n=1.46 a (2c/a) light cone = dimensionless units: Maxwells equations are scale-invariant c (2/a) PCF: Holey Silica Cladding r = 0.17717a (2c/a)

light cone = c (2/a) 2r n=1.46 a PCF: Holey Silica Cladding r = 0.22973a (2c/a) light cone

= c (2/a) 2r n=1.46 a PCF: Holey Silica Cladding r = 0.30912a (2c/a) light cone =

c (2/a) 2r n=1.46 a PCF: Holey Silica Cladding r = 0.34197a (2c/a) light cone = c

(2/a) 2r n=1.46 a PCF: Holey Silica Cladding r = 0.37193a (2c/a) light cone = c (2/a) 2r

n=1.46 a PCF: Holey Silica Cladding r = 0.4a (2c/a) light cone = c (2/a) 2r n=1.46

a PCF: Holey Silica Cladding r = 0.42557a (2c/a) light cone = c (2/a) 2r n=1.46 a PCF: Holey Silica Cladding

n=1.46 2r a r = 0.45a light cone (2c/a) gap-guided modes go here = c index-guided modes go here

(2/a) PCF: Holey Silica Cladding r = 0.45a light cone 2r n=1.46 a above air line: (2c/a) guiding in air core is possible air ht g li

e lin = c (2/a) below air line: surface states of air core Bragg Fiber Cladding Graphics are QuickTime needed decompressor toand seeathis picture. at large radius, becomes ~ planar

nhi = 4.6 Bragg fiber gaps (1d eigenproblem) nlo = 1.6 radial kr (Bloch wavevector) k 0 by conservation of angular momentum wavenumber = 0: normal incidence Omnidirectional Cladding Graphics are QuickTime

needed decompressor toand seeathis picture. omnidirectional (planar) reflection e.g. light from fluorescent sources is trapped Bragg fiber gaps (1d eigenproblem) for nhi / nlo big enough and nlo > 1 [ J. N. Winn et al, Opt. Lett. 23, 1573 (1998) ] wavenumber = 0: normal incidence

Outline What are these fibers (and why should I care)? The guiding mechanisms: index-guiding and band gaps Finding the guided modes Small corrections (with big impacts) Sequence of Computation 1 Plot all solutions of infinite cladding as vs. light cone empty spaces (gaps): guiding possibilities 2 Core introduces new states in empty spaces plot () dispersion relation 3

Compute other stuff Computing Guided (Core) Modes 1 n2 Hn = 2 Hn c constraint: H =0 z where: = +i magnetic field = Same differential equation as before, except no kt can solve the same way H(x,y) ei(z t) New considerations:

1 Boundary conditions 2 Leakage (finite-size) radiation loss 3 Interior eigenvalues Computing Guided (Core) Modes 1 computational cell Boundary conditions Only care about guided modes: exponentially decaying outside core Effect of boundary cond. decays exponentially mostly, boundaries are irrelevant! periodic (planewave), conducting, absorbing all okay

2 Leakage (finite-size) radiation loss 3 Interior eigenvalues Guided Mode in a Solid Core Graphics are QuickTime needed decompressor toand seeathis picture. small computation: only lowest- band! (~ one minute, planewave) 1.46 c/ = 1.46 neff 0.12

0.1 holey PCF light cone flux density 0.08 0.06 fundamental mode 0.04 (two polarizations) 2r 0.02 0 0.3 endlessly single mode: neff decreases with 0.4

0.5 0.6 0.7 0.8 /a 0.9 1 1.1 1.2 n=1.46 a r = 0.3a

Fixed-frequency Modes? Here, we are computing ('), but we often want (') is specified No problem! Just find root of (') ', using Newtons method: d d (Factor of 34 in time.) group velocity = power / (energy density) (a.k.a. Hellman-Feynman theorem, a.k.a. first-order perturbation theory, a.k.a. k-dot-p theory) Computing Guided (Core) Modes 1 computational cell Boundary conditions

Only care about guided modes: exponentially decaying outside core Effect of boundary cond. decays exponentially mostly, boundaries are irrelevant! periodic (planewave), conducting, absorbing all okay except when we want (small) finite-size losses 2 Leakage (finite-size) radiation loss 3 Interior eigenvalues Computing Guided (Core) Modes 1 Boundary conditions 2 Leakage (finite-size) radiation loss

Use PML absorbing boundary layer perfectly matched layer [ Berenger, J. Comp. Phys. 114, 185 (1994) ] with iterative method that works for non-Hermitian (dissipative) systems: Jacobi-Davidson, Or imaginary-distance BPM: [ Saitoh, IEEE J. Quantum Elec. 38, 927 (2002) ] in imaginary z, largest (fundamental) mode grows exponentially 3 Interior eigenvalues Computing Guided (Core) Modes d n=1.45 1 Boundary conditions 2

Leakage (finite-size) radiation loss imaginary-distance BPM [ Saitoh, IEEE J. Quantum Elec. 38, 927 (2002) ] 2 rings 3 3 rings Interior eigenvalues Computing Guided (Core) Modes 1 Boundary conditions 2 Leakage (finite-size) radiation loss 3

Interior eigenvalues [ J. Broeng et al., Opt. Lett. 25, 96 (2000) ] 2.4 fundamental & 2nd order (~ N states for N-hole cell) but most methods compute smallest (or largest ) guided modes 2.0 (2c/a) Gap-guided modes lie above continuum air 1.6

l t li h ig ne fundamental air-guided mode 1.2 bulk crystal continuum 0.8 1.11 1.27 1.43 1.59

1.75 1.91 (2/a) 2.07 2.23 2.39 Computing Guided (Core) Modes 1 Boundary conditions 2 Leakage (finite-size) radiation loss 3

Interior (of the spectrum) eigenvalues i Gap-guided modes lie above continuum (~ N states for N-hole cell) but most methods compute smallest (or largest ) Compute N lowest states first: deflation (orthogonalize to get higher states) [ see previous slide ] ii Use interior eigensolver method closest eigenvalues to (mid-gap) Jacobi-Davidson, Arnoldi with shift-and-invert, smallest eigenvalues of (A2)2 convergence often slower

iii Other methods: FDTD, etc Interior Eigenvalues by FDTD finite-difference time-domain Simulate Maxwells equations on a discrete grid, + PML boundaries + eiz z-dependence Excite with broad-spectrum dipole ( ) source signal processing complex n [ Mandelshtam, J. Chem. Phys. 107, 6756 (1997) ] Response is many sharp peaks, one peak per mode decay rate in time gives loss: Im[] = Im[] / d/d Interior Eigenvalues by FDTD

finite-difference time-domain Simulate Maxwells equations on a discrete grid, + PML boundaries + eiz z-dependence Excite with broad-spectrum dipole ( ) source Response is many sharp peaks, one peak per mode mode field profile narrow-spectrum source An Easier Problem: Bragg-fiber Modes In each concentric region, solutions are Bessel functions: c Jm (kr) + d Ym(kr) im e 2

k = 2 c angular momentum At circular interfaces match boundary conditions with 4 4 transfer matrix search for complex that satisfies: finite at r=0, outgoing at r= [ Johnson, Opt. Express 9, 748 (2001) ] Hollow Metal Waveguides, Reborn Graphics are QuickTime needed decompressor toand seeathis picture. Graphics are

QuickTime needed decompressor toand seeathis picture. metal waveguide modes OmniGuide fiber modes frequency Graphics are QuickTime needed decompressor toand seeathis picture. 1970s microwave tubes @ Bell Labs wavenumber

wavenumber modes are directly analogous to those in hollow metal waveguide An Old Friend: the TE01 mode lowest-loss mode, just as in metal (near) node at interface = strong confinement = low losses non-degenerate mode cannot be split = no birefringence or PMD r E Bushels of Bessels A General Multipole Method [ White, Opt. Express 9, 721 (2001) ]

Each cylinder has its own Bessel expansion: M only cylinders allowed field~ cmJ m +dmYm m (m is not conserved) With N cylinders, get 2NM 2NM matrix of boundary conditions Solution gives full complex , but takes O(N3) time more than 45 periods is difficult future: Fast Multipole Method should reduce to O(N log N)? Outline What are these fibers (and why should I care)? The guiding mechanisms: index-guiding and band gaps Finding the guided modes

Small corrections (with big impacts) All Imperfections are Small (or the fiber wouldnt work) Material absorption: small imaginary Nonlinearity: small ~ |E|2 Acircularity (birefringence): small boundary shift Bends: small ~ x / Rbend Roughness: small or boundary shift Weak effects, long distances: hard to compute directly use perturbation theory Perturbation Theory and Related Methods (Coupled-Mode Theory, Volume-Current Method, etc.) Given solution for ideal system compute approximate effect of small changes solves hard problems starting with easy problems & provides (semi) analytical insight Perturbation Theory

for Hermitian eigenproblems u =uu O u & u for small O given eigenvectors/values: find change Solution: expand as power series in O (1) (2) u =0+u +u + u u O u = uu

(1) & u =0+ u + (1) (first order is usually enough) Perturbation Theory for electromagnetism c H A H (1) = 2 H H 2 E =

2 E 2 2 (1) (1) = / vg e.g. absorption gives imaginary = decay! d vg = d A Quantitative Example but what about the cladding?

Gas can have low loss & nonlinearity some field penetrates! & may need to use very bad material to get high index contrast Suppressing Cladding Losses Graphics are QuickTime needed decompressor toand seeathis picture. Material absorption: small imaginary 1x10-2 Mode Losses

Bulk Cladding Losses EH11 1x10-3 Large differential loss TE01 strongly suppresses cladding absorption (like ohmic loss, for metal) 1x10-4 TE01 1x10-51.2 1.6 2 (m) 2.4 2.8

High-Power Transmission at 10.6m (no previous dielectric waveguide) Polymer losses @10.6m ~ 50,000dB/m -3.0 waveguide losses ~ 1dB/m 8 -3.5 6 [ B. Temelkuran et al., Nature 420, 650 (2002) ] -4.0 slope = -0.95 dB/m R2 = 0.99 4 -4.5 2.5

3.0 3.5 Length (meters) 2 4.0 0 5 5 6 6 7 7 8 8

9 99 10 10 10 11 11 11 12 12 12 Wavelength (m) [ figs courtesy Y. Fink et al., MIT ] Quantifying Nonlinearity Kerr nonlinearity: small ~ |E|2 ~ power P ~ 1 / lengthscale for nonlinear effects = / P = nonlinear-strength parameter determining self-phase modulation (SPM), four-wave mixing (FWM), (unlike effective area, tells where the field is, not just how big) Suppressing Cladding Nonlinearity

Graphics are QuickTime needed decompressor toand seeathis picture. Mode Nonlinearity* Cladding Nonlinearity 1x10-6 1x10-7 Will be dominated by nonlinearity of air TE01 1x10-8 ~10,000 times weaker than in silica fiber

(including factor of 10 in area) 1x10-91.2 * nonlinearity = (1) / P = 1.6 2 (m) 2.4 2.8 Acircularity & Perturbation Theory (or any shifting-boundary problem) = 1 2 2 1

= 2 1 just plug s into perturbation formulas? FAILS for high index contrast! beware field discontinuity fortunately, a simple correction exists [ S. G. Johnson et al., PRE 65, 066611 (2002) ] Acircularity & Perturbation Theory (or any shifting-boundary problem) = 1 2 2 1 = 2 1 (continuous field components) h

1 2 2 h E|| D (1) = surf. 2 2 E [ S. G. Johnson et al., PRE 65, 066611 (2002) ] Loss from Roughness/Disorder imperfection acts like a volume current r r

J ~ E0 volume-current method or Greens functions with first Born approximation Loss from Roughness/Disorder imperfection acts like a volume current r r J ~ E0 For surface roughness, including field discontinuities: r r r 1 J ~ E|| D Loss from Roughness/Disorder

uncorrelated disorder adds incoherently So, compute power P radiated by one localized source J, and loss rate ~ P * (mean disorder strength) Effect of an Incomplete Gap on uncorrelated surface roughness loss some radiation blocked reflection same reflection radiation Conventional waveguide (matching modal area) with Si/SiO2 Bragg mirrors (1D gap) 50% lower losses (in dB) same reflection

Considerations for Roughness Loss Band gap can suppress some radiation typically by at most ~ 1/2, depending on crystal Loss ~ 2 ~ 1000 times larger than for silica Loss ~ fraction of |E|2 in solid material factor of ~ 1/5 for 7-hole PCF ~ 10-5 for large-core Bragg-fiber design Hardest part is to get reliable statistics for disorder. Using perturbations to design big effects Perturbation Theory and Dispersion when two distinct modes cross & interact, unusual dispersion is produced mode 1 mode 2 no interaction/coupling

Perturbation Theory and Dispersion when two distinct modes cross & interact, unusual dispersion is produced mode 1 mode 2 coupling: anti-crossing Two Localized Modes = Very Strong Dispersion core mode localized claddingdefect mode weak coupling = rapid slope change = high dispersion

(> 500,000 ps/nm-km + dispersion-slope matching) [ T. Engeness et al., Opt. Express 11, 1175 (2003) ] Slow-light Modes = Anomalous Dispersion (Different-Symmetry) slow-light band edges at =0 TM TE =0 point has additional symmetry:

modes can be purely TE/TM polarized force different symmetry modes together [ M. Ibanescu et al., Phys. Rev. Lett. 92, 063903 (2004) ] Slow-light Modes = Anomalous Dispersion (Different-Symmetry) ultra-flat (4) backward wave slow light non-zero velocity [ M. Ibanescu et al., Phys. Rev. Lett. 92, 063903 (2004) ] Slow-light Modes = Anomalous Dispersion (Different-Symmetry)

Uses gap at =0: perfect metal [1960] or Bragg fiber or high-index PCF (n > 2.5) [ M. Ibanescu et al., Phys. Rev. Lett. 92, 063903 (2004) ] Further Reading Reviews: J. D. Joannopoulos, R. D. Meade, and J. N. Winn, Photonic Crystals: Molding the Flow of Light (Princeton Univ. Press, 1995). P. Russell, Photonic-crystal fibers, Science 299, 358 (2003). This Presentation, Free Software, Other Material: http://ab-initio.mit.edu/photons/tutorial