OSA's Digital Library

Optics Express

Optics Express

  • Editor: C. Martijn de Sterke
  • Vol. 16, Iss. 8 — Apr. 14, 2008
  • pp: 5227–5240
« Show journal navigation

Finite-difference time-domain algorithm for modeling Sagnac effect in rotating optical elements

Chao Peng, Rui Hui, Xuefeng Luo, Zhengbin Li, and Anshi Xu  »View Author Affiliations


Optics Express, Vol. 16, Issue 8, pp. 5227-5240 (2008)
http://dx.doi.org/10.1364/OE.16.005227


View Full Text Article

Acrobat PDF (323 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

Electrodynamics in rotating optical elements has attracted much interest due to its potential application to ultra-sensitive rotating sensing. And it is important to investigate the Sagnac effect in some novel photonic structures for it may lead to a variety of unusual manifestations. We propose a Finite-Difference Time-Domain (FDTD) method to model the Sagnac effect, which is based on the modified constitutive relation in rotating frame. The time-stepping expressions for the FDTD routine are derived and discussed, and the classical Sagnac phase shift along a waveguide is calculated. Further discussions about numerical dispersion, dielectric boundary condition and perfect matched layer (PML) absorbing boundary conditions in the rotating FDTD model are also presented respectively. The theoretical analysis and simulation results prove that the numerical algorithm can analyze the Sagnac effect effectively, and can be applied to general cases with various material properties and complex geometric structures. The proposed algorithm provides a promising systematic tool to study the properties of rotating optical elements, and to accurately analyze, design and optimize rotation sensitive optical devices.

© 2008 Optical Society of America

1. Introduction

Recently, electrodynamics in rotating optical elements has attracted much interest due to its potential application to ultra-sensitive gyroscope [1–10

1. U. Leonhardt and P. Piwnitski, “Ultrahigh sensitivity of slow-light gyroscope,” Phys. Rev. A 62, 055801 (2000). [CrossRef]

]. Several works have pointed out that the Sagnac effect can be observed in arrays of coupled resonator waveguides [3

3. J. Scheuer and A. Yariv, “Sagnac Effect in Coupled-Resonator Slow-Light Waveguide Structures,” Phys. Rev. Lett.. 96, 053901 (2006). [CrossRef] [PubMed]

][4

4. C. Peng, Z. Li, and A. Xu, “Optical gyroscope based on a coupled resonator with the all-optical analogous property of electromagnetically induced transparency,” Opt. Express 15, 3864–3875 (2007). [CrossRef] [PubMed]

], and it also causes mode-splitting in some photonic crystal structures [5

5. B. Z. Steinberg, “Rotating photonic crystals: A medium for compact optical gyroscopes,” Phys. Rev. E. 71, 056621 (2005). [CrossRef]

][6

6. B. Z. Steinberg and A. Boag “Splitting of microcavity degenerate modes in rotating photonic crystals-the miniature optical gyroscopes,” J. Opt. Soc. Am. B. 24, 142–151 (2006). [CrossRef]

]. In addition, the Sagnac effect in resonant micro-cavities is studied [7

7. S. Sunada and T. Harayama, “Sagnac effect in resonant microcavities,” Phys. Rev. A 74, 021801 (2006). [CrossRef]

][8], which shows that the cavity rotation causes shifts of the resonant frequencies proportional to the rotation rate if it is larger than a certain value. Research points to the fact that the Sagnac effect leads to a variety of manifestations in some novel photonic structures, which are significantly different compared with the phase shift phenomenon observed in conventional waveguides. The magnitude of the Sagnac effect depends on the light slowing ability of the respective structure [9

9. C. Peng, Z. Li, and A. Xu, “Rotation sensing based on a slow-light resonating structure with high group dispersion,” Appl. Opt. 46, 4125–4131 (2007). [CrossRef] [PubMed]

], allowing even exponentially sensitive rotation sensing [10

10. B. Z. Steinberg, J. Scheuer, and A. Boag, “Rotation-induced superstructure in slow-light waveguides with mode-degeneracy: optical gyroscopes with exponential sensitivity,” J. Opt. Soc. Am. B 24, 1216–1224 (2007). [CrossRef]

]. These phenomena may potentially be utilized as highly sensitive gyroscope, and a general approach to model and analyze such phenomena is required.

In conventional waveguides, the Sagnac phase shift can be calculated by using an integral over the light path [9

9. C. Peng, Z. Li, and A. Xu, “Rotation sensing based on a slow-light resonating structure with high group dispersion,” Appl. Opt. 46, 4125–4131 (2007). [CrossRef] [PubMed]

]. To analyze a structure’s response in rotation frame, themost intuitive way is to follow the transfer function, whilst taking into account the Sagnac phase shift along each light path [3

3. J. Scheuer and A. Yariv, “Sagnac Effect in Coupled-Resonator Slow-Light Waveguide Structures,” Phys. Rev. Lett.. 96, 053901 (2006). [CrossRef] [PubMed]

][4

4. C. Peng, Z. Li, and A. Xu, “Optical gyroscope based on a coupled resonator with the all-optical analogous property of electromagnetically induced transparency,” Opt. Express 15, 3864–3875 (2007). [CrossRef] [PubMed]

][9

9. C. Peng, Z. Li, and A. Xu, “Rotation sensing based on a slow-light resonating structure with high group dispersion,” Appl. Opt. 46, 4125–4131 (2007). [CrossRef] [PubMed]

]. However, this method depends on some phenomenological parameters, such as the reflection coefficient and the attenuation factor, which are difficult to determine in the designing phase. Additionally, in some scattering mechanism based photonic structures (such as photonic crystal), there is no exact light path, therefore the transfer function method can not describe this kind of Sagnac effect completely.

On the other hand, an analytic model for analyzing the Sagnac effect in micro-cavity structures was proposed [5–8

5. B. Z. Steinberg, “Rotating photonic crystals: A medium for compact optical gyroscopes,” Phys. Rev. E. 71, 056621 (2005). [CrossRef]

]. For single resonant cavity, a wave-dynamical approach can be applied, and a tight-binding approximation assuming weak-coupling between micro-cavities can be used to deal with the case of multiple resonant cavities. This method works well for simple cavities with good symmetries, but it may encounter difficulties for the cavities with complex structures, for which the eigenstate function may not be resolved easily. And because of the approximation, this analytic model may also lead to some other limitations when applied to general cases.

Furthermore, the two-dimensional Green’s function theory was proposed for electrodynamics of a rotating medium [12

12. B. Z. Steinberg, A. Shamir, and A. Boag, “Two-dimensional Green’s function theory for the electrodynamics of a rotating medium,” Phys. Rev. E 74, 016608 (2006). [CrossRef]

], which represents the response to a point source. It provides a systematic tool for the general study of the spectral properties of rotating systems. However, as a spectral method, it is not convenient for analyzing transient processes. Currently, the three-dimensional theory is still being researched on.

Numerical algorithms can thus provide a vital tool to accurately analyze, design and optimize a rotation sensitive optical device. In this paper, we develop a Finite-Difference Time-Domain (FDTD) method for the simulation of the Sagnac effect which is based on the modified constitutive relation in rotating frame. This is an ab initio, time domain method which can simulate the Sagnac effect from real physical parameters, such as material property and geometric structure. It can easily be applied to 2D or 3D situations and complex structures, for studying the device’s behavior in a rotating frame.

This paper is organized as follows: In Section 2, we present our theory of the FDTD method. We derive and discuss the 2D situation, and obtain the time-stepping expressions for the FDTD routine; In Section 3, the classical Sagnac phase shift along a conventional silicon waveguide is calculated with the light path integral and the FDTD algorithm respectively, in order to verify our theory; Then in Section.4, we discuss the issues about numerical dispersion, dielectric boundary condition and perfect matched layer (PML) absorbing boundary conditions in the rotating FDTD model respectively. Finally in Section 5, a summary of this paper is provided.

2. Theory of the FDTD algorithm in rotating frame

Assume the medium rotates slowly around a fixed axis at an angular velocity Ω, which satisfies the slow velocity assumption as |ΩRc|, R is the rotating radius and c is the speed of light in vacuum. As some early works postulated [13

13. T. Shiozawa, “Phenomenological and electron-theoretical study of the electrodynamics of rotating systems,” Proc. IEEE 61, 1694–1702 (1973). [CrossRef]

][14

14. J. L. Anderson and J. W. Ryon, “Electromagnetic radiation in accelerated systems,” Phys. Rev. 181, 1765–1775 (1969). [CrossRef]

], the basic physical laws governing the electromagnetic fields are invariant under any coordinate transformation, including a no inertial one. The transformation from the stationary frame to a rotating frame is only manifested by an appropriate change of the constitutive relation. If the rotation is rather slow, no other relativistic effects can be observed. Then up to the first order in velocity the constitutive relations in the rotation frame become [6

6. B. Z. Steinberg and A. Boag “Splitting of microcavity degenerate modes in rotating photonic crystals-the miniature optical gyroscopes,” J. Opt. Soc. Am. B. 24, 142–151 (2006). [CrossRef]

][13

13. T. Shiozawa, “Phenomenological and electron-theoretical study of the electrodynamics of rotating systems,” Proc. IEEE 61, 1694–1702 (1973). [CrossRef]

][14

14. J. L. Anderson and J. W. Ryon, “Electromagnetic radiation in accelerated systems,” Phys. Rev. 181, 1765–1775 (1969). [CrossRef]

][15

15. J. Van Bladel, “Relativity and Engineering”, Springer, Berlin, (1984)

]:

D=εEc2Ω×r×H
B=μH+c2Ω×r×E
(2.1)

And Maxwell’s equations in rotation frame keep their basic forms as:

Dt=×HJJ=Js+σE
Bt=×EMM=Ms+σmH
(2.2)

Where r=xx̂+zẑ+zẑ, Ω⃗=Ω. The derivation can be started with Eq. 2.1 and Eq. 2.2. As Ω is assumed varying slowly, Ω=const, Ω/∂t=0. And the axis is fixed in rotation frame, which leads to ∂r⃗/∂t=0. Thus, we get:

Hxt=1μ[EyzEzy(Msx+σmHx)c2xΩEzt]
(2.3)
Hyt=1μ[EzxExz(Msy+σmHy)c2yΩEzt]
(2.4)
Hzt=1μ[ExyEyx(Msz+σmHz)+c2yΩEzt+c2xΩExt]
(2.5)
Ext=1ε[HzyHyz(Jsx+σEx)+c2xΩHzt]
(2.6)
Eyt=1ε[HxyHzx(Jsy+σEy)+c2yΩHzt]
(2.7)
Ezt=1ε[HyxHxy(Jsz+σEz)c2yΩHytc2xΩHxt]
(2.8)

From these modified Maxwell’s equations above, it is noticed that the electromagnetic fields can still be divided into two independent groups: Ex, Ey, Hz for TE modes and Hx, Hy, Ez for TM modes. Compared with the ordinary form of Maxwell’s Equations, some extra terms which represent the Sagnac effect are induced, whose magnitudes are in the order of c -2 and proportional to rotation velocity. Without loss of generality, we focus the discussion on TE mode, and ∂Ex/∂t, ∂Ey/∂t, ∂Hz/∂t can be resolved with discarding the c -4 terms, then it leads to:

Ext+σεEx=1εHzy+c2xΩεμ(ExyEyxσmHz)
(2.9)
Eyt+σεEy=1εHzx+c2yΩεμ(ExyEyxσmHz)
(2.10)
Hzt+σmμHz=1μExy1μEyx+c2yΩεμ(HzxσEy)
+c2xΩεμ(HzyσEx)
(2.11)

Here, Hz can be split to two components Hzx and Hzy as Berenger’s solution for the perfect matched layer (PML) [16

16. J. P. Berenger, “A perfectly matched layer for the absorption of electromagnetic,” J. Computational Physics114185–200 (1994).

], therefore, Eq. 2.11 can be written as:

Hzxt+σmμHzx=1μEyx+c2yΩεμ(HzxσEy)
(2.12)
Hzyt+σmμHzy=1μExy+c2xΩεμ(HzyσEx)
(2.13)
Fig. 1. 2D Yee lattice

Equation 2.9–2.13 are ready for discretization now. This discretization is based on the Yee lattice and leapfrog scheme, using exponential stepping, as described by Katz, et al.[18

18. D. S. Katz, E. T. Thiele, and A. Taflove, “Validation and extension to three dimensions of the Berenger PM-Labsorbing boundary condition for FD-TD meshes,” Microwave and Guided Wave Letters, IEEE , 4268–270, (1994). [CrossRef]

], and following Yee’s notations [17

17. A. Taflove, “Computational Electrodynamics: The Finite-Difference Time-Domain Method,” Artech House, Boston, (1995).

], we obtain:

Exi,j+12n+12=exp(σi,j+12Δtεi,j+12)Ex|i,j+12n121exp(σi,j+12Δtεi,j+12)Δyσi,j+12Hz|i,j+12n
+1exp(σi,j+12Δtεi,j+12)σi,j+12c2xΩμi,j+12(ExyEyxσmHz)i,j+12n
(2.14)
Eyi+12,jn+12=exp(σi+12,jΔtεi+12,j)Ey|i+12,jn121exp(σi+12,jΔtεi+12,j)Δxσi+12,jHz|i+12,2jn
1exp(σi+12,jΔtεi+12,j)σi+12,jc2yΩμi+12,j(ExyEyxσmHz)i+12,jn
(2.15)
Hzxi,jn+1=exp(σm;i,jΔtμi,j)Hzx|i,jn1exp(σm;i,jΔtμi,j)Δxσm;i,jEy|i,jn+12
1exp(σm;i,jΔtμi,j)σm;i,jc2yΩεi,j(HzxσEy)i,jn+12
(2.16)
Hzyi,jn+1=exp(σm;i,jΔtμi,j)Hzy|i,jn1exp(σm;i,jΔtμi,j)Δyσm;i,jEx|i,jn+12
+1exp(σm;i,jΔtμi,j)σm;i,jc2xΩεi,j(HzyσEx)i,jn+12
(2.17)

Equation 2.14–2.17 are the time-stepping expressions for FDTD routine. They are valid in the whole computational region, including the PML region. For the free propagation case, if (σ,σm)=0, the exponential term should be treated carefully, as in

limσ0=1exp(σΔtε)σ=Δtε,limσm01exp(σmΔtmu)σm=Δtμ
(2.18)

However, there is a problem when applying the Yee lattice to the FDTD algorithm we proposed, that is not all of the required field components are stored in the grid. Such as Hz|n i,j+1/2, which is required by Eq. 2.14, but it is not stored at point (i, j+1/2). A way to resolve the problem is to approximately express it as (Hz|n i,j+Hz|n i,j+1)/2, which are all readily available. This strategy can also be applied to Ex|n+1/2 i,j and Ey|n+1/2 i,j.

At any given time step, the centered finite-difference (central-difference) expressions are used for the spatial and temporal derivatives that are second-order accurate in space and time increment [17

17. A. Taflove, “Computational Electrodynamics: The Finite-Difference Time-Domain Method,” Artech House, Boston, (1995).

]. For instance, ∂Hz/∂y|n i,j, which is required by Eq. 2.17, can be expressed as (Hz|n i,j+1-Hz|n i,j-1)/2Δy. But it doesn’t work for the cells at the boundary because they can not be “centered”. We also need an extrapolation method to estimate their values, which are based on the spatial derivatives in internal cells. At least two adjacent internal cells are required for the simplest extrapolation, and more cells may be needed for some sophisticated extrapolation methods.

To extend the analysis to the 3D situation, substitute Eq. 2.8 to Eq. 2.3 and 2.4, and discard the c -4 terms, then the temporal differentials ∂Hx/∂t and ∂Hy/∂t can be expressed by spatial differentials. We can repeat this procedure and obtain six equations with temporal differentials on one side and spatial differentials on the other side. Then by discretizing these equations, the time-stepping equations for FDTD routine can be obtained. The discussion about the extrapolation technique can be also extended to the 3D case.

3. Modeling 2D rotating waveguide by the FDTD algorithm

The Sagnac effect shows that, when an electromagnetic wave propagates in a moving medium, it accumulates additional phase shift compared with a wave propagating in a stationary medium [11

11. E. J. Post, “Sagnac Effect,” Rev. Mod. Phys. 39, 475–493(1967). [CrossRef]

]. This phase shift depends on the scalar product of the wave propagating direction and the velocity vector of the medium, as in

Δϕ=2πv·Δlcλ
(3.1)

A well investigated configuration is that of two waves counter propagating along a circular path in rotating medium. The total phase shift can be calculated with an integral over the light path as Δϕ=4πcλlv·Δl . It is well known that this total phase shift, a) does not depend on the shape of the light path; b) does not depend on the location of the center of rotation; c) does not depend on the presence of commoving refracting medium in the light path [11

11. E. J. Post, “Sagnac Effect,” Rev. Mod. Phys. 39, 475–493(1967). [CrossRef]

]. Actually, R. Wang et al [19

19. R. Wang, Y. Zheng, and A. Yao “Generalized Sagnac Effect,” Phys. Rev. Lett. 93, 143901 (2004). [CrossRef] [PubMed]

] point out that the Sagnac effect does not only exist in circular motion, but any moving path contributes to the total phase shift. So Eq. 3.1 is the basic differential form to describe the Sagnac effect along a light path.

