中 華 大 學

中 華 大 學 碩 士 論 文

題目： 利用 FDTD 模擬微雷射共振腔之研究

Numerical Study of Microcavity by FDTD

系 所 別：電機工程學系碩士班 學號姓名：M09501047 許益誠 指導教授：高 川 原 博士

中 華 民 國 九十八 年 七 月

摘要 摘要 摘要 摘要

由於奈米製程技術的進步與成熟，使得雷射共振腔的尺寸能夠縮小到和光波 長的數量級相當，得以實現微共振腔雷射。精確的數值模擬可以幫助了解詳盡的 微共振腔特性，並可以在元件製作前優化設計的參數。

本論文是藉由時變的數值方法去模擬電磁波與線性和非線性增益介質的交 互作用，同樣的也探討色散介質與電磁波的交互作用。FDTD 是經由轉換時變的 馬克斯威爾方程式來獲得的一套數值方法，我們將運用這樣的技術來模擬微雷射 共振腔的數值模型。模擬的結果獲得了微雷射共振腔的增益頻譜，同時也可得其 Q-factor 和 finesse。在未來這樣的數值方法將可以發展高功率的光電元件應用在 光通訊的領域。

ABSTRACT

Recent advances in nanotechnology-based manufacturing have made it possible to manufacture optical resonators having physical dimensions of the order of the optical wavelength. Accurate numerical simulations can provide a detailed understanding of the characters of optical microcavities and allow for design optimization before devices are fabricated.

This purpose is to develop time-domain numerical algorithms for modeling electromagnetic wave interactions with linear and nonlinear optical gain media, as well as dispersions. FDTD techniques is used to study these phenomena and simulate the numerical modeling of microcavity lasers directly from the time-dependent Maxwell's equations. Simulations obtained the gain spectrum, Q-factor and finesse of microcavity resonator. The computational method will be served as a tool to develop useful components in future high-intensity photonic integrated circuits applied in optical communications.

誌謝 誌謝 誌謝 誌謝

在光電實驗室中度過了充實的研究生之旅，此時此刻感到特別高興，終於完 成了學生時期的最後一個階段並可以順利的進入社會和投入職場。

誠蒙恩師高川原老師悉心教誨與指導，藉由這樣的過程讓我以更多元的思考 來面對和分析問題，使我能夠在碩士的階段裡獲得顯著的成長，這一切的歷練將 讓我受用無窮，在此僅致上我最誠摯的謝意。同時也感謝溫勝發老師、鄭邵家老 師平時在課業上的教導與課後的關心。

感謝學長張民宗、黃豪昇、歐俊偉、張燕華，感謝同學廖偉帆、周家旭、黃 仁傑、郭修志、蘇瑞誠、于大容、林煜勝、廖苑辰、楊舜博以及學弟們對我的關 心與協助，讓我在這段日子裡能夠以非常愉快的心情來學習。

最後僅以本論文獻給我的家人和關心我的所有朋友和同學們。

目錄 目錄 目錄 目錄

摘要 摘要 摘要

摘要 ... I

ABSTRACT...II

誌謝

誌謝 誌謝

誌謝 ...III 目錄

目錄 目錄

目錄 ... IV 圖目錄

圖目錄 圖目錄

圖目錄 ... VI 第一章

第一章 第一章

第一章 序論序論序論序論... 1

1.1 簡介

簡介簡介... 1 簡介

1.2

論文大綱論文大綱論文大綱論文大綱 ... 2

第二章第二章 第二章第二章 微共振腔與雷射概述微共振腔與雷射概述微共振腔與雷射概述微共振腔與雷射概述... 3

2.1

雷射概述雷射概述雷射概述雷射概述 ... 3

2.1.1

雷射原子和雷射雷射原子和雷射泵雷射原子和雷射雷射原子和雷射泵泵泵浦浦浦浦 ... 3

2.1.2

雷射放大率雷射放大率雷射放大率 ... 3 雷射放大率

2.2

微共振腔概述微共振腔概述微共振腔概述微共振腔概述 ... 4

2.2.1

品質因子和品質因子和品質因子和 FINESSE ... 5 品質因子和

2.2.2

共振腔內電場的放大率共振腔內電場的放大率共振腔內電場的放大率 ... 7 共振腔內電場的放大率

2.2.3

光子壽命光子壽命光子壽命 ... 7 光子壽命 第三章 第三章 第三章 第三章 數值方法數值方法數值方法數值方法... 10

3.1

有限差分法有限差分法有限差分法有限差分法 ... 10

3.2

有限差分時域法有限差分時域法有限差分時域法有限差分時域法 ... 12

3.3 PLRC ... 18

3.4 ADE ... 23

3.4.1 Linear Lorentz Polarization... 25

3.4.2 Nonlinear Polarization ... 26

3.5 CPML ... 33

3.6 CPML

低頻吸收的參數控制低頻吸收的參數控制低頻吸收的參數控制低頻吸收的參數控制 ... 40

第四章第四章 第四章第四章 數值模擬與分析數值模擬與分析數值模擬與分析數值模擬與分析 ... 43

4.1

數值分析的比較數值分析的比較數值分析的比較數值分析的比較 ... 43

4.2

線性射散介質的模擬線性射散介質的模擬線性射散介質的模擬線性射散介質的模擬... 44

4.3

非線性色散介質非線性色散介質非線性色散介質非線性色散介質 ... 46

4.4

增益介質增益介質增益介質增益介質 ... 48

4.5

雷射共振腔雷射共振腔雷射共振腔雷射共振腔 ... 51

第五章 第五章 第五章 第五章 結論與未來展望結論與未來展望結論與未來展望結論與未來展望 ... 55

參考文獻 參考文獻 參考文獻 參考文獻 ... 56

圖目錄 圖目錄 圖目錄 圖目錄

圖 2.1

(

A

) 單一入射面的反射和 (

B

) 多層入射面的反射 ... 5

圖 2.2

(

A

) 微共振腔內模態的分佈 (

B

) 模態的頻譜... 6

圖 2.3

共振腔內光子轉換的過程 ... 8

圖 2.4

POPULATION

INVERSION 介於高低能階的分佈圖... 4

圖 2.5

雷射的光放大簡易圖... 4

圖 3.1 描述的是以 P 點為基準，所對應的 forward、backward、central differences ………..11

圖 3.2

以一個平面的分佈來描述離散後時間與空間相對應的關係。 ... 13

圖 3.3

