compute-methods_cornell by liaoxiuli

VIEWS: 17 PAGES: 69

									  Computational Photonics:
Frequency and Time Domain Methods


         Steven G. Johnson
      MIT Applied Mathematics
          Nano-photonic media (l-scale)
                     strange waveguides


                                              & microcavities
[B. Norris, UMN]    [Assefa & Kolodziejski,
                             MIT]
                                                                3d
                                                                structures

                                [Mangan,
                                Corning]




   synthetic materials
                                                            optical phenomena
                                  hollow-core fibers
                 Photonic Crystals
                   periodic electromagnetic media
1887                     1987
          1-D                    2-D             3-D




        periodic in            periodic in       periodic in
       one direction         two directions   three directions




   can have a band gap: optical ―insulators‖
                Electronic and Photonic Crystals
                atoms in diamond structure      dielectric spheres, diamond lattice
Medium
Periodic




                                                  photon frequency
Band Diagram
 Bloch waves:


                 electron energy




                                   wavevector                        wavevector
   Electronic & Photonic Modeling
       Electronic                  Photonic
• strongly interacting     • non-interacting (or weakly),
  —entanglement, Coulomb     —simple approximations
  —tricky approximations       (finite resolution)
                             —any desired accuracy

• lengthscale dependent    • scale-invariant
    (from Planck’s h)        —e.g. size 10  l 10
                               (except materials may change)
Computational Photonics Problems
• Time-domain simulation
      — start with current J(x,t)
      — run ―numerical experiment‖ to simulate E(x, t), H(x, t)

• Frequency-domain linear response
      — start with harmonic current J(x, t) = e–it J(x)
      — solve for steady-state harmonic fields E(x), H(x)
      — involves solving linear equation Ax=b

• Frequency-domain eigensolver
      — solve for source-free harmonic eigenfields
             E(x), H(x) ~ e–it
      — involves solving eigenequation Ax=2x
 Numerical Methods: Basis Choices
  finite difference                                 finite elements
                                                                      in irregular ―elements,‖
        discretize                                                    approximate unknowns
        unknowns                                                      by low-degree polynomial
      on regular grid

df   f ( x  x)  f ( x  x)
                               O(x 2 )
dx              x

                                                            boundary-element methods
spectral methods                                                       discretize only the
                                                                       boundaries between
                                   complete basis of                   homogeneous media
 +                                 smooth functions
                                   (e.g. Fourier series)
                                                                               …solve
        +                                                                      integral equation
               + …..                                                   via Green’s functions
   Numerical Methods: Basis Choices
   finite difference                        finite elements
                                                              in irregular ―elements,‖
         discretize                                           approximate unknowns
         unknowns                                             by low-degree polynomial
       on regular grid

                                                    boundary-element methods
  spectral methods
                                                               discretize only the
                                                               boundaries between
                         complete basis of                     homogeneous media
   +                     smooth functions
                         (e.g. Fourier series)
       +                                                               …solve
                                                                       integral equation
            + …..
                                                               via Green’s functions

Much easier to analyze, implement,                Potentially much more efficient,
generalize, parallelize, optimize, …               especially for high resolution
Computational Photonics Problems                          Numerical Methods: Basis Choices

                                                            finite difference
• Time-domain simulation
    — start with current J(x,t)
    — run ―numerical experiment‖ to                                                                      spectral methods
          simulate E(x, t), H(x, t)
                                                         df   f ( x  x)  f ( x  x)
                                                         dx
                                                            
                                                                         x
                                                                                         O(x 2 )      +
• Frequency-domain linear response
 — start with harmonic current J(x, t) = e–it J(x)                                                         +
                                                    
 — solve for steady-state harmonic fields E(x), H(x)                                                            + …..
 — involves solving linear equation Ax=b
                                                                 finite elements
• Frequency-domain eigensolver                                                                 in irregular ―elements,‖
                                                                                               approximate unknowns
 — solve for source-free harmonic eigenfields
                                                                                               by low-degree polynomial
             E(x), H(x) ~ e–it
 — involves solving eigenequation Ax=2x

                                                                                                     boundary-element methods
                                                                                                            discretize only the
                                                                                                            boundaries between
                                                                                                            homogeneous media
              FDTD
Finite-Difference Time-Domain methods

         Divide both space and time into discrete grids
                — spatial resolution ∆x
                — temporal resolution ∆t

         Very general: arbitrary geometries, materials,
         nonlinearities, dispersion, sources, …
           — any photonics calculation, in principle