Based on the FDTD model, it is expected that Eq. 3.1 can be verified. Our simulation model is set up for a 2D silicon dielectric waveguide, as Fig. 2 shows. A continuous point source is put into the waveguide, and the waveguide’s width is tailored so that only single mode can propagate. In fact, multi-modes may still exist near the source, because those modes, which cannot propagate, do not decay as fast. However, theoretically, Eq. 3.1 should still hold in spite of the multi-modes.

Unlike most of other FDTD simulations, the phase information is much more interesting for studying the Sagnac effect, compared with the field amplitude. For the position at (i, j), the value of the field component at every time step can be recorded as f(nΔt), and the spectrum F(ω) for this time domain signal can be calculated by Fourier transformation, as F(ω)=∫NΔt 0 f(t)e -jωt dt⋍ΣN-1 n=0 f(nΔt)e -jωnΔtΔt. The phase angle of the complex frequency represents the phase information. Thus the phase difference between two positions A and B can be obtained from the respective spectrum immediately, as Δϕ=arg(FA(ω))-arg(FB(ω)). Indeed, how to extract the unknown frequencies and phases from discrete time signals, sometimes called “harmonic inversion”, is a well investigated problem and various algorithms are proposed to resolve it [20

20. V. A. Mandelshtam and H. S. Taylor, “Harmonic inversion of time signals,” J. Chem. Phys. 107, 6756–6769 (1997). [CrossRef]

].