描述電場與磁場在時間與空間上對應的關係。 ... 14

圖 3.4

L

ORENTZ MEDIA

的線性色散響應。 ... 19

圖 3.5

PLRC-FDTD 模擬色散介質的流程圖。... 23

圖 3.6

非線性色散響應對時間所做的圖。 ... 27

圖 3.7

ADE-FDTD 模擬色散介質的流程圖。 ... 30

圖 3.8

低頻吸收的模擬空間圖。... 41

圖 3.9

以

α

為變數的反射關係圖。... 41

圖 4.1 Lorentz media 一維的架構圖………..43

圖 4.2

PLRC-FDTD 與 ADE-FDTD 模擬的比較。 ... 44

圖 4.3

分析波在色散介質傳遞的模型圖。 ... 45

圖 4.4

ADE-FDTD 模擬電場在 L

ORENTZ MEDIA

傳遞的過程。 ... 45

圖 4.5

敘述色散的機制 ... 47

圖 4.6

非色散介質中脈衝傳遞的變化過程。... 48

圖 4.7

波傳遞於非色散介質所產生的

PRECURSOR

。 ... 48

圖 4.8

脈衝傳遞 l 後所獲得的增益頻譜。 ... 50

圖 4.9

β l

在頻域分佈的圖形。 ... 51 圖 4.10 原子放射的型態。 ... 52 圖 4.11 顯示了共振腔每微米所獲得的增益。... 54

第一章 第一章 第一章

第一章 序 序 序 序論 論 論 論

1.1 簡介 簡介 簡介 簡介

由於近代奈米光子學 [1] 的興起以及材料與製程技術的進步，使得雷射共 振腔的尺寸能夠縮小到與光波長的數量級相當。為了以低成本的方式來優化設計 的參數，將運用數值模擬來分析這樣的系統。雷射共振腔的介質分析包括了電磁 波與介質交互作用後所產生的極化現象，而這樣的現象有色散效應和非線性效應 兩種。我們將以 FDTD(finite-difference time-domain) [2]為主要的數值模擬方法。

FDTD 在 1966 年由 K.S.Yee 提出。此種運算法是經由馬克斯威爾方程式來獲 得時域的差分方程式。由於電場與磁場的正交性，所以在 FDTD 的架構裡每一 個電場(或磁場)都將有磁場(或電場)分佈在周圍，這樣空間上交錯排列的分佈如 同網格一般。運用此種離散方式將時變的馬克斯威爾旋度方程式轉換為一組差分 方程式，在時間軸上逐步推進並且求解空間的電磁場。FDTD 在空間上能對每一 個離散後的單元設定其對應的介質參數，如此一來將可以模擬複雜形狀或非均勻 介質物體的散射與輻射。

電磁波與介質的四個交互作用關鍵特性 [3]：

(1) 線性色散：介質介電常數(磁導係數)為頻率的函數，並且是在低強度的電磁 場下。

(2) 非線性：介質的變化與局部的電場、磁場強度有關。

(3) 非線性色散：介質其非線性的強度變化跟電磁波的頻率有關。

(4) 增益介質：增益介質能將外加能量轉換至介質中的電磁波，而電磁波的振幅 將會在介質內逐步的增加。這樣的增益可以是色散介質非線性介質。