H     1             E 1        J
       E            H 
 t                 t         

                 dielectric function (x) = n2(x)
                 The Yee Discretization (1966)
                a cubic ―voxel‖: ∆x  ∆y  ∆z                      (i, j+1)
                                                 (i+1, j+1, k+1)

                             Hz
                                                                    Ey
(i, j, k+1)                                                                     Hz

                                            Hx
     Ez              Hy                                                (i, j)        (i+1, j)
                                                 (i+1, j+1, k)
                                                                                Ex
                                            Ey
    (i, j, k)                 (i+1, j, k)
                      Ex

 Staggered grid in space:
        — every field component is stored on a different grid
          The Yee Discretization (1966)
                                                   (i, j+1)
       H     1
              E 
        t                                         Ey
                                                                Hz
H z                   1 Ey Ex 
                               
 t      1
       i , j 
                1        x   y                  (i, j)           (i+1, j)
         2      2                                               Ex
                         1              1         1                  1     
              E (i  1, j  )  Ey (i, j  ) Ex (i  , j  1)  Ex (i  , j) 
          1  y           2              2        2                  2
                                                                         
                         x                             y                
                                                                           
                                                + O(∆x2) + O(∆y2)
               all derivatives become center differences…
      The Yee Discretization (1966)
         all derivatives become center differences…
               including derivatives in time
                                               1         1
                                          H(n  )  H(n  )
H                  1                          2         2
                       E            
 t   t nt                  t nt
                                                 t
                                                            + O(∆t2)
Explicit time-stepping:
                                                  x
         stability requires t 
                                              # dimensions

      (vs. implicit time steps: invert large matrix at each step)
      FDTD Discretization Upshot

• For stability, space and time resolutions are proportional
       — doubling resolution in 3d requires
              at least 16 = 24 times the work!

• But at least the error goes quadratically with resolution
                     …right?
                     …not necessarily!
          Difficulty with a grid:
representing discontinuous materials?




                             ―staircasing‖

     … how does this affect accuracy?
Field Discontinuity Degrades
     Order of Accuracy
E
    TE polarization (E in plane: discontinuous)
                                                  a

                                                               QuickTime™ and a
                                                      TIFF (Uncompressed) decompressor
                                                         are neede d to see this picture.




             a
Sub-pixel smoothing
                 Can eliminate
                 discontinuity
               by ―grayscaling‖
            — assign some average 
                 to each pixel


       ?
   = discretizing a smoothed structure
   — that means we are changing geometry
   — can actually add to error
Past sub-pixel smoothing methods
      can make error worse!
                                                                 & convergence is
Three previous smoothing methods
                                                                 still only linear
                                    QuickTime™ and a
                           TIFF (Uncompressed) decompressor
                              are neede d to see this picture.




                                                                     [ Dey, 1999 ]
                                                                   [ Kaneda, 1997 ]
                                                                 [ Mohammadi, 2005]
 A Criterion for Accurate Smoothing
1st-order errors                       1    2 
                        ~   E||  ( ) D 
                                    2
     from
 smoothing                                 

We want the smoothing errors to be zero to 1st order
      — minimizes error and 2nd-order convergent!
              Use a tensor : 
                                                           
                                                           E||
                   (in principal axes:)
                                                         
                                                       1 
    [ Meade et al., 1993 ]                       1       E
Consistently the Lowest Error
                                                      a

quadratic accuracy!                      quadratic!                QuickTime™ and a
                                                          TIFF (Uncompressed) decompressor
                                                             are neede d to see this picture.




                       a




        [ Farjadpour et al., Opt. Lett. 2006 ]
A qualitatively different case: corners
    still ~lowest error, but not quadratic
                                                              QuickTime™ and a
                                                     TIFF (Uncompressed) decompressor
                                                        are neede d to see this picture.




                                              zero-perturbation
                                                   criterion
                                                 not satisfied
                                             due to E divergence
                                                   at corner
                                               — analytically,
                                                error ~ ∆x1.404
Yes, but what can you do with FDTD?
 Some common tasks:

   • Frequency-domain response:
          — put in harmonic source and wait for steady-state

   • Transmission/reflection spectra:
          — get entire spectrum from a single simulation
                 (Fourier transform of impulse response)

   • Eigenmodes and resonant modes:
          — get all modes from a single simulation
                 (some tricky signal processing)
     Transmission Spectra in FDTD
                          =1