Fig. 2. The field diagram for Hz
The simulation is setup for a silicon waveguide. We use 2-order Berenger PML boundary condition. The simulation condition is set to:
The number of grid in x direction: Nx=100
The number of grid in y direction: Ny=50
The number of grid for PML layer: Nb=80
The space differentials: Δxy=λ/100
The courant stability factor: S=0.5
The wavelength of the point source: 1.550µm
The position of the point source: x=5, y=25
Half-width of the waveguide: halfwidth=10
Relative permittivity of the waveguide: ε=12

Generally, the velocity can be resolved into two orthogonal directions, we notate Vx=-yΩ and Vy=xΩ. Because the Sagnac phase shift is proportional to the scalar product of the wave propagating direction and the velocity, thus Vx≠0 leads to an additional phase shift, and it is a constant when Vy varies. It is also well known that the waveguide dispersion can in no way influence the magnitude of the Sagnac effect. Different width and refractive index of the waveguide may cause different propagation modes, but the Sagnac phase shift should be unchanging.

The Sagnac phase shift along the waveguide can be written as Δϕsagnacϕrotatingϕstationary. We calculate it for a series of positions along the center of the waveguide, and the results are plotted in Fig. 3(a)–3(c). The simulation results agree well with Eq. 3.1, which proves that the Sagnac phase shift is a constant regardless of the velocity in the perpendicular direction to the waveguide, the width of the waveguide and the refractive index of the material. These exactly confirm the theory’s assumptions. Therefore the proposed FDTD algorithm can effectively model the classical Sagnac effect along the waveguide. Thus it is believed that: this algorithm can be applied to many other photonic structures with various material properties and complex geometric structures, to study their behavior in rotating frame.