為了模擬上述的不同介質，主要的做法有 PLRC [4, 5, 6] (piecewise - linear

FDTD 兩種，這兩種數值方法最大的差異性是在於處理捲積時遞迴運算的概念。

PLRC - FDTD 是在時域由片段線性來獲得以遞迴關係表達的捲積，ADE - FDTD 是在頻域描述捲積，經由傅利葉反轉換進而達到遞迴運算的微分方程。為了準確 模擬微共腔雷射，我們將以模擬線性色散介質的結果來比較這兩種方法，因為模 擬線性色散介質所運用的遞迴是由線性計算的方式獲得，適合用來做準確度的比 較。

接著將模擬微共振腔雷射，由上述比較的結果，選定適合的模擬方法，模擬 我們設計的參數，藉由模擬的結果可以了解到，各個模態所對應的增益大小為 何，由此就可以分析此種微共振腔雷射的特性。

1.2 論文大綱 論文大綱 論文大綱 論文大綱

第一章說明論文主題；第二章是微共振腔的基本原理；第三章介紹 FDTD、

PLRC 和 ADE 離散的推導過程及概念，並分析了解參數設計的要領；第四章分 別討論線性色散介質、非線性色散介質、增益介質的建構過程和如何設定參數與 模擬。在這一章的最後以模擬雷射共振腔為主，並計算其頻譜。第五章為結論與 未來的展望。

第二章 第二章 第二章

第二章 微共振腔 微共振腔 微共振腔 微共振腔與雷射 與雷射 與雷射 與雷射概述 概述 概述 概述

2.1 雷射概述 雷射概述 雷射概述 雷射概述

由參考文獻 [9] 可知，雷射廣義的來說是一種能產生光放大的元件，就好 比真空管放大器能放大音頻的電子信號一樣。在這裡光的頻寬必須是足夠大的，

因為不同類型的雷射可以放大的放射線範圍從由紅外線至可見光的範圍發展到 紫外線到 X-ray 的範圍。雷射已到達多樣式的型態變換，例如運用不同的雷射介 質，不同的原子系統，不同的激發技術。雷射所放射或放大的光束有著顯著的特 性，就像方向性，頻譜的純粹性，和強度的特性。這些特性確實構成雷射巨大的 變化，並且在應用上主導著雷射的型式。

2.1.1 雷射原子和雷射 雷射原子和雷射 雷射原子和雷射 雷射原子和雷射泵 泵 泵浦 泵 浦 浦 浦

以一個簡單的方式來說，我們將運用原子來說明雷射的介質。一個泵浦的過

程是必須要在原子被激發至量子力學所描述的高能階後才能作用的。而典型的雷 射介質是可以藉由許多方式而獲得泵浦的過程。

雷射的行為，在泵浦的過程不單單只是要原子被激發至高能階外，還有一個 限制就 population inversion 如圖 2.1 所示，其中描述的是雷射介質中高能階的原 子數大過低能階的原子數。而這樣的必要條件是由具多變性的雷射介質來獲得。

2.1.2 雷射放大率 雷射放大率 雷射放大率 雷射放大率

一旦獲得了 population inversion 的雷射介質之後，電磁的傳播將會在一個較

窄的頻寬內被放大，如圖 2.2 所示，而這樣的放大率頻寬的範圍是由原子的頻譜 上，高原子密度高能階的頻率到較低原子密度低能階的頻率所構成的範圍。

圖 2. 1 population inversion 介於高低能階的分佈圖。

圖 2. 2 雷射的光放大簡易圖。

2.2 微共振腔概述 微共振腔概述 微共振腔概述 微共振腔概述

根據文獻 [10] 可知，微共振腔是一種尺寸接近或低於光波長的光共振器。

而微米尺寸的共振器是運用兩種不同的方式來限定光波，第一種是運用兩個介電 Inversed

laser media Input light

beam

Output light beam population

energy

Laser action

Population inversion

High level

Low level

質在邊界各別建構單一的反射面，藉由這樣的方式來控制光的反射，第二種是以 總反射率滿足光波長共振的週期性結構來達到目的，這樣結構就好比布拉格反射 器或光子晶體一般。

圖 2.3 (a) 單一入射面的反射和 (b) 多層入射面的反射。

為了能夠優化微共振腔設計的參數，我們必須要先了解微共振腔的特性。這 一章首先我們要定義的是微共振腔的反射為 R，往返一次的最短路徑長為 L，並 且微共振腔擁有自己的光共振模態頻譜，這樣的圖是由光波長的間格和其他特性 所決定的。對於微共振腔模態圖，還完整的表達了半波長的共振模態，在模態的 轉換上提供了不同的空間配置。然而在傳統的模態區分下，將會漏掉相同基礎下 原本就存在的模態，進而影響其精準度。

2.2.1 品質因子和 品質因子和 品質因子和 FINESSE 品質因子和

根據文獻[10]，品質因子(Q-factor)在光共振腔裡扮演的角色是將共振頻寬的

(A) (b)

增大率參數化，其簡單的定義是共振腔頻率

ω

_c相對於 FWHM(full width at half of maximum)

δω

_c的比值：

c c

Q ω

=

δω

. (2.1)

而品質因子最主要要表達的是光能量在共振腔內衰減比例的度量。

(a) (b)

圖 2.4 (a) 微共振腔內模態的分佈 (b) 模態的頻譜。

FINESSE 在共振腔的定義，就相當於相鄰兩共振模態所對應的頻寬(free spectral range)與 FWHM 的比值。

1

c c

F R

R

ω π

δω

=∆ =

−

.

(2.2)

其中 R 為共振腔內往返一次的總反射率，而 FINESS 所呈現的是共振腔濾波特性 的度量。

T

δω

_c ∆

ω

_c

ω

2.2.2 共振腔內電場的 共振腔內電場的 共振腔內電場的放大 共振腔內電場的 放大 放大率 放大 率 率 率

光在共振時所產生的強度增加率可由下列式子獲得[10]

1 1

Intracavity

Incident

I F

I

⁻

R

⁼

π R .

(2. 3)

當反射鏡的損耗控制了 FINESSE，強度的增加率在行進波型共振腔內將會成均 勻的分佈。然而在駐波型微共振腔內，這樣的增加率將要由相對空間位置上波干 涉的形式來決定。

2.2.3 光子壽命 光子壽命 光子壽命 光子壽命

根據文獻 [11] 共振腔中與品質因子和 FINESSE 最有相關性的就是光子壽 命(photon lifetime)。光子壽命是一個描述共振腔建立或衰敗能量的時間常數，用 這樣的參數來描述共振腔是非常有用而且簡單。

圖 2.5 描述的是一個共振腔，而內部有一個類似裝載光子的包裹在前後兩端 被反射鏡限制的範圍內彈跳。假設在時間

t

= 時，包裹內有0

N 個光子，所對應

_p

的能量是

h N υ

_p。光子的數量在經過一次的往返後只剩下初始光子的 S (

S

=

R

=

R R

_{1 2})倍，

R 、

₁

R 為各別共振腔邊界的反射率，因此損耗的光子為：

₂

往返一次的損耗=

[

¹−

S N ]

p. (2.4)

在這裡 S 相當於光子倖存的因子。對於共振腔內光子數單位時間變換的比例則由 負的損耗數量除以一次往返的時間

τ

_RT。

( ) ( ) [

¹

]

p p p RT p p

RT RT

N dN N t N t S N

t dt

τ

τ τ

∆ + − −

→ = = −

∆ . (2.5) 或者是

p p

p

dN N

dt

^−

τ

. (2.6)

1

RT

p

S

τ

=

τ

=

− photon lifetime . (2.7)

由(2.6)式可以獲得一個簡單的解

( ) ( )

^{0 exp /}

p p p

N t

=

N

^−

t τ

^ . (2.8)

圖 2.5 共振腔內光子轉換的過程。

d

N

p

R

1

R

₂

2 p

R N

1 2 p

R R N

藉由上述的內容便可以達到對微共振腔初步的了解，運用上述的關係就可以簡單 的區別共振腔的特性，只要搭配精確的製程與良好的參數便可以獲得所需的光共 振元件。

第三 第三 第三

第三章 章 章 章 數值方法 數值方法 數值方法 數值方法

3.1 有限差分法 有限差分法 有限差分法 有限差分法

求得離散化的微分方程式是電磁數值模擬的第一步，所以一開始我們將需要 一個函數

^{f x 如圖} ( )

3.1 所示。接著我們可以運用 PB 兩點所構成的曲線來近似 P 點的斜率、導數、或切線，進而求得 P 點的 forward-difference 公式 [12]，

( ) (

0

) ( )

0 0

f x x f x f x

x

+ ∆ −

′
∆ . (3.1)

或者是以 AP 曲線來獲得 backward-difference 公式，

( ) ( )

0

(

0

)

0

f x f x x f x

x

− − ∆

′
∆ . (3.2)

或者是由 AB 曲線來獲得 central-difference 公式，

( ) (

0

) (

0

)

0 2

f x x f x x f x

x

+ ∆ − − ∆

′
∆ . (3.3)

或者是

( ) (

0

) (

0

)

0

/ 2 / 2

f x x f x x

f x

x

+ ∆ − − ∆

′
∆ . (3.4)

我們也可以推算出 P 點的二階導數

( ) ( ) ( )

( ) ( ) ( ) ( )

0 0

0

0 0 0 0

/ 2 / 2

1 .

f x x f x x

f x

x

f x x f x f x f x x

x x x

′ + ∆ − ′ − ∆

′′ =

∆

+ ∆ − − − ∆

 

=  − 

∆  ∆ ∆ 

(3.5)

或者是

( ) ( ) ( ) ( ) ( )

0 0 0

0 2

2

f x x f x f x x f x

x

+ ∆ − + − ∆

′′ ∆

. (3.6)

一個範圍內對一個點以近似的方式將其導數離散成為一個差分方程式，我們 都稱這樣的結果為有限差分近似法。

圖 3.1 描述的是以 P 點為基準，所對應的 forward、backward、central differences。

運用上述的過程，以有限差分近似法來獲得該函數的導數是相當直接的。然

( )

f x

x

A

P

B

x

0

x

0− ∆

x x

₀+ ∆

x

而還有一個比較普遍的做法是運用泰勒級數。所以接下來我們將對函數

^{f x} ( )

^做

泰勒級數展開

(

0

) ( )

0

( )

0

( )

²

( )

0

( )

³

( )

0

1 1

2! 3!

f x

+ ∆

x

=

f x

+ ∆

xf

′

x

+ ∆

x f

′′

x

+ ∆

x f

′′′

x

+⋅⋅⋅ , (3.7)

和

(

0

) ( )

0

( )

0

( )

²

( )

0

( )

³

( )

0

1 1

2! 3!

f x

− ∆

x

=

f x

− ∆

xf

′

x

+ ∆

x f

′′

x

− ∆

x f

′′′

x

+⋅⋅⋅, (3.8)

接著將對

f x 做泰勒級數展開的(3.7)及(3.8)式做相加的動作， ( )

(

0

) (

0

)

2

( ) ( )

0 ²

( )

0

( )

⁴

f x

+ ∆

x

+

f x

− ∆

x

=

f x

+ ∆

x f

′′

x

+

O

∆

x

. (3.9)

( )

⁴

O

∆

x

是截斷級數後的誤差。我們說這樣的誤差是四階 x∆ 或簡化的誤差。因此

( )

⁴

O

∆

x

所得到的值比

( )

^∆

^x

⁴還小，所以允許這一項忽略不計，則移項整理後可得 (3.6)式，接著再以(3.7)式減掉(3.8)式，而且忽略含有

( )

^∆

^x

³的項，這樣亦可以得(3.3) 式。

3.2 有限差分時域法 有限差分時域法 有限差分時域法 有限差分時域法

馬克斯威爾方程式是描述電磁現象的一組基本方程式，

H D J

t

∇ × =∂ +

∂

   

. (3.10)

E B t

∇ × = −∂

∂

  

. (3.11)

.

.

D E

B H

J E

ε

µ

σ

=

=

=

 

 

 

.

(3.12a,b,c)

(3.12a,b,c)式中的

ε

是介電常數，

µ

是磁導係數，

σ

是電導率。由於我們模擬的 電磁場具有行進與偏振的方向性，所以先做幾個基本的定義與假設。首先我們決 定入射的電場為 y 方向偏振，磁場是 z 方向偏振，而波行進方向是朝 x 方向前進。

接著定義空間與時間離散後的參數如圖 3.2 所描述的分配方式。藉由上述的過程 可知 FDTD 是由馬克斯威爾的微分方程式出發進行差分離散而獲得的一組時域 推進公式。

0,1, 2, 3....

0,1, 2, 3....

x i t i

t n t n

= ∆ =

= ∆ =

圖 3.2 以一個平面的分佈來描述離散後時間與空間相對應的關係。

這時我們可以將 P 點位置的電場或磁場定義如下

t

x

∆

x

∆

t

P

, 1,

i n

i

−

n i

+1,

n

, 1

i n

+

, 1

i n

−

(

^,

) (

^,

)

z z

E i x n t

∆ ∆ =

E i n

 

. (3.13)

(

^,

) (

^,

)

y y

H i x n t

∆ ∆ =

H i n

 

. (3.14)

接下來要探討的是磁場與電場之間的關係。由於為了在離散後的電磁場能滿足馬 克斯威爾方程式所表達的正交性與連續性，Yee 對電磁場的空間與時間的分佈做 了做一些幾何上的分配，這樣的分配猶如網格一般，每一個電磁分量稱之為 Yee 元胞如圖 3.3 所描述的狀況

圖 3.3 描述電場與磁場在時間與空間上對應的關係。

接著我們將可以對上述的電磁場運用到馬克斯威爾方程式，並且以

central-difference 的差分方法來獲得我們所要的 FDTD。我們要先展開極座標的

n

x E

y



H

z



H

z



H

z



H

z



i

¹

i

+2

i

+ 1 3

i

+2

n

1

n

+ 2

1

n

+

3

n

+ 2

馬克斯威爾方程式，獲得所有分量的關係式後，再依據之前的定義選取適合的偏 微方程來獲得我們所要的 FDTD。

E H µ

^∂

t

∇ × = −

∂

 

. (3.15)

 

(

Ez Ey

) (

Ex Ez

) (

Ey Ex

)

E x y z

y z z x x y

∂ ∂ ∂ ∂ ∂ ∂

∇ × = − + − + −

∂ ∂ ∂ ∂ ∂ ∂

 

, (3.16)

先將(3.15)式中磁場對時間的微分離散化，並把所對應的時間用上標

n

來表示

1 1

2 2 1

n n

H H

n

t µ E

+ −

− = − ∇ ×

∆

 



,

1 1

2 2

n n

t

n

H H E

µ

+ − ∆

= − ∇ ×

  

. (3.17)

接著再將(3.16)式代入(3.17)式

 

1 1

2 2

( ) ( )

( )

n n

n n

y x

z z

n n

n n

y x

E E

E E

x y

y z z x

H H t

E E

z x y

µ

+ −

 ∂ ∂ ∂ ∂ 

− + −

 

∂ ∂ ∂ ∂

∆  

= −  

∂ ∂

 

+ −

 ∂ ∂ 

 

 



. (3.18)

再將各個不同偏振方向的馬克斯威爾方程式分開表示

1 1

2 2 ( )

n n

n n z y

x x

E E H H t

y z

µ

+ − ∆  ∂ ∂ 

= −  − 

∂ ∂

 

 

. (3.19a)

1 1

2 2 ( )

n n

n n

x z

y y

E E

H H t

z x

µ

+ − ∆  ∂ ∂ 

= −  − 

∂ ∂

  . (3.19b)

1 1

2 2 ( )

n n

n n y x

z z

E E

H H t

x y

µ

+ − ∆  ∂ ∂ 

= −  − 

∂ ∂

 

  . (3.19c)

法拉第定律各個偏振方向所對應的離散結果如(3.19a, b, c)式所示，接著是展開安 培定律的式子

( , ) ( , )

( , )

D r t

H r t E

t E r t

t E σ

ε σ

∇ × =∂ +

∂

= ∂ +

∂

  

 

. (3.20)

 

(

Hz Hy

) (

Hx Hz

) (

Hy Hx

)

H x y z

y z z x x y

∂ ∂ ∂ ∂ ∂ ∂

∇ × = − + − + −

∂ ∂ ∂ ∂ ∂ ∂

 

, (3.21)

先把(3.20)式中電場的部分做離散的動作

1 1 1

1 2

2

n n n n

E E E E

n

t H

σ

ε ε

− −

− − −

= − + ∇ ×

∆

   



. (3.22)

接著再將(3.21)式代入(3.22)式

 

1 1

1 1

2 2

2 2

1 1

1 1

2 2

( ) ( )

1 2

( )

n n

n n

y x

z z

n n n n

n n

y x

H H

H H

x y

y z z x

E E E E

t H H

z x y

σ

ε ε

− −

− −

− −

− −

 

∂ ∂

∂ ∂

 

− + −

 ∂ ∂ ∂ ∂ 

− +  

= − +  

∆  

∂ ∂

 + − 

 ∂ ∂ 

 

   



 

1 1

1 1

2 2

2 2

1 1

1 1

2 2

( ) ( )

( ) ( ) 1

2 2

( )

n n

n n

y x

z z

n n n n

n n

y x

H H

H H

x y

y z z x

E E E E

t t

H H

z x y

σ σ

ε ε ε

− −

− −

− −

− −

 

∂ ∂

∂ ∂

 

− + −

 ∂ ∂ ∂ ∂ 

 

+ = − +  

∆ ∆  

∂ ∂

 + − 

 ∂ ∂ 

 

   



 

1 1

1 1

2 2

2 2

1

1 1

2 2

1 ( ) ( )

( )

1 1 2

1 1

2 2 ( )

n n

n n

y x

z z

n n

n n

y x

H H

H H

x y

y z z x

E t E

H H

t t z

x y

σ ε

σ ε σ

ε ε

− −

− −

−

− −

 

∂ ∂

∂ ∂

 

− + −

−  ∂ ∂ ∂ ∂ 

= ∆ +  

 

+ + ∂ ∂

∆ ∆  + − 

 ∂ ∂ 

 

 



 

1 1

1 1

2 2

2 2

1

1 1

2 2

( ) ( )

(1 ) 2

1 1

2 2 ( )

n n

n n

y x

z z

n n

n n

y x

H H

H H

t t x y

y z z x

E E

t t

H H

z x y

σ

ε ε

σ σ

ε ε

− −

− −

−

− −

 

∂ ∂

∂ ∂

 

∆ ∆ − + −

−  ∂ ∂ ∂ ∂ 

= +  

∆ ∆  

+ +  + ∂ −∂ 

 ∂ ∂ 

 

 



1 1 2 2 1

(1 )

2 ( )

1 1

2 2

n n

n n z y

x x

t t

H H

E E

t t y z

σ

ε ε

σ σ

ε ε

− −

−

∆ ∆  

−  ∂ ∂ 

= +  − 

∆ ∆  ∂ ∂ 

+ +  

. (3.23a)

1 1

2 2

1

(1 )

2 ( )

1 1

2 2

n n

n n x z

y y

t t

H H

E E

t t z x

σ

ε ε

σ σ

ε ε

− −

−

∆ ∆  

−  ∂ ∂ 

= +  − 

∆ ∆  ∂ ∂ 

+ +  

. (3.23b)

1 1

2 2

1

(1 )

2 ( )

1 1

2 2

n n

n n y x

z z

t t

H H

E E

t t x y

σ

ε ε

σ σ

ε ε

− −

−

∆ ∆  

−  ∂ ∂ 

= +  − 

∆ ∆  ∂ ∂ 

+ +  

. (3.23c)

安培定律各別偏振方向的離散結果如(3.23a, b, c)式所示，接著我們再令

E

_x、

E

_z、

H

x和

H 為零，之後再做離散的動作，就可以獲得電磁場的 FDTD。

_y

3.3 PLRC

根據參考文獻 [4, 5, 6] 描述，PLRC ( piecewise-linear recursive-convolution ) - FDTD 模擬色散介質時，最適合運用的介質參數是介電常數，因此我們將從電 通密度與電場的關係式來建構線性色散介質的方程式，

( )

0

( )

0

( ) ( )

0 t

D t E t E t d

τ

ε ε

_∞

ε τ χ τ τ

=

= +

∫

−

  

. (3.24)

接著我們將使用

^D

^{n =}

^{D n t} (

^{∆ 和}

) ^{E t}

^{n( )=}

^{E n t}

⁽ ∆ 來作為描述時間離散的符號，因此⁾ (3.24)式可寫為

( ) ( )

0 0

0

n n n t

D E E n t d

τ

ε ε ε τ χ τ τ

∆

∞

=

= +

∫

∆ −

  

. (3.25)

( ) ^t

χ

為 Lorentz media 與電場之間的響應，其時間的分佈如圖 3.4 所描述的

圖 3.4 Lorentz media 的線性色散響應。

我們將根據圖 3.4 獲得適當的數學式來描述線性色散響應，其式子如下

( )

2^{p 2p}2 ^ptsin

(

²_{p p}²

)

p p

t ε ω e

^δ

t

χ ω δ

ω δ

∆ −

= −

−

. (3.26)

根據文獻 [ 10,11] 可以將(3.25)式中電場

^{E t} ( )

^{在 t}∆ 的範圍內用常數來表達，

( ) ( )

1

n n

n

E E

E t E t n t

t

 + − 

= + ⋅ − ∆

 ∆ 

 

 

 

, (3.27)

(3.27)式這樣的表達方式即為片段線性的近似法，我們將藉由這樣的近似來獲得 對時間離散的電場。為了描述(3.25)式中的捲積，我們將以反向的時間偏移將(3.27) 式改寫為

( ) ( )

1

n m n m

n m

E E

E n t E m t

τ t τ

− − −

−  − 

∆ − = + ⋅ − ∆

 ∆ 

 

 

 

, (3.28)

將(3.28)式代入(3.25)式我們將獲得

( )

1 1

0 0

0

n n n n m m n m n m m

m

D ε ε E ε E χ E E ξ

− − − − −

∞

=

 

= +

∑

 + − 

    

. (3.29)

其中

( )

( )

( )

( )

( )

1 1

1

m t m t

m m

m t m t

d m t d

χ χ τ τ ξ t τ χ τ τ

+ ∆ + ∆

∆ ∆

= = − ∆

∫

^； ∆

∫

. (3.30a,b)

由於(3.26)式含有弦波函數，無法輕易的獲得

χ

^m和

ξ

^m，因此我們將藉由複數的 方式來重新描述(3.26)式

( ) t j

p

e

^{( αp jβp)t}

U t ( )

χ

 = −

γ

^{− +}

. (3.31)

2 2 2

p p p p p p p p/ p

α

=

δ

；

β

=

ω

−

δ

；

γ

= ∆

ε ω β

. (3.32a,b,c)

接著我們將(3.31)式和(3.32a, b, c)式一同代入(3.30a, b)將可以獲得

( ) ( )

1 ^{p j p t p j pm t}

m p p

p p

j e e

j

α β α β

χ γ

α β

− + ∆ − + ∆

−  

= −  − 

 , (3.33a)

(

^p/

)

²

{ ( )

^{1 ( p j p) t 1}

}

^{( p j p)m t}

m

p p p

p p

j t

j t e e

j

α β α β

ξ γ α β

α β

− + ∆ − + ∆

∆  

=  − ∆ +  −

−



, (3.33b)

其中“︶”是代表複數。接著我們將(3.29)式代入 Ampere’s law

1 1/2

n n

n

D D

H t

+

+ −

∇ × =

∆

 

 

. (3.34)

藉由(3.34)式我們將可以獲得電場在 FDTD 中更新的離散公式

( )

1 0 1/ 2

0

0 0 0 0

1 1

0 0

0

/

1 .

n n n

n n m m n m n m m

m

E E t H

E E E

ε ε ξ

ε ξ χ ε ξ χ

χ ξ

ε ξ χ

+ ∞ +

∞ ∞

− − − − −

∞ =

 −   ∆ 

=  + ∇ ×

− + − +

   

   

+ − + 

∑

 ∆ + − ∆ 

   

   (3.35)

為了要獲得∆

χ

^m和∆

ξ

^m我們將要先計算下列式子

( ) ( )

1 ^p j ^p t 1 ^p j ^p t

m m m m

p p

e

^{α β} p p

e

^{α β}

χ

 ⁺ =

χ

 ^{− + ∆}

ξ

 ⁺ =

ξ

 ^{− + ∆}

； . (3.36a, b)

接著我們將獲得

χ

^m_p 和

ξ

^_p^m實數的值

( ) ( )

m Re m m m

p p p p

χ

=

χ



ξ ξ



； =Re . (3.7a,b)

獲得了以上的色散參數後我們將會獲得複數的遞迴公式如下

( )

1 1

0

n n m n m n m

n m m

p p p

m

E E E

ψ χ ξ

− − − − −

=

 

=

∑

^ ∆ + ^ −^ ∆^ 

 

. (3.38)

其中

1 1

m m m m m m

p p p p p p

χ χ χ

⁺

ξ ξ ξ

⁺

∆ =  −  ∆ =  − 

； . (3.39a, b)

根據參考文獻 [4] 我們將可以運用遞迴的方式將(3.39)式改寫為

(

^{0 0}

)

^{n 0 n 1 1 ( p j p) t}

n n

p p p

E

p

E

p

e

^{α β}

ψ

= ∆

χ

− ∆

ξ

+ ∆

ξ

⁻ +

ψ

^{− − + ∆}

 

 

  

. (3.40)

接著將獲得實數

ψ

( )

n Re n

p p

ψ

=

ψ



. (3.41)

將上述的參數代入(3.35)式，我們將獲得電場的更新式

1 0 0 1/2

0 0 0 0

0 0

/

1 .

n n n

n

E ε ξ E t ε H

ε ξ χ ε ξ χ

ε ξ χ ψ

+ ∞ +

∞ ∞

∞

 −   ∆ 

=  + ∇ ×

− + − +

   

 

+  

− +

 

   

 (3.42)

圖 3.5 PLRC-FDTD 模擬色散介質的流程圖。

遞迴運算

(

^{0 0}

)

^{n 0 n1 1 ( p j p) t}

n n

p p p

E

p

E

p

e

^{α β}

ψ

= ∆

χ

− ∆

ξ

+ ∆

ξ

⁻ +

ψ

^{− − + ∆}

 

 

  

更新電場

1 0 0 1/2

0 0 0 0

0 0

/

1 .

n n n

n

E ε ξ E t ε H

ε ξ χ ε ξ χ

ε ξ χ ψ

+ ∞ +

∞ ∞

∞

 −   ∆ 

=  + ∇ ×

− + − +

   

 

+  

− +

 

   



儲存遞迴運算所需的電場。直接做變數的 對換。其方式如下

1 1

; .

n n n n

E

⁻ =

E E

=

E

⁺

   

計算磁場

 

1 1

2 2

( ) ( )

( )

n n

n n

y x

z z

n n

n n

y x

E E

E E

x y

y z z x

H H t

E E

z x y

µ

+ −

 ∂ ∂ ∂ ∂ 

− + −

 

∂ ∂ ∂ ∂

∆  

= −  

∂ ∂

 

+ −

 ∂ ∂ 

 

 



Time

Step=

1

n

+

3.4 ADE

ADE(auxiliary differential equation)運算法 [7, 8] 是將介質的極化向量藉由 附屬微分方程式來描述，離散後的微分方程將搭配 FDTD 的運算來模擬我們所 要的介質。這樣的做法不但能模擬非線性的介質外，對即將要模擬的色散介質在 運算上更有節省運算資源的好處，而且還能保持運算的精確度。本論文所討論的 介質其極化的響應以數學來說是與電場相依的捲積。因為要計算捲積的關係，所 以在運算時需要儲存從一開始到現階段的所有色散響應。如此一來將會消耗相當 大量的運算資源。ADE-FDTD 的運算法在第一步就會將原本的積分方程式轉換 成一個等價的微分方程式，這樣一來就不需要先前儲存的動作，換句話說就是運 用這樣的做法將會節省大量的運算資源。而且 ADE 的運算法在更新電場時並不 會改變 FDTD 原本的運算結構，在修正和調整時都非常方便。

首先我們來定義色散介質的極化向量，其與電通量的關係如下

D

=

ε ε

0 _∞

E

+

P

  

.

L NL

P

=

P

+

P

  

.

P

L



：線性色散極化向量，

P

NL



：非線性色散極化向量，

ω

p：pole frequency，

δ

p：damping coefficient。

其中∆

ε

_p是

ε

_{∞, p}與

ε

_{s p,} 的差值，

ε

_{∞, p}是在頻率無窮大時的介電常數，而

ε

_{s p,} 是再穩 定或零頻率時的介電常數。

3.4.1 Linear Lorentz Polarization

有了線性色散介質的響應後，我們將可以運用與電場相依的捲積運算而獲得 線性色散的極化向量

( )

0

( ) ( )

0 L t

P t E t d

τ

ε τ χ τ τ

=

=

∫

−

 

. (3.43)

為了要將(3.43)式轉換為微分方程式，第一步將要把所對應的響應(3.26)式做傅利 葉轉換



( )

2

2 2 2

p p p

p

j

p

χ ω ε ω

ω ωδ ω

= ∆

+ − . (3.44)

接著再與傅利葉轉換後的電場相乘便可以獲得傅利葉轉換後的捲積公式，因此我 們將獲得頻域相對應線性色散的極化向量



( )

^

( )

2 0

2 2

2

p p Lorentz

p p

P E

j ε ε ω

ω ω

ω ωδ ω

= ∆

+ − , (3.45)

再將(3.45)轉換如下

Lorentz

( ) (

²p 2 p ²

)

0 p ²p^

^{( )}

P ω ω

+

j ωδ

−

ω

=

ε ε ω

∆

E ω

, (3.46)

接著再對(3.46)式做傅利葉反轉換

2

2 2

2 0

2 ^{Lorentz Lorentz} .

Lorentz

p p p p

P P

P E

t t

ω

+

δ

^∂ +^∂ =

ε ε ω

∆

∂ ∂

 

 

(3.47)

(3.47)式就是我們所要的線性色散介質所對應的 ADE。

3.4.2 Nonlinear Polarization

由於非線性射散介質有兩種，所以我們將會把非線性色散極化向量區分成兩 部分，但在這之前我們還是先從極化向量與電場的關係式開始

( )

0 0

^{( )}

³

( ) ( ) ( )

²

P

NL

t ε χ E t

^∞

g t t E t dt

−∞ ′ ′ ′

=

∫

−

 

. (3.48)

( )

³

χ

0 ：third order nonlinearity，

g t ：causal response function. ( )

其中

g t 為一歸一化的函數， ( ) ∫

_−∞^∞

^{g t} ( )

^=1。而且

^{g t ( )}

⁼

^{αδ ( ) ( t}

^{+ 1−}

^{α ) g}

^Raman

^{( ) t 。}

藉由 0<

α

< 將1

g t 以百分比的方式區分為二，其中 ( ) g

Raman

( ) t 對時間的分佈如圖

3.6 我們將獲得一個適當的數學函數來描述這樣的響應

( )

¹² 2^{22 /2}

(

1

)

1 2

sin /

t Raman

g t τ τ e

^τ

t

τ τ τ

 +  −

=  

  . (3.49)

(3.49)式中

τ

₁相當於

g

Raman

( ) t 的週期， τ

₂相當於衰減控制的因子，這樣的響應稱 為 Raman effect。而另一個部分，我們稱為 Kerr effect [13]，所以當極化向量要描 述 Kerr effect 時就如(3.48)式一般，

( ) ^{( )} ( ) ( ) ( )

( ) ( )

3 2

0 0 3 3 0 0

P

Kerr

t E t t t E t dt E t

ε χ αδ

αε χ

∞

−∞ ′ ′ ′

= −

=

 

∫



. (3.50)

圖 3.6 非線性色散響應對時間所做的圖

。

因為(3.50)式對^

^E

³

( ) ^t

來說是線性的方程式，所以可以對時間直接離散化

( )

³

( )

³

1 1

0 0

n n

P

Kerr⁺ =

αε χ E

⁺ . (3.51)

相對描述 Raman effect 的極化向量是必須要做捲積的，因此就必須要先將(3.52) 式轉換為微分方程。

( )

0 0

^{( )}

³

( ) ( )

²

( )

Raman Raman

P t ε χ E t

^∞

g t t E t dt

−∞ ′ ′ ′

=

∫

−

 

. (3.52)

首先先定義函數

S t ( )

=

χ

Raman

^{( )}

³

( ) t

∗

E

²

( ) t

，再對

S t 做傅利葉轉換， ( )



( )

Raman

^{( )}

³

( )

^2

( )

S ω

=

χ ω E ω

, (3.53)



^{( )}

³

( ) ( )

0

^{( )}

^{3 2}

2 2

1 2

Raman Raman

Raman

j

Raman

α χ ω

χ ω

ω ωδ ω

= −

+ − , (3.54)

其中

2 2

1 2

2 2 2 1 2

; 1 /

Raman Raman

τ τ

ω δ τ

τ τ

= + = 。



( ) ( ) ^{( )}



( )

3 2

0 2

2 2

1 2

Raman

Raman Raman

S E

j α χ ω

ω ω

ω ωδ ω

= −

+ − , (3.55)



( ) (

_Raman² 2 _Raman ²

) (

1

)

0

^{( )}

³ _Raman^{2 2}

( )

S ω ω

+

j ωδ

−

ω

= −

α χ ω E ω

, (3.56)

( ) ^{( )}

2

2 3 2 2

2 0

2 1

Raman Raman Raman

S S

S E

t t

ω

+

δ

^∂ +^∂ = −

α χ ω

∂ ∂ . (3.57)

(3.57)式即為 Raman effect 的 ADE。接著我們再將(3.47)、(3.57)式做對時間離散 的動作

( )

( )

2 2

1 1

2 2 0

2 1

1 1

1

p p

n n n

Lorentz Lorentz Lorentz

p p

p p n

p

t t

P P P

t t

t E t

ω δ

δ δ

ε ε ω δ

+  − ∆   ∆ −  −

= ∆ +  + ∆ + 

 ∆ ∆ 

+  

 ∆ + 

 

. (3.58)

( )

( ) ^{( )} ( )

( )

2 2

1 1

3 2 2 0 2

2 1

1 1

1 .

1

Raman

n n Raman n

Raman Raman

Raman n

Raman

t t

S S S

t t

t E t

ω δ

δ δ

α χ ω δ

+  − ∆   ∆ −  −

=  + 

∆ + ∆ +

   

 

 − ∆ 

+  

 ∆ + 

 

(3.59)

( )

0

( ) ( )

Raman

P t

=

ε E t S t

 

. (3.60)

1 1 1

0

n n n

Raman

P

⁺ =

ε E

⁺

S

⁺ . (3.61)

接著再將離散後的公式代入電通密度的公式

1 1 1 1 1

0

n n n n n

Lorentz Kerr Raman

D

⁺ =

ε ε

_∞

E

⁺ +

P

⁺ +

P

⁺ +

P

⁺

    

,

( )

1 1 1 1 1

/ 0

n n n n n

Lorentz Kerr Raman

E

⁺ =

D

⁺ −

P

⁺ −

P

⁺ −

P

⁺

ε ε

_∞

    

. (3.62)

由於 Kerr effect 含有

^E

^3

( ) ^t

，所以必須要對(3.62)式做迭代，藉由迭代的方式來逼 近我們所需要的解，因此將(3.62)式改寫為

( ) ( )

1 1

1

3 2 1

0 0 0 0

n n

m Lorentz

m n

D P

E

E S

ε ε αε χ ε

+ +

+

+

∞

= −

+ +

 



 . (3.63)

1

2

m m

m

E E

E

+ +

=

 



. (3.64)

(3.64)中的 m=0,1,2,3,…..，而起始值

E

⁰ =

E

ⁿ

 

。藉由(3.63)式的方式來對電場更 新，就可以運用 ADE-FDTD 來模擬線性與非線性的色散介質。

圖 3.7 ADE-FDTD 模擬色散介質的流程圖。

遞迴運算

( )

( )

2 2

1 1

2 2 0

2 1

1 1

1

p p

n n n

Lorentz Lorentz Lorentz

p p

p p n

p

t t

P P P

t t

t E t

ω δ

δ δ

ε ε ω δ

+  − ∆   ∆ −  −

= ∆ +  + ∆ + 

 ∆ ∆ 

+  

∆ +

 

 

.

( )

( ) ^{( )} ( ) ( )

2 2

1 1

3 2 2 0 2

2 1

1 1

1 .

1

Raman

n n Raman n

Raman Raman

Raman n

Raman

t t

S S S

t t

t E t

ω δ

δ δ

α χ ω δ

+  − ∆   ∆ −  −

=  + 

∆ + ∆ +

   

 

 − ∆ 

+  

 ∆ + 

 

計算電通密度

1 1/ 2

n n n /

D

⁺ =

D

+ ∇ ×

H

⁺ ∆

t

   

更新電場

( )

( )

1 1

1

3 2 1

0 0 0 0

n n

m Lorentz

m n

D P

E

E S

ε ε αε χ ε

+ +

+

+

∞

= −

+ +

 





計算磁場



1 1



2 2

( )

( )

( )

n n z y

n n

n n x z

n n

y x

E E

x y z

E E

H H t y

z x

E E

z x y

µ

+ −

 ∂ ∂ 

 − 

∂ ∂

 

 ∂ ∂ 

∆  

= −  + − 

∂ ∂

 

 ∂ ∂ 

+ −

 

∂ ∂

 

 

 



儲存遞迴函數：

P

ⁿ⁻¹=

P

ⁿ ;

P

ⁿ=

P

ⁿ⁺¹;

S

ⁿ⁻¹=

S

ⁿ ;

S

ⁿ=

S

ⁿ⁺¹。

介紹完了色散介質的 ADE-FDTD 推導過程，接著要推導的是增益介質的 ADE-FDTD。這裡與上述有一個比較不一樣的做法，那就是運用極化電流密度來 描述被極化的介質，跟上述的極化向量不同，因為如果運用極化向量無法穩定的 運算傳遞於增益介質的電磁場，所以這一段所將運用的極化響應將會是電導率，

其所對應的時間函數如(3.65)式，雖然會有隨時間衰退的趨勢，但與(3.43)式或 (3.49)式不同處，是衰退的速度更慢，

( )

⁰

(

0

)

^/2

2

cos ^{t T}

t t e

T

σ

=

σ ω

⁻ . (3.65)

σ

0：peak value，

T ：damping coefficient，

2

ω

0：transition frequency。

接下來將與色散介質的做法相同，將要對(3.65)式做傅利葉轉換



( )



( )



( )

( )

( )

0 2

2 2 2 2

0 2 2 2

1 1

1 2

1

s

J j T

E I T j T T

I

ω σ ω

σ ω ω ω ω ω

 + 

 

= =

+ + −

 

+  

, (3.66)

(

1+

ω

0²

T

2²

)

^

J ( ) ω

+2

j T ω

2−

ω

²

T

2²=

σ

0

(

1+

j T E ω

2

)

^

( ) ω

, (3.67)

再對(3.67)式做傅利葉反轉換

(

1 0² 2²

)

2 2

J

2^{2 2}2

J

0 0 2

E

T J T T E T

t t t

ω

^{∂ ∂}

σ σ

^∂

+ + + = +

∂ ∂ ∂

  

 

. (3.68)

由於(3.68)式是二次微分的方程式，所以我們定義

F J t

≡∂

∂

 

, (3.69)

將(3.69)式代入(3.68)式

(

1 0² 2²

)

2 2 2²

F

0 0 2

E

T J T F T E T

t t

ω

^∂

σ σ

^∂

+ + + = +

∂ ∂

 

  

. (3.70)

(3.70)式即是我們所要的 ADE，接著將(3.70)式對時間做離散的動作，在這之前 如同先前定義的電磁波是 x 方向傳遞，電場是 y 方向偏振，磁場是 z 方向偏振，

經過整理後可獲得以下的式子

( ) ( )

( ) ( )

1 1

1 2 2

1 2

3 4

1 1

2 2

n n

n n

z z y

n n

y

F i A H i H i A E i

A J i A F i

+ +

+     

=   + −  − +

   

 

+ + . (3.71)

( ) (

0 2

)

1

4

ts i t

2

T

A x

σ β

∆ ∆ +

= ∆ , (3.72)

( )

0 2

8

s i t

,

A ε σ

β

= ∆ (3.73)

(

0² 2²

) ( ) (

0 2

)

3

4 2 1 2

,

t T s i t T

A ε ω σ

β

 

∆  + + ∆ + 

= − (3.74)