a                     = 12


example: a 90° bend, 2d strip waveguide




    transmitted power = energy flux here:
       Transmission Spectra in FDTD
                           =1

   a                    = 12



Gaussian-pulse
current source J
                      Fourier-transform the fields at each x:
                       E         E(t)e
                                        it
                                               dt   E(nt)e   int
                                                                         t
                                                   n
          1
   P() 
          2
                    Re[ E   *
                                   H  ]dx
      Transmission Spectra in FDTD
   must always do two simulations: one for normalization

electric                                                                          P0()
field Ez:

                   QuickTime™ and a                    QuickTime™ and a
             YUV420 codec decompressor           YUV420 codec decompressor
            are needed to see this picture.     are needed to see this picture.




                                P()
                                              transmission = P() / P0()
        Reflection Spectra in FDTD
                            =1

                         = 12


         1
PR () 
         2
                          0               0
                 Re[(E  E )*  (H   H  )]dx


    for reflection, subtract incident fields
           (from normalization run)

         1
  P() 
         2
                   Re[ E   *
                                  H  ]dx
    Transmission/Reflection Spectra
    QuickTime™ see this
    are needed toand a picture.
    Graphics decompressor

      1

    0.9
                                   =1
    0.8

a   0.7                      = 12
    0.6

    0.5      T                        1–T–R
    0.4

    0.3

    0.2

    0.1          R
      0
       0.1   0.11    0.12   0.13    0.14   0.15   0.16   0.17   0.18   0.19   0.2

                                   a / 2πc = a / l
         Dimensionless Units

Maxwell’s equations are scale invariant
      — most useful quantities are dimensionless ratios
           like a / l, for a characteristic lengthscale a
      — same ratio, same ,  = same solution
           regardless of whether a = 1µm or 1km

Our (typical) approach:
  pick characteristic lengthscale a
      – measure distance in units of a
      – measure time in units of a/c
      – measure  in units of 2πc/a = a / l
      – ....
           Absorbing Boundaries:
     Perfectly Matched Layers
―perfect‖ absorber: PML   Artificial absorbing material
                          overlapping the computation

                          Theoretically reflectionless

                           … but PML is no longer perfect
                                with finite resolution,
                          so ―gradually turn on‖ absorption
                              over finite-thickness PML
Computational Photonics Problems                          Numerical Methods: Basis Choices

                                                            finite difference
• Time-domain simulation
    — start with current J(x,t)
    — run ―numerical experiment‖ to                                                                      spectral methods
          simulate E(x, t), H(x, t)
                                                         df   f ( x  x)  f ( x  x)
                                                         dx
                                                            
                                                                         x
                                                                                         O(x 2 )      +
• Frequency-domain linear response
 — start with harmonic current J(x, t) = e–it J(x)                                                         +
                                                    
 — solve for steady-state harmonic fields E(x), H(x)                                                            + …..
 — involves solving linear equation Ax=b
                                                                 finite elements
• Frequency-domain eigensolver                                                                 in irregular ―elements,‖
                                                                                               approximate unknowns
 — solve for source-free harmonic eigenfields
                                                                                               by low-degree polynomial
             E(x), H(x) ~ e–it
 — involves solving eigenequation Ax=2x

                                                                                                     boundary-element methods
                                                                                                            discretize only the
                                                                                                            boundaries between
                                                                                                            homogeneous media
         A Maxwell Eigenproblem
          1         
 E        H i H                            First task:
          c t       c                      get rid of this mess
          1         0   
 H        E  J  i E
          c t           c

dielectric function (x) = n2(x)


           1                      2
                                                   + constraint
      H    H
              c                                H  0
    eigen-operator            eigen-value   eigen-state
Electronic & Photonic Eigenproblems
         Electronic                       Photonic
 2 2               1         2
 
           
        V   E     H    H
 2m                         c 
  nonlinear eigenproblem          simple linear eigenproblem
 (V depends on e density ||2)        (for linear materials)


                                 —many well-known
                                   computational techniques

     Hermitian = real E & , … Periodicity = Bloch’s theorem…
A 2d Model System

                  dielectric ―atom‖
                     =12 (e.g. Si)


                    square lattice,
                       period a
          a
              a
                               E
                     TM
                               H
            Periodic Eigenproblems
 if eigen-operator is periodic, then Bloch-Floquet theorem applies:

                                        
                                       i k x t   
       can choose:    H(x ,t)  e                    Hk (x )

                              planewave
                                                    periodic ―envelope‖