Fig. 3. The Sagnac phase shift vs. distance. The reference position is x=15 and y=25. The solid line is the theoretic result from Eq. 3.1, and the marked points are simulation results from the FDTD program.
(a) halfwidth=10;ε=12;Vy=-5,0,5m/s(b) halfwidth=8,10,12;ε=12;Vy=0m/s(c) halfwidth=10;ε=10,12,14;Vy=0m/s

An extrapolation method is applied to obtain the field components which are not available in the computer’s memory. The accuracy of the model certainly depends on the extrapolation accuracy. A two point linear extrapolation technique is used in our simulation. This is a very simple implementation, and other sophisticated extrapolation algorithms can be applied to achieve higher accuracy.

It is necessary to point out that, computational mesh setting and the PML layer setting are subtle for the FDTD algorithm. The numerical dispersion which is caused by the discrete gridding may bring some negative impacts to the Sagnac phase shift. We present an analysis about what is the contribution of rotation to the numerical dispersion of FDTD algorithm in Section. 4.1. On the other hand, the PML absorption boundary insures theoretically no wave reflection [16

16. J. P. Berenger, “A perfectly matched layer for the absorption of electromagnetic,” J. Computational Physics114185–200 (1994).

][17

17. A. Taflove, “Computational Electrodynamics: The Finite-Difference Time-Domain Method,” Artech House, Boston, (1995).

], but the additional term which represents the Sagnac effect may break this law. To what extent the absorbing boundary conditions are appropriate for the new set of equations Further discussions about the PML boundary are presented in Section.4.3.

4. Analysis and discussion about the FDTD model in rotating frame

4.1. Numerical dispersion and stability

The Sagnac effect in its basic form is manifested as a phase accumulation. However, the numerical dispersion which is caused by the discrete gridding will also induce non-physical phase error, and it may bring some negative impacts to the FDTD model. Therefore, it is necessary to investigate the numerical dispersion and stability applicable to the FDTD model in rotating frame.

The discussion can be started with the modified Maxwell’s equations Eq. 2.9–2.13. For simplicity, assume the FDTD modeling space is filled with homogeneous material having no variation of ε,µ with position in the grid, and there are no magnetic or electric loss, then the system yields

Ext=1εHzy+Vyn2(ExyEyx)
(4.1)
Eyt=1εHzxVxn2(ExyEyx)
(4.2)
Hzt=1μExy1μEyx+Vxn2Hzx+Vyn2Hzy
(4.3)

Here Vx=-yΩ and Vy=xΩ, n=εrμr is the refractive index of the material. Based on the equations above, the finite-difference expressions for this TE case can be obtained. Such as Eq. 4.1, it becomes:

Ex|i,j+12n+12Ex|i,j+12n12Δt=1εHzi,j+1nHzi,jnΔy
+Vyn2(Ex|i,j+1nEx|i,jnΔyEy|i+12,j+12nEx|i12,j+12nΔx)
(4.4)

As mentioned in Section.2, some field components may not be available in computer’s memory because of the spatial and temporal grid discretization under certain algorithm. An extrapolation method is suggested to obtain these field components. If a fine enough extrapolation is applied, Eq. 4.4 is valid with good enough accuracy.

Next, we consider a plane, monochromatic, traveling-wave solution in rotating frame, as

F|I,Jn=F0exp[j(ωΔtk~xIΔxk~yJΔy)]whereF=Ex,Ey,Hz
(4.5)

where ω is the wave angular frequency, x and y are the x- and y-components of the numerical wavevector. It should be noticed that the numerical wavevector here contains the contribution from the Sagnac effect, therefore it is slightly different with the wavevector in the stationary frame.

The basic procedure for the numerical dispersion analysis involves substitution of the plain wave solution Eq. 4.5 into the finite-difference expressions Eq. 4.4, after simplification, it leads that

Ex0=ΔtHz0εΔy·sin(k~yΔy2)sin(ωΔt2)ΔtVyEx0n2Δy·sin(k~yΔy2)sin(ωΔt2)+ΔtVyEy0n2Δx·sin(k~xΔy2)sin(ωΔt2)
(4.6)

Eq. 4.2,4.3 can be treated similarly,

Ey0=ΔtHz0εΔx·sin(k~xΔx2)sin(ωΔt2)+ΔtVxEx0n2Δy·sin(k~yΔy2)sin(ωΔt2)ΔtVxEy0n2Δx·sin(k~xΔx2)sin(ωΔt2)
(4.7)
Hz0=ΔtEx0μΔy·sin(k~yΔy2)sin(ωΔt2)+ΔtEy0μΔx·sin(k~xΔx2)sin(ωΔt2)
ΔtVxHz0n2Δx·sin(k~xΔx2)sin(ωΔt2)ΔtVyHz0n2Δy·sin(k~yΔy2)sin(ωΔt2)
(4.8)

Upon substituting E x0 of Eq. 4.6, E y0 of Eq. 4.7 into Eq. 4.8, we obtain

[ncΔtsin(ωΔt2)]2=[1GΔxsin(k~xΔx2)]2+[1GΔysin(k~yΔy2)]2
whereG=1+ΔtVxn2Δx·sin(k~xΔx2)sin(ωΔt2)+ΔtVyn2Δy·sin(k~yΔy2)sin(ωΔt2)
(4.9)

To obtain the ideal physical dispersion relation in rotating frame, we consider a physical plane wave propagating in a homogenous lossless rotating medium. Because of the rotation, the Sagnac effect induces an extra phase shift, then the plane wave can be written as:

F|I,Jn=F0exp[j(ωΔtkxIΔxkyJΔyΔϕ)]
(4.10)

Here kx, ky are the wavevectors in the stationary frame without taking accout the Sagnac effect; Δϕ presents the Sagnac effect shift, as Eq. 3.1 shows, which depends on the scalar product of the wave propagating direction and the velocity vector of the medium, then it can be rewritten as:

Δϕ=2πλc(kk·V)(kk·Δr)
=1n2ω(kxVx+kyVy)(kxIΔx+kyJΔy)
(4.11)

From Eq. 4.10 and 4.11, the wavevector which contains the contribution from the Sagnac effect in Eq. 4.5 yields

k~x=G0kx,k~y=G0kywhereG0=(1+1n2ωkxVx+1n2ωkyVy)
(4.12)

In contrast to the numerical dispersion relation of Eq. 4.9, the ideal dispersion relation can be obtained from (/c)2=(kx)2+(ky)2, and Eq. 4.12, it is

(nωc)2=(k~xG0)2+(k~yG0)2
(4.13)