Corollary 1: k is conserved, i.e. no scattering of Bloch wave

Corollary 2: H k given by finite unit cell,
             so  are discrete n(k)
     A More Familiar Eigenproblem
                         =1

 a                             = 12                             y
                                                                     x
                                       band diagram / dispersion relation
                                   
find the normal modes                          light cone
                                        (all non-guided modes)
   of the waveguide:
                        i(kxt)
H(y,t)  H k (y)e
  (propagation constant k
         a.k.a. )
                                                                         k
    Solving the Maxwell Eigenproblem
                                                                   1                          n 2
 Finite cell  discrete eigenvalues n                  ik    ik  H n                      Hn
                                                                                             c   2


  Want to solve for n(k),
& plot vs. ―all‖ k for ―all‖ n,
       QuickTime™ see this
       are needed toand a picture.
       Graphics decompressor
                                                          constraint:     ik H  0
                                                 
        1


                                                         where:         H(x,y) ei(kx – t)
       0.9
       0.8
       0.7
       0.6
       0.5




                                                      
       0.4
                        Photonic Band Gap
       0.3
       0.2
                                      TM bands
       0.1
        0




                     1           Limit range of k: irreducible Brillouin zone

                     2          Limit degrees of freedom: expand H in finite basis

                     3           Efficiently solve eigenproblem: iterative methods
Solving the Maxwell Eigenproblem: 1
        1   Limit range of k: irreducible Brillouin zone
             —Bloch’s theorem: solutions are periodic in k
                                                 M

                                                      2     ky
     first Brillouin zone                G        X    a
                                                                  kx
= minimum |k| ―primitive cell‖

                                            
                     irreducible Brillouin zone: reduced by symmetry

        2   Limit degrees of freedom: expand H in finite basis

        3   Efficiently solve eigenproblem: iterative methods
Solving the Maxwell Eigenproblem: 2a
      1   Limit range of k: irreducible Brillouin zone
      2   Limit degrees of freedom: expand H in finite basis (N)

                   N
  H  H(xt )   hm bm (x t )               ˆ H  2 H
                                     solve: A
                   m1


          finite matrix problem:   Ah   Bh  2



   f g   f* g                 ˆ
                         Am  bm A b              Bm  bm b

      3   Efficiently solve eigenproblem: iterative methods
Solving the Maxwell Eigenproblem: 2b
                1   Limit range of k: irreducible Brillouin zone
                2   Limit degrees of freedom: expand H in finite basis
                      — must satisfy constraint: (  ik) H  0

     Planewave (FFT) basis                               Finite-element basis
                                                                            constraint, boundary conditions:
  H(x t )   HG e            iGxt
                                                                            Nédélec elements
                     G                                                       [ Nédélec, Numerische Math.

                HG  G  k  0
                                                                                   35, 315 (1980) ]
  constraint:
                                                                                nonuniform mesh,
 uniform ―grid,‖ periodic boundaries,
                                                                             more arbitrary boundaries,
       simple code, O(N log N)               [ figure: Peyrilloux et al.,
                                                 J. Lightwave Tech.
                                                  21, 536 (2003) ]
                                                                            complex code & mesh, O(N)

              3   Efficiently solve eigenproblem: iterative methods
Solving the Maxwell Eigenproblem: 3a
     1   Limit range of k: irreducible Brillouin zone
     2   Limit degrees of freedom: expand H in finite basis
     3   Efficiently solve eigenproblem: iterative methods

                        Ah   Bh  2


    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)
Solving the Maxwell Eigenproblem: 3b
     1   Limit range of k: irreducible Brillouin zone
     2   Limit degrees of freedom: expand H in finite basis
     3   Efficiently solve eigenproblem: iterative methods

                        Ah   Bh  2


    Many iterative methods:
           — Arnoldi, Lanczos, Davidson, Jacobi-Davidson, …,
              Rayleigh-quotient minimization
Solving the Maxwell Eigenproblem: 3c
               1   Limit range of k: irreducible Brillouin zone
               2   Limit degrees of freedom: expand H in finite basis
               3   Efficiently solve eigenproblem: iterative methods

                                  Ah   Bh  2


           Many iterative methods:
                  — Arnoldi, Lanczos, Davidson, Jacobi-Davidson, …,
                     Rayleigh-quotient minimization

               for Hermitian matrices, smallest eigenvalue 0 minimizes:
                           h' Ah