Equation 4.9 and Eq. 4.13 show that the two dispersion relations are identical in the limit as Δxyt approach zero. Qualitatively, this suggests that if fine enough gridding is chosen, the impact of numerical dispersion can be reduced to any desired degree. The numerical wavevector in rotating medium is proportional to the numerical wavevector in stationary frame with the factor G, which represents the Sagnac effect numerically. So that the Sagnac phase shift will not be submerged in the numerical phase error due to the algorithm’s nature.

For the FDTD model in stationary frame, the numerical waves in a 2D/3D Yee lattice have a propagation velocity which is dependent upon the direction of wave propagation. This is also true for the FDTD algorithm in rotating frame. From the expression of G, it is apparent that the numerical wave experiences an anisotropic Sagnac phase shift. These errors represent a fundamental limitation of all grid-based Maxwell’s equations’ algorithms, and can be trouble-some when modeling electronically large structures. Further work to improve the computation accuracy is required.

4.2. Dielectric boundary conditions in rotating frame

Another important question is to what extent the boundary conditions at a dielectric material boundary are preserved in rotating frame. For non-rotating cases, the boundary conditions can be derived directly from the integral form of Maxwell’s equations focusing on a small region at an interface of two media. As we discussed in Section.2, theMaxwell’s equations governing the electromagnetic fields are invariant in rotating frame, but an appropriate change of the constitutive relation is induce to represent the Sagnac effect. Because the integral form of Maxwell’s equations is independent of the constitutive relation, the basic forms of dielectric boundary conditions are also invariant in rotating frame, as

n·(D2D1)=n·(B2B1)=0
(4.14)
n×(E2E1)=0n×(H2H1)=K
(4.15)

However, if it is necessary, the modified constitutive relation should be applied to the boundary conditions to calculate a rotating model. For instance,

n·D=n·(εE)c2n·(Ω×r×H)
(4.16)
n·B=n·(μH)+c2n·(Ω×r×E)
(4.17)

Fortunately, in 2D geometries rotation case with Ω⃗=Ω, the boundary conditions can be further simplified. Form Eq. 4.15, it is notice that the continuity of the tangential E and H fields still preserved across the dielectric interface. If Ω⃗ is in the z-direction, for a small region at the media interface,

(Ω×r×H)·dS=0;(Ω×r×E)·dS=0;
(4.18)

so that Eq. 4.14 can be written as:

n·(εE2εE1)=n·(μH2μH2)=0
(4.19)

They are the same as the boundary conditions in the stationary frame. The 2D geometries rotation does not affect boundary conditions at the dielectric interface. With this result, those well tested legacy codes for modeling stationary structures can also be extended to dealing with rotating structures, with the new set of the time-stepping equations as Eq. 2.14–2.17.

4.3. Absorbing boundary conditions (ABC) in rotating FDTD model

Recall that the ABCs used in the rotating FDTD model were originally developed for non-rotating cases. It is also necessary to consider whether they are appropriate for the new set of equations of rotating medium.

To discuss the PML absorbers in rotating medium analytically, consider the time-harmonic forms of Eq. 2.3–2.8 in a rotating medium, where no electronic or magnetic source exist, it leads to:

jωεsyE~x=H~zy+jωc2VyH~z
(4.20)
jωεsxE~y=H~zxjωc2VxH~z
(4.21)
jωμsx*Hzx~=Ey~xjωc2VxEy~
(4.22)
jωμsy*Hzy~=Ex~yjωc2VyEx~
(4.23)

Here the notations are simplified by introducing the variables

sw=1+σwjωε;sw*=1+σw*jωμ;wherew=x,y
(4.24)

To evaluate the ABCs, the basic procedure is to derive the plain-wave solution in the rotating PML medium and calculate the field reflection across the media interface. We consider a TE-polarized wave incident from the lossless medium half-space (Region 1, x<0), onto the material half-space having the eletric conductivity σ and magnetic loss σ* (Region 2, x>0). Substituting the expressions for ∂Ẽy/∂x and ∂Ẽx/∂y form Eq. 4.20,4.21 into Eq. 4.22,4.23, we obtain the wave-equation in rotating medium:

1sxsx*2H~zx2+1sysy*2H~zy2+jωc2Vxsxsx*H~zx+jωc2Vysysy*H~zy+ω2εμH~z=0
(4.25)

this wave-equation support the solution:

H~z=H0exp(jsxsx*kx~x)·exp(jsysy*ky~y)
(4.26)

with the dispersion relation

k~x2+k~X22ωc2Vxk~xsxsx*2ωc2Vyk~ysysy*=(nωc)2
(4.27)

Then, from Eq. 4.28, 4.20, 4.21, all fields in Region 2 can be obtained:

H2z~=H0τ·exp(js2xs2x*k2x~x)·exp(jsysy*k2y~y)
(4.28)
E2x~=H0τ[k2y~ωε2s2x*s2y+Vyε2s2yc2]·exp(js2xs2x*k2x~x)·exp(js2ys2y*k2y~y)
(4.29)
E2y~=H0τ[k2x~ωε2s2x*s2xVxε2s2xc2]·exp(js2xs2x*k2x~x)·exp(js2ys2y*k2y~y)
(4.30)

We know that the PML formulation represents a generalization of normally modeled physical media. If σ x=σ y=0 and σ*x=σ*y=0, it reduces to a lossless homogenous medium. In this case, the above dispersion relation Eq. 4.27 consists with Eq. 4.13. Therefore, in Region 1, the total field are given by