―variational
 theorem‖            min
                    2
                    0
                                             minimize by preconditioned
                                              conjugate-gradient (or…)
                        h h' Bh
                        are of 2dsee
                        QuickTime™ Model
             Band Diagram needed toand athis picture. System
                        Graphics decompressor

   a                    (radius 0.2a rods, =12)
                                                            1
                                                           0.9




                             frequency  (2πc/a) = a / l
                                                           0.8
                                                           0.7
                                                           0.6
                                                           0.5
                                                           0.4
                                                                    Photonic Band Gap
                                                           0.3
                                                           0.2
                                                                                  TM bands
                                                           0.1
                                                            0
irreducible Brillouin zone                                      G       X        M           G
                 M
                                                                    E             gap for
   k             X                                           TM
         G                                                          H           n > ~1.75:1
  The Iteration Scheme is Important
   (minimizing function of 104–108+ variables!)
                      h' Ah
                min
                2
                0            f (h)
                   h h' Bh

Steepest-descent: minimize (h + af) over a … repeat

Conjugate-gradient: minimize (h + ad)
  — d is f + (stuff): conjugate to previous search dirs
Preconditioned steepest descent: minimize (h + ad)
      — d = (approximate A-1) f ~ Newton’s method
Preconditioned conjugate-gradient: minimize (h + ad)
      — d is (approximate A-1) [f + (stuff)]
   The Iteration Scheme is Important
     QuickTime™ see this
     are needed toand a picture.
     Graphics decompressor

            (minimizing function of ~40,000 variables)
      1000000
          100000
           10000
            1000 Ñ
                 J            E       EE EE
             100
                        J
                        Ñ
                                                 E EEEE
                                                        EE EE
                                                          E EEEEEE
                                                                           no preconditioning
% error




                           J
                           Ñ Ñ                                 EEEE
                                                                 EEEE
              10             J ÑÑ
                               J Ñ
                                                                   EEEEEEE
                                                                    E EEEEE
                                                                     EEEEEE
                                                                      E EE EE
                                J ÑÑ                                        EE
                                                                             EEE
                                                                             EEE
                                                                              EE
                                                                               EE
                                                                                EE
               1                 J J ÑÑÑÑ                                        EE
                                                                                 EE
                                                                                  EE
                                                                                   EE
                                                                                   EE
                                                                                    EE
                                          ÑÑÑ
                                    J JJ ÑÑ                                         EE
                                                                                     EE
                                                                                      EE
                                                                                      EE
                                               ÑÑÑÑÑÑ
                                        JJJJJ ÑÑÑÑÑÑ
                                                 ÑÑÑÑÑÑ
                                                      ÑÑÑÑÑ
                                                        ÑÑÑÑ                           EE
                                                                                        EEE
                                                                                        EE
                                                                                         EE
              0.1                            JJ           ÑÑÑ
                                                          ÑÑÑ
                                                            ÑÑÑ
                                                             ÑÑÑ                         EEEE
                                                                                           EE
                                                               ÑÑÑ
                                                                ÑÑÑ
                                                                ÑÑÑ
                                                                  ÑÑ
                                                                  ÑÑ                        EE
                                                                                            EE
                                                                                             EE
                                                                                              EE
                                               J                    ÑÑ
                                                                    ÑÑ
                                                                    ÑÑ
                                                                     ÑÑ                       EE
                                                                                               EE
             0.01                               J                    ÑÑÑÑ                       EE
                                                                                                EE
                                                                                                 EE
                                                                                                 EE
                                                 J                     Ñ
                                                                       ÑÑÑ
                                                                        Ñ                         EE
                                                                                                   EE
                                                                         Ñ
                                                                         ÑÑ
                                                                          Ñ                        EE
                                                                                                    E
                                                  J                       ÑÑ
                                                                           Ñ
                                                                           Ñ                         EE
                                                                                                      E
            0.001                                                           Ñ
                                                                            Ñ                         EE
                          preconditioned J        J                         ÑÑ
                                                                             Ñ
                                                                             ÑÑ
                                                                              Ñ
                                                                                                       EE
                                                                                                        E
                                                                                                        EE
                                                                                                         EE
                                                                                                          E
                                                                                                          E
           0.0001
                        conjugate-gradient J        J        no conjugate-gradient                         E
                                                                                                           EEE
                                                                                                            EE
                                                                                                            EE
                                                                                                             E
          0.00001                                    J
                                                      J
     0.000001
                    1                     10                       100                      1000

                                               # iterations
The Boundary Conditions are Tricky
                       E|| is continuous

                          E is discontinuous
                        (D = E is continuous)



                   Any single scalar  fails:

              ?    (mean D) ≠ (any ) (mean E)
                   Use a tensor :
                                           
                                             E||
                                           
                                         1 
                                   1       E
                 The -averaging is Important
      QuickTime™ see this
      are needed toand a picture.
      Graphics decompressor


          100


                     H
                         H       backwards averaging
           10                H H H
                                        H
                     J
                               B             H H H
                                                         H   H H
                                                                     correct averaging
                                    B                            H
% error




                         J   J           B        B                   changes order
                                                           averaging
                                                      B no B
                                 J B    B
            1                       J    J
                                             B
                                              B                 B     of convergence
                                                             B
                                              J
                                                  J   J
                                                                  B
                                                                      from ∆x to ∆x2
                         tensor averaging                J   J
          0.1
                                                                 J J



      0.01
                10                                                 100     (similar effects
                         resolution (pixels/period)                         in other E&M
                                                                         numerics & analyses)
         Gap, Schmap?     QuickTime™ see this
                          are needed toand a picture.
                          Graphics decompressor


                           1
a                         0.9
                          0.8
                          0.7




            frequency 
                          0.6
                          0.5
                          0.4
                                           Photonic Band Gap
                          0.3
                          0.2
                                                         TM bands
                          0.1
                           0
                               G            X           M           G



But, what can we do with the gap?
         Intentional ―defects‖ are good

         microcavities                    waveguides (―wires‖)




3D Pho to nic C rysta l with De fe c ts
             Intentional ―defects‖ in 2d
                                      QuickTime™ see this
                                      are needed toand a picture.
                                      Graphics decompressor

        (Same computation, with supercell = many primitive cells)
        QuickTime™
GraphicsGraphics a toand a picture.
QuickTime™ see this picture.
are needed toand decompressor
        are needed see this
         decompressor




a
Microcavity Blues
    For cavities (point defects)
    frequency-domain has its drawbacks:

    • Best methods compute lowest- bands,
      but Nd supercells have Nd modes
      below the cavity mode — expensive

    • Best methods are for Hermitian operators,
      but losses requires non-Hermitian
        Time-Domain Eigensolvers
             (finite-difference time-domain = FDTD)

                       Simulate Maxwell’s equations on a discrete grid,
                              + absorbing boundaries (leakage loss)

                        • Excite with broad-spectrum dipole ( ) source

                                                 

               signal processing                      Response is many
complex n                                               sharp peaks,
                      [ Mandelshtam,
             J. Chem. Phys. 107, 6756 (1997) ]        one peak per mode


decay rate in time gives loss
              Signal Processing is Tricky
                                     signal processing
                                                             complex n
                                                 ?
                                                       spectrum
        a common approach: least-squares fit of picture.
QuickTime™ see this
are needed toand a picture.
Graphics decompressor              QuickTime™ see this
                                   are needed toand a
                                   Graphics decompressor

0.8                                                    450
                                                                   E
 0.6 E E                                               400
        E
 0.4   EE
        E
                                                       350                  fit to:
      E E
 0.2 E EE    E
             E
            EE EEEE
              E E EE E EE
                                                       300      EE
                                                                                 A
                E EE E EE EEE EEEEEEEEEEEEEEEE
   0 E E E E E EEEEEEEEEEEEEEEEEEEEEEEEEEE
                     E E EEE EEEEEEEEEEEEEEE
                                                                           (   0 ) 2  G 2
                E E E                                  250
           E EE E
         E E   EE
               E
-0.2 E EE                                              200
-0.4
          E
          E
          E
          E
           E
                                                 FFT   150
                                                               E       E
      E
-0.6                                                   100       E E
     E
-0.8                                                   50    E
                                                                E   E
                                                                     E
     E
  -1 E                                                  0  EE E
                                                             E        EE E
                                                                        E E E EE EEE E EE
                                                                           E E EE E EEE EEE
     0 1 2 3 4 5 6 7 8 9 10                               0 0.5 1 1.5 2 2.5 3 3.5 4

          Decaying signal (t)                                Lorentzian peak ()
                         Fits and Uncertainty
                                                completely
       problem: have to run long enough to and athis picture. decay