H1z~=H0[1+Γexp(2jk1x)]·exp(jk1xx)exp(jk1yy)
(4.31)
E1x~=H0[1+Γexp(2jk1x)][k1yωε1+Vyε1c2]·exp(2jk1x)]·exp(jk1xx)exp(jk1yy)
(4.32)
E1y~=H0[1+Γexp(2jk1x)][k1xωε1Vxε1c2]·exp(2jk1x)]·exp(jk1xx)exp(jk1yy)
(4.33)

where Γ and τ is the H-field reflection and transmission coefficients. As we discussed in Section. 4.2, the dielectric boundary conditions still hold for 2D rotating medium, thus the tangential E and H fields must be preserved across the x=0 interface. For non-rotating cases, the PML-ABC requires sy=s*y=1, sx=s*x,ε 1=ε 2,µ 1=µ 2, or equivalently σ y=0=σ*y,σ x/ε 1=σ*x/µ 1. Apply them to Eq. 4.28–4.33, it is noticed that the phase-matching condition k x1=k x2, k y1=ky2 still hold for the rotating case. Further, we derive the H-field reflection and transmission coefficients:

Γ=(k1xωε1k2xωε2sx*sx+Vxε2s2xc2Vxε1c2)(k1xωε1+k2xωε2sx*sxVxε2s2xc2+Vxε1c2)1
(4.34)

From above eqation, we notice that the reflectionless condition Γ=0 does hold strictly with the conventional PML-ABC. For all incident angles, it becomes

Γ=Vx2c(1s21)
(4.35)

It is apparent that the reflection coefficient is in the order of ν/c. For a slowly-rotating case, the reflection should be considerably small in quantity. Therefore, the PML-ABC which originally developed for non-rotating case still has good accuracy for the rotating medium.

5. Conclusion

In this paper, a FDTD algorithm is proposed for numerical modeling the Sagnac effect in rotating optical elements. This algorithm is based on the modified constitutive relation in rotating frame. We derive and discuss the 2D case, and obtain the time-stepping expressions for the FDTD routine. We also suggest an extrapolation technique to approximate the field components which are not available in the computer’s memory. These discussions can also be extended to the 3D case.

Further, a series of deeper issues about numerical dispersion, dielectric boundary condition and PML absorbing boundary conditions in the rotating FDTD model are discussed respectively. To analysis to what extend the FDTD model for rotating optical elements is valid.

Acknowledgments

This work is supported by NSFC with contract No. 60772003 and the National High-Tech Research and Development Program of China (863 Program) with contract No. 2006AA12Z310.

References and links

1.

U. Leonhardt and P. Piwnitski, “Ultrahigh sensitivity of slow-light gyroscope,” Phys. Rev. A 62, 055801 (2000). [CrossRef]

2.

A. B. Matsko, A. A. Savchenkov, V. S. Ilchenko, and L. Maleki, “Optical gyroscope with whispering gallery mode optical cavities,” Opt. Commun. 233, 107–112 (2004). [CrossRef]

3.

J. Scheuer and A. Yariv, “Sagnac Effect in Coupled-Resonator Slow-Light Waveguide Structures,” Phys. Rev. Lett.. 96, 053901 (2006). [CrossRef] [PubMed]

4.

C. Peng, Z. Li, and A. Xu, “Optical gyroscope based on a coupled resonator with the all-optical analogous property of electromagnetically induced transparency,” Opt. Express 15, 3864–3875 (2007). [CrossRef] [PubMed]

5.

B. Z. Steinberg, “Rotating photonic crystals: A medium for compact optical gyroscopes,” Phys. Rev. E. 71, 056621 (2005). [CrossRef]

6.

B. Z. Steinberg and A. Boag “Splitting of microcavity degenerate modes in rotating photonic crystals-the miniature optical gyroscopes,” J. Opt. Soc. Am. B. 24, 142–151 (2006). [CrossRef]

7.

S. Sunada and T. Harayama, “Sagnac effect in resonant microcavities,” Phys. Rev. A 74, 021801 (2006). [CrossRef]

8.

S. Sunada and T. Harayama, “Design of resonant microcavities: application to optical gyroscopes,” Opt. Express 15, 16245–16254 (2007). [CrossRef] [PubMed]

9.

C. Peng, Z. Li, and A. Xu, “Rotation sensing based on a slow-light resonating structure with high group dispersion,” Appl. Opt. 46, 4125–4131 (2007). [CrossRef] [PubMed]

10.

B. Z. Steinberg, J. Scheuer, and A. Boag, “Rotation-induced superstructure in slow-light waveguides with mode-degeneracy: optical gyroscopes with exponential sensitivity,” J. Opt. Soc. Am. B 24, 1216–1224 (2007). [CrossRef]

11.

E. J. Post, “Sagnac Effect,” Rev. Mod. Phys. 39, 475–493(1967). [CrossRef]

12.

B. Z. Steinberg, A. Shamir, and A. Boag, “Two-dimensional Green’s function theory for the electrodynamics of a rotating medium,” Phys. Rev. E 74, 016608 (2006). [CrossRef]

13.

T. Shiozawa, “Phenomenological and electron-theoretical study of the electrodynamics of rotating systems,” Proc. IEEE 61, 1694–1702 (1973). [CrossRef]

14.

J. L. Anderson and J. W. Ryon, “Electromagnetic radiation in accelerated systems,” Phys. Rev. 181, 1765–1775 (1969). [CrossRef]

15.

J. Van Bladel, “Relativity and Engineering”, Springer, Berlin, (1984)

16.

J. P. Berenger, “A perfectly matched layer for the absorption of electromagnetic,” J. Computational Physics114185–200 (1994).

17.