are needed toand a picture.
QuickTime™ see this
Graphics decompressor              QuickTime™
                                   Graphics decompressor
                                   are needed to see

  1  E                                             40000
    E E
     E
0.8 E E
     E EE     E
              E E                                  35000
              E E
0.6 EE
        E
              E E
             E E
                  E      E
                         E
                          E   E E
                              E E
                             E E
                                                                             actual
       EE    EE EE        E   E E E E
                             E E E E               30000
0.4                      EE           E
                                      E       E
                             EE EE EE EE E    E
    EE E E   E E EE                    E EE EE     25000
0.2                      EE  E E E E E E E E EE
                                               E
  0 EEEEEEEEE             EEEEEEEEEEE              20000
                                             E E
-0.2            EE E E    E E E E E E E E EE E                           E
      E EE EE                        E EE EE E
                            E EE EE E E E
                                                   15000
                          E             E E
-0.4
           EE   EE EE           E E
                            E E E E
                          E E       E E E E
                                                   10000
-0.6  E EE
            E   E E
                 E EE
                    E
                           E E
                           E                                 signal
         E E    E E
                 E
        E E
-0.8 E E E
         E
            E    E                                 5000
                                                             portion
     E E
  -1 E                                               0E E E E E             E E E E E
     0 1 2      3    4    5   6   7   8   9   10      0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 1.4 1.5

       Portion of decaying signal (t)               Unresolved Lorentzian peak ()

           There is a better way, which gets complex  to > 10 digits
        Unreliability of Fitting Process
          Resolving two overlapping peaks is
     near-impossible 6-parameter nonlinear fit
QuickTime™ see this
are needed toand a picture.
Graphics decompressor

       (too many local minima to converge reliably)
1200
                              EE
1000
                                   sum of two peaks
                                                               There is a better
800                               E
                              E                                way, which gets
600                                                               complex 
                                                                for both peaks
400         = 1+0.033i   E        E                             to > 10 digits
200
                          E           E  = 1.03+0.025i
                         E             E
                       EE               EE
  0 EEEEEEEEEEEEEEEEEEE                   EEEEEE
                                                EEEEEEEEEEEE
  0.5 0.6    0.7 0.8 0.9      1       1.1  1.2 1.3 1.4 1.5
            Sum of two Lorentzian peaks ()
 Quantum-inspired signal processing (NMR spectroscopy):
Filter-Diagonalization Method (FDM)
               [ Mandelshtam, J. Chem. Phys. 107, 6756 (1997) ]



Given time series yn, write:        y n  y(nt)   ak e                i k nt

                                                                  k

     …find complex amplitudes ak & frequencies k
         by a simple linear-algebra problem!

                
 Idea: pretend y(t) is autocorrelation of a quantum system:
  ˆ  i  
  H                         time-∆t evolution-operator:           ˆ  eiHt /
                                                                  U
                                                                         ˆ

         t
       say:                               ˆ n  (0)
              y n   (0)  (nt)   (0) U
                                       Method
Filter-Diagonalization107, 6756 (1997) ] (FDM)
         [ Mandelshtam, J. Chem. Phys.



                             ˆ
 y n   (0)  (nt)   (0) U n  (0)        ˆ  eiHt /
                                              U
                                                     ˆ




We want to diagonalize U: eigenvalues of U are ei∆t
                              
        …expand U in basis of |(n∆t)>:

                ˆ                 ˆ ˆˆ
U m,n   (mt) U  (nt)   (0) U mUU n  (0)  y m n 1


Umn given by yn’s — just diagonalize known matrix!
   Filter-Diagonalization Summary
         [ Mandelshtam, J. Chem. Phys. 107, 6756 (1997) ]



Umn given by yn’s — just diagonalize known matrix!
  A few omitted steps:
   —Generalized eigenvalue problem (basis not orthogonal)
   —Filter yn’s (Fourier transform):
         small bandwidth = smaller matrix (less singular)


                       • resolves many peaks at once
                       • # peaks not known a priori
                       • resolve overlapping peaks
                       • resolution >> Fourier uncertainty
             Do try this at home
FDTD simulation:
     http://ab-initio.mit.edu/meep/

Bloch-mode eigensolver:
      http://ab-initio.mit.edu/mpb/

Filter-diagonalization:
        http://ab-initio.mit.edu/harminv/

     Photonic-crystal tutorials (+ THIS TALK):
           http://ab-initio.mit.edu/
                 /photons/tutorial/
     Meep (FDTD) MPB (Eigensolver)
• Arbitrary(x) — including
  dispersive, loss/gain,              • Arbitrary periodic(x) —
  and nonlinear [(2) and (3)]       anisotropic, magneto-optic, …
                                      (lossless, linear materials)
• ArbitraryJ(x,t)
• PML/periodic/metal bound.                            • 1d/2d/3d

• 1d/2d/3d/cylindrical               • band diagrams, group velocities
                                            perturbation theory, …
• power spectra • eigenmodes

        Free/open-source          • fully scriptable interface
        software (GNU)
                                  • built-in multivariate optimization,
            • MPI parallelism         integration, root-finding, …
• exploit mirror symmetries       • field output (standard HDF5 format)
                  Unix Philosophy
   combine small, well-designed tools, via files

Input text file       MPB/Meep           standard formats
                                           (text + HDF5)
Disadvantage:
— have to learn several programs
                                      Visualization / Analysis
Advantages:                                  software
— flexibility                          (Matlab, Mayavi [vtk],
— batch processing, shell scripting   command-line tools, …)
— ease of development
                   Unix Philosophy
    combine small, well-designed tools, via files

 Input text file      MPB/Meep           standard formats
                                           (text + HDF5)

    GNU Guile scripting interpreter
        (Scheme language)
                                      Visualization / Analysis
                                             software
Embed a full scripting language:       (Matlab, Mayavi [vtk],
— parameter sweeps                    command-line tools, …)
— complex parameterized geometries
— optimization, integration, etc.
— programmable J(x, t), etc.
— … Turing complete
         A Simple Example (MPB)
                    =1

a                         = 12                    y
             find the normal modes n(k)               x
                   of the waveguide:
                                     i(kxt)
             H(y,t)  H k (y)e

            Need to specify:
            • computational cell size/resolution
          • geometry, i.e. (y)
            • what k values
            • how many modes (n = 1, 2, … ?)
  A File Format Made of Parentheses
Need to specify:
 • computational cell size/resolution
(set! geometry-lattice (make lattice (size no-size 10 no-size)
(set! resolution 32)


 • geometry, i.e. (y)
 • what k values
 • how many modes (n = 1, 2, … ?)                  1 pixel




                                                             10 (320 pixels)
                                 =1

  a                                     = 12                                  y
                                                                                   x
  A File Format Made of Parentheses
Need to specify:
 • computational cell size/resolution
 • geometry, i.e. (y)
(set! geometry                                  (choose   units of a)
 (list
  (make block (size infinity 1 infinity)
         (center 0 0 0)
         (material (make dielectric (epsilon 12))))))
 • what k values
 • how many modes (n = 1, 2, … ?)

                               =1

  a                                  = 12                              y
                                                                            x
  A File Format Made of Parentheses
Need to specify:
 • computational cell size/resolution
 • geometry, i.e. (y)
                                                   (units of 2π/a)
 • what k values
(set! k-points
 (interpolate 10 (list (vector3 0 0 0) (vector3 2 0 0))))

                        (built-in function)
 • how many modes (n = 1, 2, … ?)

                                  =1

   a                                      = 12                      y
                                                                         x
  A File Format Made of Parentheses
Need to specify:
 • computational cell size/resolution
 • geometry, i.e. (y)                  …Then run:
 • what k values                        (run)
 • how many modes (n = 1, 2, … ?)
(set! num-bands 5)                      or only TM polarization:
                                        (run-tm)

                                        or only TM, even modes:
                                        (run-tm-yeven)

                          =1

  a                             = 12                    y
                                                             x
    Simple Example (MPB) Results
                       =1

a                            = 12             y
                 find the normal modes n(k)       x
                       of the waveguide:




    red = even
    blue = odd
             Do try this at home
FDTD simulation:
     http://ab-initio.mit.edu/meep/

Bloch-mode eigensolver:
      http://ab-initio.mit.edu/mpb/

Filter-diagonalization:
        http://ab-initio.mit.edu/harminv/

     Photonic-crystal tutorials (+ THIS TALK):
           http://ab-initio.mit.edu/
                 /photons/tutorial/

								
To top