A. Taflove, “Computational Electrodynamics: The Finite-Difference Time-Domain Method,” Artech House, Boston, (1995).

18.

D. S. Katz, E. T. Thiele, and A. Taflove, “Validation and extension to three dimensions of the Berenger PM-Labsorbing boundary condition for FD-TD meshes,” Microwave and Guided Wave Letters, IEEE , 4268–270, (1994). [CrossRef]

19.

R. Wang, Y. Zheng, and A. Yao “Generalized Sagnac Effect,” Phys. Rev. Lett. 93, 143901 (2004). [CrossRef] [PubMed]

20.

V. A. Mandelshtam and H. S. Taylor, “Harmonic inversion of time signals,” J. Chem. Phys. 107, 6756–6769 (1997). [CrossRef]

OCIS Codes
(000.4430) General : Numerical approximation and analysis
(060.2800) Fiber optics and optical communications : Gyroscopes
(120.5790) Instrumentation, measurement, and metrology : Sagnac effect
(140.3370) Lasers and laser optics : Laser gyroscopes

ToC Category:
Physical Optics

History
Original Manuscript: January 25, 2008
Revised Manuscript: March 21, 2008
Manuscript Accepted: March 22, 2008
Published: April 1, 2008

Citation
Chao Peng, Rui Hui, Xuefeng Luo, Zhengbin Li, and Anshi Xu, "Finite-difference time-domain algorithm for modeling Sagnac effect in rotating optical elements," Opt. Express 16, 5227-5240 (2008)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-16-8-5227


Sort:  Year  |  Journal  |  Reset  

References

  1. U. Leonhardt and P. Piwnitski, "Ultrahigh sensitivity of slow-light gyroscope," Phys. Rev. A 62, 055801 (2000). [CrossRef]
  2. A. B. Matsko, A. A. Savchenkov, V. S. Ilchenko and L. Maleki, "Optical gyroscope with whispering gallery mode optical cavities," Opt. Commun. 233, 107-112 (2004). [CrossRef]
  3. J. Scheuer and A. Yariv, "Sagnac Effect in Coupled-Resonator Slow-Light Waveguide Structures," Phys. Rev. Lett.  96, 053901 (2006). [CrossRef] [PubMed]
  4. C. Peng, Z. Li, and A. Xu, "Optical gyroscope based on a coupled resonator with the all-optical analogous property of electromagnetically induced transparency," Opt. Express 15, 3864-3875 (2007). [CrossRef] [PubMed]
  5. B. Z. Steinberg, "Rotating photonic crystals: A medium for compact optical gyroscopes," Phys. Rev. E.  71, 056621 (2005). [CrossRef]
  6. B. Z. Steinberg and A. Boag "Splitting of microcavity degenerate modes in rotating photonic crystals-the miniature optical gyroscopes," J. Opt. Soc. Am. B.  24, 142-151 (2006). [CrossRef]
  7. S. Sunada and T. Harayama, "Sagnac effect in resonant microcavities," Phys. Rev. A 74, 021801 (2006). [CrossRef]
  8. S. Sunada and T. Harayama, "Design of resonant microcavities: application to optical gyroscopes," Opt. Express 15, 16245-16254 (2007). [CrossRef] [PubMed]
  9. C. Peng, Z. Li, and A. Xu, "Rotation sensing based on a slow-light resonating structure with high group dispersion," Appl. Opt. 46, 4125-4131 (2007). [CrossRef] [PubMed]
  10. B. Z. Steinberg, J. Scheuer, and A. Boag, "Rotation-induced superstructure in slow-light waveguides with modedegeneracy: optical gyroscopes with exponential sensitivity," J. Opt. Soc. Am. B 24, 1216-1224 (2007). [CrossRef]
  11. E. J. Post, "Sagnac Effect," Rev. Mod. Phys. 39, 475-493 (1967). [CrossRef]
  12. B. Z. Steinberg, A. Shamir, and A. Boag, "Two-dimensional Green’s function theory for the electrodynamics of a rotating medium," Phys. Rev. E 74, 016608 (2006). [CrossRef]
  13. T. Shiozawa, "Phenomenological and electron-theoretical study of the electrodynamics of rotating systems," Proc. IEEE 61, 1694-1702 (1973). [CrossRef]
  14. J. L. Anderson and J. W. Ryon, "Electromagnetic radiation in accelerated systems," Phys. Rev. 181, 1765-1775 (1969) [CrossRef]
  15. J. Van Bladel, Relativity and Engineering (Springer, Berlin, 1984)
  16. J. P. Berenger, "A perfectly matched layer for the absorption of electromagnetic," J. Computational Physics 114, 185-200 (1994).
  17. A. Taflove, Computational Electrodynamics: The Finite-Difference Time-Domain Method (Artech House, Boston, 1995).
  18. D. S. Katz, E. T. Thiele, and A. Taflove, "Validation and extension to three dimensions of the Berenger PMLabsorbing boundary condition for FD-TD meshes," Microwave and Guided Wave Letters, IEEE,  4, 268-270 (1994). [CrossRef]
  19. R. Wang, Y. Zheng, and A. Yao "Generalized Sagnac Effect," Phys. Rev. Lett. 93, 143901 (2004). [CrossRef] [PubMed]
  20. V. A. Mandelshtam and H. S. Taylor, "Harmonic inversion of time signals," J. Chem. Phys. 107, 6756-6769 (1997). [CrossRef]

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

Figures

Fig. 1. Fig. 2. Fig. 3.
 

« Previous Article  |  Next Article »

OSA is a member of CrossRef.

CrossCheck Deposited