以時域有限差分方法求解三維馬克斯威爾方程組，探討手機于頭部之電磁波之比吸收率之分布

(1)

國立臺灣大學工學院工程科學及海洋工程學系碩士論文

Department of Engineering Science and Ocean Engineering College of Engineering

National Taiwan University Master Thesis

以時域有限差分方法求解三維馬克斯威爾方程組，探討手機于頭部之電磁波之比吸收率之分布

Calculation of Mobile Phone Induced Specific Absorption Rate Prediction in Phantom Head by a 3D FDTD code for

Maxwell’s equations

林緯皓 Wei-Hao Lin

指導教授：許文翰博士

Advisor: Tony Wen-Hann Sheu, Ph.D.

中華民國 107 年 7 月

July, 2018

(2)

(3)

誌謝　　　

光陰匆匆，霎那間，兩年的碩士班生涯也接近了尾聲。兩年的汗水、淚水與歷練，逐一記錄在這本碩士論文裡。

首先，我要將這兩年的一切成果與榮耀，歸功於 - 我的父親 - 林宏鈞先生，與我的母親 - 鄧桂霞女士。從我過去的不懂事到大學畢業，到當兵入伍，到退伍後出國長居，到考上台大，到我取得碩士學位。兒時的叛逆無知，年少輕狂，回首一望彷彿昨日。然而，父母從沒放棄過教養我、並導正我的價值觀，今日，我才有機會見到自己的蛻變。小時候，母親常常開玩笑的問我，我長大後會養爸媽嗎，

當時不免擔心自己能有多少能力讓爸媽過好生活，現在終於可以很肯定的回答這個問題，讓爸媽未來的退休生活遊山玩水，衣食無虞，共享天倫之樂。

同時，我也要特別感謝我的指導教授 - 許文翰恩師。記得當初來到台大，

對許文翰老師的研究相當有熱忱，唯恐基礎不佳，老師擔憂我是否能負荷，所幸，三顧茅廬之下，老師看見了我的決心，願意破例超收我這個大學成績不佳，

又看起來混混的學生。然而，這兩年的時光裡，許文翰老師，就像一個嚴父的角色，為了讓底子不好的我能夠步上軌道，給予我比其他同儕更多的關心與壓力，

還能有幸跟老師一同前往青島做研討會報告，有鑑於老師願意給予我這樣的環境與充足的資源下，最初懵懂的我，才能快速的在碩一下學期上了軌道，順利的找到方向，最終完成了我的碩士論文，發生在我身上如此大的改變，除了感念老師對我的知遇之恩以及教養之恩，師恩浩蕩，著實無以回報，萬幸能在 TWSIAM 研討會中獲得了一個傑出碩士海報論文獎，為師爭光。惟學生畢業前，僅僅能將此獎項，獻給我的恩師-許文翰老師，以表學生的感念之情。

謝謝共同指導老師 - 李佳翰老師，願意不吝共同指導我，在每周的 meeting 上，給予我在研究上專業領域的幫助與建議，也感謝李佳翰老師的緣故，讓我有機會能和張宏鈞老師學習，謝謝李佳翰老師兩年來的照顧。

謝謝電信所 - 張宏鈞老師，在百忙之中能夠抽空給予我研究上的幫助與教學，

更在口試期間擔任學生的口試委員，給予學生相當多的建議與肯定！讓我在研究上有更廣闊的視野。

謝謝我的妹妹林筱函，在我在台北求學的這段期間，下班之後都有好好陪伴爸媽！希望我經歷的故事與轉變，也能夠給予我妹妹在人生道路上一些方向與借鏡！

你是我最優秀的妹妹。

(4)

在這個使我蛻變良多的環境裡，首先我要先感謝很照顧我的高仕超學長，在我剛進實驗室人生地不熟的時候，鼓勵我好好念書，也分享自己的故事，告訴我不要受限於自己的過去，並開創自己的人生！超哥除了是我的學長，也是我人生的前輩，常常受你照顧，一日學長終生學長！

接著我要感謝同樣很照顧我的柳冠碩學長，在實驗室與研究上，各種大大小小的事情中，幾乎透過請教你，都可以得到最有效率且最正確的答案，使得我也想成為 Rex 學長這樣能夠獨當一面解決問題的人！所以後來才能漸漸的憑自己的力量完成我的研究，謝謝 Rex 學長一路上的幫助與鼓勵！

同樣要特別感謝從 SCCS 實驗室畢業的張育瑋學長，一路上無條件的付出，

給我專業知識上的協助，以及未來我的人生規劃上也給我很多很大的建議！在我口試期間一個禮拜，在百忙之中還願意撥空為我的口試進行預備練習，您真的是我人生中，遇過最無私付出，又有想法，並且有決心與毅力去做好每一件事的學長兼益友兼恩師！未來人生的道路，除了當一名好學弟之外，我一定會成為一名有能力回饋您的恩情的好朋友的。

另外，也要感謝實驗室畢業的林樂學長，無私的分享給我未來就業的方向與建議，希望未來出去工作我有能力能夠幫上學長的忙。

特別感謝楊世銓同學，在我剛來到這個環境，就把同學們都介紹給我認識，

讓我很快的融入工科 B01 的大家庭，把我加入了你們班的群組裡，彷彿我就像同班同學一樣，另外在實驗室裡，你也常常勉勵我，給我許多研究上的資源與建議，

一起生出很厲害的研究！其中最重要的是，你在實驗室準備的躺椅，結果都是我在用，多少個我和老鼠一起度過的夜晚，有你購置的這張躺椅在，如今，我才能在今年跟你一起順利畢業，以後人生的道路要一起加油，繼續戰鬥下去！我的好朋友！

感謝張育維同學在我碩一的時候，無數個夜晚跟我一起熬夜奮鬥到四點多，

很懷念那段我們熬夜跟研究奮鬥加上兩支啤酒的日子，真是一段前無古人後無來者的時光。你的努力與研究的天分，還有你對人生的態度，都激勵我堅持著往前走下去，你也是我學習的對象，謝謝有你跟楊世銓同學在這裡的照顧，才有現在的我。

同時，感謝周書聖同學，能夠認識你這樣有夢想的朋友，是相當難得的，碩一上也熱心的分享關於電磁學的細節觀念給我，使我如沐春風，後來，還記得聽你深聊過關於夢想與現實的故事，很讓我感動，每個禮拜二的光機電組 Meeting，

兩年中，與你一起努力的過程也歷歷在目，未來畢業以後，也要常保持聯絡，謝謝你！

(5)

謝謝實驗室去德國念書的王智達學長，您在實驗室的期間，也給了我許多鼓勵與建議，達哥的心路歷程，也給了我很多很多遇到困難時之正面的想法。

感謝碩一上時翁凡學長與高千渝學姊把我引薦進入 SCCS LAB 中！沒有你們的出現，我也不會在這麼好的環境中學習與成長。

謝謝吳政道同學的教導一些關於 coding 方面的專業知識，讓我受益良多。謝謝聶廷叡同學、王志明同學、羅鈺棋同學加入了 SCCS LAB，為實驗室增添了許多新的氣象。

感謝我的室友物理所的林合俊同學，在我無數個深夜晚歸，但是並沒有因此而覺得我是不良室友，相反的我們反而很照顧彼此，也常常請教你一些關於國外學校與就業方面的事，相當感謝跟你住在同一間寢室的時光！良好的室友合俊！

讓我有良好的休息環境，迎戰每一個明天的到來。

特別感謝 B01 的蔡宜芳同學，沒有你的 Latex 格式照應，我可能要花兩倍的時間才能完成我的碩士論文，以及修資結的時候給我的資源與觀念上的幫助，讓我對資料結構有更進一步的了解！感謝 B01 的江明謙同學，遠在西雅圖工作的你的故事，讓我看見會玩又會念書的表率，我會努力實現跟你的約定，美國見！謝謝 B01 的魏丞鴻同學，在工科所畢業一肩擔起了負責人身分，讓大家能夠無憂的去著手自己的碩士論文。感謝光電所李家豪同學，從碩一上認識到現在一起修課、

一起去社團認識新朋友、一起畢業找工作！在兩年研究所中，許多個煩悶的日子，

常常找你一起吃消夜聊天消除苦悶，謝謝你照應，才讓我度過那些苦難與低潮的時光。感謝電子所徐瑋澤同學，跟我一起在念書念累了相約健身房，還有一起修課的時光，從你身上也學到很多，謝謝你！

謝謝光電所謝承哲同學，在就業面試時給予我相當多的建議與幫助！未來再一起努力互相照應！謝謝化工系的曾耀賢同學、中文系的翁翊烋同學與會計系的吳妍庭同學，能遇見有相同興趣的你們，在我碩班兩年的生涯中，能夠有一個情感的抒發去處！在社團認識你們，是我在台大兩年當中最大的收穫。

僅將此碩士論文獻給我的家人 - 林宏鈞先生、鄧桂霞女士與妹妹筱函。

(6)

摘要

本論文是在非交錯網格上發展一三維時域有限差分法 (FDTD)，以求解馬克斯威爾方程。本文的方法是在時域內，在滿足高斯定律 (即電場和磁場零散度條件) 的架構下求解法拉第定律和安培定律。

本文所提出的數值方法能在時間上和空間上保有相當好的理論收斂斜率，且能有效地減少實解相速度與數值相速度之間的誤差，而得以顯著地降低了因時域有限差分所造成的數值色散誤差以及各向異性誤差。本研究證實了所提出的數值方法在具辛結構與色散關係上皆具有良好的保持性，尤其在針對經長時間馬克斯威爾方程的數值模擬後，

其效果尤為顯著。

本文進而將此數值方法針對人體頭部進行預測及數值分析其暴露在手機輻射 (RF) 下之特定比吸收率 (Specific Absorbtion Rate) 的電磁場與 SAR 場的在頭部各器官組織的分布情形。

人體在使用手機進行通話時，通常將手機聽筒貼置在左耳或右耳上，使得頭部將與手機直接貼觸，直接接受由手機天線發射出的低強度射頻電磁場 (RF-EMF) 曝曬。然而，電磁曝曬的強度，將與手機種類以及手機輻射功率和作用頻段相關聯。本文選用複雜幾何之 Apple iPhone4-like 模型，並與複雜幾何頭部組織進行電磁曝曬分析模擬，使用顯式非交錯 (或稱並列) 網格方法進行模擬計算。此方法相當適合使用多圖形處理單元 (GPUs) 平行計算，透過增加更多圖形處理單元減少計算時間，以換取計算空間之網格密度。

由於馬克斯威爾方程組屬於完全可積之方程，因此，我們採用具辛結構之 Runge-Kutta 方法來逼近時間導數項，並且保持馬克斯威爾方程組能量守恆的性質；同時透過最小化數值色散關係式與色散關係式之間的差，以減少數值色散誤差。結果顯示，所模擬行動電話的數值

(7)

結果與實驗測量值相當接近，顯示本文所使用之數值方法，可以準確的預測出低頻射頻場對人腦的影響。

關鍵字：馬克斯威爾方程; 時域有限差分法; 非交錯網格; 人體特定比吸收率 (SAR) ; 實解與數值相速度; 手機射頻電磁場分析; 數值色散關係式

(8)

Abstract

An explicit finite-difference scheme for solving the three-dimensional Maxwell’s equations in non-staggered grids is presented in time domain. Our aim is to solve the Faraday’s law and Ampère’s law in time domain under the constraint of Gauss’law.

The numerical method presented in this paper can maintain a fairly good theoretical convergence slope in time and in space. It can effectively reduce the error between the actual solution phase velocity and the numerical phase velocity by dispersion relation analysis. with the concept of phase velocity preserving, this numerical method can significantly reduce the numerical dispersion error and anisotropy error. This study confirms that the proposed numerical method can retain on symplectic structure and dispersion relation- ship.

Exposure to mobile (or cell) phone radiation will be numerically investi- gated in human head by solving the Maxwell’s equation. Our aim is to get the distribution of the electrical field in the calculation of Specific Absorption Rate (SAR). Cell phone handset is normally placed over left/right ear. Expo- sure to low-intensity Radio Frequency–Electro/Magnetic Fields (RF-EMF) from cell phone is therefore a well localized issue. Moreover, the accompa- nied electrical field takes its highest magnitude in brain regions closest to cell phone antenna. The degree of exposure depends on the type of cell phone being used. As a result, Apple iPhone4 and a phantom head are chosen in this three-dimensional simulation of Maxwell’s equations. For performing a computationally effective simulation of Maxwell’s equations, calculation of Maxwell’s solution will be performed in non-staggered (or collocated) grids using the explicit finite difference scheme. As the result of the employed

(9)

explicit discretization scheme in non-staggered grids, this simulation can be suitably executed in parallel on Graphic Process Units (GPUs), thereby re- ducing a dramatic amount of computing time. Maxwell’s equations belong to a class of completely integrable equations. Symplectic Runge-Kutta tem- poral scheme is therefore adopted to approximate time derivative terms so as to be able to preserve the embedded Hamiltonians and invariants embedded in Maxwell’s equations. In addition, the introduced numerical dispersion error is reduced by minimizing the difference between of numerical and exact dispersion relation equations. As a result, the emitted low-frequency radio frequency fields can be accurately predicted.

Keywords : M axwell^′s equations; non−staggered grids; F inite difference time domain methods; N on− staggered grid; Specific absorption rate;

dispersion relation equation; Exact and numerical phase velocities; M obile radio f requency electromagnetic f ield analysis; N umerical dispersion relation equations.

(10)

圖目錄

Fig 2.1 CPML 吸收邊界之示意圖。 . . . 18

Fig 2.2 求解電磁波傳遞問題的方程及計算空間之示意圖。 . . . 18

Fig 3.1 交錯式之 Yee 網格系統。 . . . 26

Fig 3.2 本文所採用之非交錯式網格系統。 . . . 26

Fig 3.3 非交錯網格內部節點之示意圖及編號。 . . . 26

Fig 3.4 座標系示意圖。 . . . 27

Fig 4.1 平面波在三維情況下傳遞時，定義出天頂角 (Zenith angle) θ，方 位角 (Azimuth angle) ϕ。 . . . 34

Fig 4.2 固定 Cr 為 0.2 時，在不同積分範圍下之色散關係之實解與數值 解的比較。 . . . 34

Fig 4.3 固定積分範圍 ^−π₂ ∼ ^π₂ 時，在不同 Cr 值情況下之色散關係之實 解與數值解的比較。 . . . 35

Fig 6.1 吾人所討論之實際問題描述。 . . . 40

Fig 6.2 計算空間中節點 P 與三角面三點 ABC 之關係圖。 . . . 41

Fig 6.3 真實手機的各個部件。 . . . 42

Fig 6.4 計算空間中，手機各個部件的建模。 . . . 43

Fig 6.5 高斯脈衝函數之時域波型。 . . . 46

Fig 6.6 程式計算之流程圖。 . . . 48

Fig 6.7 iPhone4 316L 不鏽鋼邊框之模型。 . . . 51

Fig 6.8 iPhone4 UMTS/GSM 頻段之參考電路圖。 . . . 52

Fig 6.9 人體頭部之結構圖。 . . . 54

Fig 7.1 Apple iPhone4-like 的 316L 不鏽鋼金屬外框天線輻射部的結構，于 T=100∆t 到 T=500∆t 在五個不同時間點的電場傳播情形。 . . . 57

Fig 7.2 頭部組織內手機輻射之電場強度分布之模擬結果。 . . . 59

(13)

Fig 7.3 大腦及腦幹組織之電場強度分布情形。 . . . 60

Fig 7.4 小腦組織之電場強度分布情形。 . . . 61

Fig 7.5 顱骨之電場強度分布情形。 . . . 62

Fig 7.6 肌肉組織之電場強度分布情形。 . . . 63

Fig 7.7 皮膚組織之電場強度分布情形。 . . . 64

Fig 7.8 由左耳到右耳，所考慮之人體組織層的截面尺寸、導電率大小以及電場強度分布。 . . . 65

Fig 7.9 人體組織最大 SAR 發生處在總模擬時間為 5000∆t 的電場大小 值|E|。 . . . 66

Fig 7.10 移去顱骨後，頭部組織內手機輻射之電場強度分布之模擬結果。 . . . 68

Fig 7.11 移去顱骨後，大腦及腦幹組織之電場強度分布情形。 . . . 69

Fig 7.12 移去顱骨後，小腦組織之電場強度分布情形。 . . . 70

Fig 7.13 移去顱骨後，肌肉組織之電場強度分布情形。 . . . 71

Fig 7.14 移去顱骨後，皮膚組織之電場強度分布情形。 . . . 72

Fig 7.15 移去顱骨後，由左耳到右耳，所考慮之人體組織層的截面尺寸、導電率大小以及電場強度分布。 . . . 73

Fig 7.16 大腦內，SAR 最大值發生位置，在 5000∆t 的時間範圍內，將原本的模擬結果和移除顱骨後的模擬結果進行比較。 . . . 74

Fig 7.17 小腦內，SAR 最大值發生位置，在 5000∆t 的時間範圍內，將原本的模擬結果和移除顱骨後的模擬結果進行比較。 . . . 74

Fig 7.18 大腦及腦幹組織於各單位格點於計算時間為 5000∆t，比較有無顱骨時的電場強度 |E| 之分布結果。 . . . 75

Fig 7.19 小腦組織於各單位格點於計算時間為 5000∆t，比較有無顱骨時的電場強度 |E| 之分布結果。 . . . 76

Fig 8.1 在 log 尺度下，頭部組織 SAR 分布之模擬結果。 . . . 78

(14)

表目錄

Table 4.1 經長時間 (即 T 為 30(s)) 計算後，在 Cr = 0.2 及 Cr = 0.05 情

況下計算誤差以及所需 CPU TIME (s) 之比較。 . . . 33

Table 4.2 取天頂角 (Zenith angle) θ 為 0^◦和 90^◦、30^◦和 60^◦以及 45^◦，配 合方位角 (Azimuth angle)ϕ 為 0^◦ 和 90^◦、6^◦ 和 84^◦、9^◦ 和 81^◦、12^◦ 和 78^◦、22.5^◦以及 67.5^◦、30^◦和 60^◦、36^◦和 54^◦、6^◦和 84^◦以及 45^◦ 所求係數 a₁ ∼ a3 的分布情況。當天頂角 (Zenith angle)θ 為 90^◦ 時， 可得到二維情況下的係數；當天頂角 (Zenith angle)θ 為 90^◦ 和方位 角 (Azimuth angle)ϕ 為 0^◦ 或 90^◦時，可得到一維情況下的係數。 . . 33

Table 4.3 在三維情況下，本文所提出的 FDTD 數值方法與 Yee 方法的穩定性範圍比較。 . . . 33

Table 5.1 比較本文所使用的 PRK-DRP FDTD 與 ADI-FDTD [1] 在 h = 0.01 情況下計算電場 E 與磁場 H 的最大誤差和時間收斂斜率。 . . . 37

Table 5.2 本文所使用的 PRK-DRP FDTD 當 t = 1，並選取 CFL number =0.2 時，計算電場 E 與磁場 H 的最大誤差和空間收斂斜率。 . . . . 37

Table 6.1 Apple iPhone4-like phone 所參考的材料性質. . . 49

Table 6.2 316L 不鏽鋼的電磁材料特性.[2, 3] . . . 50

Table 6.3 本研究所使用的人體器官材料特性。 . . . 53

Table 7.1 Apple iPhone4 (model A1332) 產品的真實量測值 [4]。 . . . 57

Table 7.2 頭部各器官組織其電場強度、能量比、該器官組織內部最大電場強度以及該器官組織的於頭部之體積比。 . . . 58

Table 7.3 在 5000∆t 的計算時間內，不同器官組織的 SAR_peak之模擬結果。 . . . 65

Table 7.4 不同器官組織的 SARpeak、SAR1g 和 SAR10g 模擬結果。 . . . 65

Table 7.5 移去顱骨後，頭部各器官組織其電場強度、能量比、該器官組織內部最大電場強度以及該器官組織的於頭部之體積比。 . . . 67

(15)

符號說明

符號對照表

E 電場強度 (Electric field intensity) V/m D 電通密度 (Electric flux density) C/m² H 磁場強度 (Magnetic field intensity) A/m B 磁通密度 (Magnetic flux density) W/m² J 電流密度 (Electric current density) A/m²

σ 電導率 (Conductivity) S/m

ϵ₀ 自由空間中的介電常數 (Permittivity of free space) F/m ϵ_r 相對介電常數 (Relative permittivity) － µ0 自由空間中的磁導係數 (Permeability of free space) H/m µ_r 相對磁導係數 (Relative permeability) － k₀ 自由空間中的波數 (Wavenumber in free space) 1/m β 行進波數 (Propagation constant) 1/m ω 角頻率 (Angular frequency) rad/s

λ 波長 (Wavelength) m

c 光速 (Speed of light in free space) m/s

(16)

第一章序論

1.1 前言

自 1873 年馬克斯威爾 (Maxwell) 建立電磁場基本方程組以來，電磁波理論與應用已經發展了一段相當長的時間，並已深入各個領域，應用範圍包括無線電波傳播、通訊、光纖、射頻電路、天線、微波、電動力學、地下電磁探測、電磁兼容、及電磁隱身應用等等。

實際上，電磁波在真實空間中的傳播過程，是具相當複雜之物理傳播特性，

例如複雜結構的天線輻射、各種目標物的散射，以及波導和微帶天線中的傳播、

城市複雜的地形、海面、亦對電磁波傳遞深具影響。

深入探討電磁波的物理現象，以及發展出各項針對電磁波的應用，對人類的科技具有十分深遠之影響性。因為環境或散射物體的幾何複雜性，因此通常必須經由數值模擬，才能得到於實際測量中所無法得知的電磁波傳遞現象，因此，結合科學計算、計算機硬體技術、和電磁理論等相關的知識的計算電磁學 (Computational Electromagnetics)，儼然是電磁領域中一個相當重要的研究領域。

在科學計算領域 (Scientific Computing)，數值分析是近幾年受益於計算機技術的快速發展，已演變成一個相當有效的科學研究手段。

對於計算電磁領域中的應用，係以電磁理論作為基礎，並透過高性能計算技術為手段，利用所開發之數值方法，以獲得準確度與計算效率和穩定性，最終在進行實際測量與量產前，可獲得符合成本與精確度的模擬結果。

計算電磁學可以根據計算域的不同，可分為頻域法與時域法兩大類。頻域法主要有頻域有限差分 (Frequency Dependent Finite Difference, FDFD) 方法、矩量法 (Moment Method, MoM) 和有限元素法 (Finite Element Method, FEM) 等；而時域法主要有時域有限差分 (Finite Difference Time Domain, FDTD) 方法、時域有限體積 (Finite Volume Time Domain, FVTD) 法、和時域積分方程 (Time Domain Integral Equation, TDIE) 法。其中，時域有限差分法 (FDTD) 是一種對馬克斯威爾方程組進行直接求解的數值方法，它既可以分析時諧場，又能計算寬頻帶的瞬態場，同時也能處理複雜外型目標的電磁散射或者非均勻、各向異性色散介質的電磁問題，

對電磁散射以及近場總場的分析具有其優勢，目前是處理各類介質電磁散射與近

(17)

場問題的主要方法。

自從 1973 年，全世界第一支手機-Dyda TAC 8000X [黑金剛] 的誕生，經歷無線蜂巢行動通訊的發明，時至今日在全民手機時代，實際上，也代表電磁波早已遍布在我們生活周遭環境中。因此，開始有人關心，手機所產生的高功率輻射電磁能，對健康上的影響，已是需要進行研究的主題。然而，一些生物實驗已經指出，手機輻射或許是導致實驗樣本致病的一個原因，因此，了解手機設計、材質選擇、手機使用方式與距離，便成為一個值得探討的課題。

(18)

1.2 文獻回顧

1966 年 Yee[5] 首次提出了時域有限差分 (Finite Difference Time Domain, FDTD) 方法，並將其應用於電磁波時域之計算，爾後這一研究取得非常迅速的發展與廣泛的應用。時域有限差分方法係直接由馬克斯威爾旋度方程出發，通過對電場 E 和磁場 H 的各分量在空間和時間上採取交替的離散方式，即每一個 E (或 H) 場分量周圍有四個 H (或 E) 場分量的環繞，此一電場與磁場交錯的方式被稱為 (Yee’s cell)。使用這種離散方式，將含時間變量的馬克斯威爾旋度方程轉化為一組時域之離散差分方程，然後沿時間軸逐步推進地求解空間電磁場，此一方法清楚地顯示了電磁波傳遞的過程，亦便於分析與設計。在 1975 年 Taflove 等用 FDTD 計算非均勻介質在正弦波入射時的穩態時諧場電磁散射，並提出了數值穩定的分析方法。Taflove 于 1988 年詳細地探討了 Yee 算則的數值色散誤差。

當時由 Yee 所提出的 Yee’s cell 僅用來模擬電磁波在空間上的傳遞；為了在有限計算區域內模擬無限空間中的電磁問題，必須在計算區域的截斷邊界上設置吸收邊界條件，自從 FDTD 方法誕出以來，對於吸收邊界的研究一直是個熱門且重要的研究課題。吸收邊界方法最初于 1969 年由 Taylor 所提出與插值方法有關的吸收邊界條件，後來並廣泛地被採用的 Mur 吸收邊界 (Absorbing Boundary Condition，ABC) [6]，及至於上世紀 90 年代所發展的完美匹配層 (Perfectly Matched Layer，PML) [7, 8]，其吸收的效果越來越好。

近年來，人們以 FDTD 方法應用於色散介質，並已進行了大量的研究。目前，FDTD 適用於截斷色散介質的最常用及有效的吸收邊界條件稱為卷積完全匹配層 (Convolutional Perfect Match Layer, CPML) [9]。

於射頻方面的應用，如完美電導體 (Perfect Electric Conductor) 此類的不連續性材質間，材質所造成之不連續性 [10][11][12][13]，在時域有限差分模擬上，已經獲得突破。此後，由 David M.Sheen 和 Sami M. Ali 等人 [14] 利用 Yee 所開發在交錯網格系統下的時域有限差分方法，針對微帶天線技術進行模擬，透過在 FDTD 網格計算域中設置一個極薄的介質基板上，底面以金屬薄層接地，其上則刻出一特定形狀之金屬貼片，並利用簡化微帶線，對金屬貼片進行饋電，其模擬結果和實驗值曲線相當的接近。

同時，在傳輸線方面的應用 [15]，其線內與外部環境的不連續性也在計算電磁學的領域獲得突破，並且在模擬時，可作為天線的瞬態源項。

上述的射頻應用，可知電磁輻射已被大量的應用，，利用時域有限差分法進行模擬，可以節省研究人力與實驗經費。

隨著這些電磁輻射的科技開發，日常生活中開始充滿各種來源的低頻電磁輻射，電磁輻射造成人體潛在的危險，例如手機和基地台天線所發出的電磁場，已

(19)

被大量的關注。目前，比吸收率 (Specific Absorption Rate, SAR) 是國際上公認的有關電磁輻射計量標準。日常環境中，人體輻射曝曬所受到頻率範圍 1 MMhz 至 1 Ghz 之間，並使用了 30 Mhz 時的參考輻射功率 0.4 W，針對身體各組織進行了一系列模擬觀察各部位的 SAR 分布的模擬研究 [16]。然而，由於手機在使用上，

會直接靠近人體頭部進行通話，針對不同手機通信之工作頻段 (通常介於 900 Mhz 至 1800 Mhz) 與不同手機之發射功率，經常通過實驗 [17] 與模擬 [18] 兩種方式，

以觀察或預測出手機對於人體頭部之電磁比吸收率的分布。因此，SAR 作為輻射的安全標準，在學術與業界亦被大量的用來作為評估與設計手機的安全規範。

(20)

1.3 研究動機

就物理意義而言，電磁波傳遞行為之時間與空間有一定的相關性，即角頻率 (Angular frequency) 與波數 (Wave number) 的關係；傳統的 FDTD 方法一般採用蛙跳格式 (即空間與時間離散階為二階準確之中央差分格式)，嚴格看來，其所使用的數值方法與物理意義上顯然並無實質的連結，為了彌補這些方法在物理意義上的不足，本文嘗試以電磁波之色散關係為方法開發的核心，將在時間與空間離散上取得直接的連結，並利用具守恆性質的時間離散，在空間上搭配具最佳色散關係保持之非交錯網格的空間離散，以期電磁波傳遞行為在時間與空間有著更正確的緊密結合，旨在確保每一個模擬的時間過程中，皆能得到具物理意義的最佳色散關係解。

本文所使用之數值方法，具備時間及空間相關聯之物理特性。並且利用模擬預測電磁輻射的實際問題，得以符合物理意義地去佐證本文所使用之方法的實用性，透過非交錯網格系統進行程式的實作，這樣的結構有助於高速平行計算的發展，進而使得原本交錯網格系統下，因平行計算程式的不易編寫及較複雜物理量儲存特性所造成之昂貴計算量，透過本文在非交錯網格架構下所開發之數值方法，

有助於未來在計算電磁領域上能有更多元的發展性。

(21)

1.4 研究目標

目前時域有限差分法大多沿用 Yee [5] 所提出之交錯式網格 (Yee’s cell) 來求解馬克斯威爾方程組，並能適當地模擬電磁波的行為。但大部的時域有限差分方法其時間離散與空間離散卻是彼此獨立，即僅僅針對數值求解安培定律即法拉第定律。然而，在電磁波傳遞的物理行為上，時間與空間必須是具相依的關係；換言之，與波數與角頻率之間的關係 (亦稱為色散關係，Dispersion relation) 在模擬的過程中是需要被考慮的。

我們以理想的馬克斯威爾方程為核心的概念，使用具辛結構 (symplecticity) 的時間離散方法，以確保滿足方程中所當具有的守恆性質，並以能量守恆的觀點檢驗時間格式是否滿足守恆性質。往外延伸至于色散介質中傳播的電磁波方程，進而處理色散介質所造成的極化電流項 (Polarization current density)。對於波傳遞之時間與空間的相依關係，我們以色散關係保持 (Dispersion relation preserving) 式作為核心，並且引入平面波的概念於方程式中，如此，可取得時間與空間合乎物理的正確聯結，發展出一具最佳色散關係保持性之馬克斯威爾方程的離散方法。

為證明此一數值方法的實用性，本論文將進一步討論手機輻射對人體頭部的吸收情形，並且觀察人體組織器官與手機間的距離、人體組織以及手機的材料性質，是否將影響組織器官的電磁輻射吸收，考慮一般手機通話時的使用方式，並透過結果與實際測量值進行對照，通過本論文探討人體電磁特定比吸收率的模擬結果，並加以思考手機輻射是否會對人體造成不良的影響。

(22)

1.5 論文大綱

本文本諸理想的馬克斯威爾方程與其色散關係方程式的核心概念，發展出一具最佳數值色散關係式之非交錯格式的空間離散方法。首先，於第一章中介紹計算電磁波的相關背景與發展之文獻回顧，以及本研究的動機與目標。本文的第二章介紹三維的馬克斯威爾方程組，以及在時域中，開發時間與空間具物理相關特徵的顯式有限差分離散格式。第三章中，為了有效的執行有限差分法，此章節簡要的描述了非交錯網格中所使用的數值模型。然而，所提出的 FDTD (時域有限差分法) 數值模型滿足了具辛結構以及數值色散關係方程的性質。此外，也針對被卷積完全匹配層 (CPML) [9] 所包圍的三維計算空間中進行了簡要的描述。第四章中，詳細敘述了吾人所開發之數值方法，將時間與空間之離散的方法進行推導。

在第五章，透過實解驗證，並且與 ADI-FDTD [1] 交替方向隱式的時域有限差分方法進行時間及空間上收斂斜率之比較，同時驗證了本文所開發之方法準確性，以及程式的可信度。第六章中，詳細介紹了正在研究的問題，其中包括貼置在模擬人頭部耳上的 Apple iPhone4-like phone 的模型建模，以及建模方法。並且回顧了大腦組織的主要組成。在第七章中，詳細討論了透過求解馬克斯威爾方程組模擬手機對人腦之輻射問題的模擬結果。最後，我們將在第八章中給出結論以及未來的展望。

(23)

第二章馬克斯威爾方程組

2.1 法拉第/安培/高斯方程組及其推導

馬克斯威爾方程組係由法拉第定律、安培定律及高斯定律所組成。在時域中，

分別可以表示為如下的安培定律及法拉第定律表示式，

∇ × ⃗E = −∂ ⃗B

∂t − ⃗Jm, (2.1)

∇ × ⃗H = ∂ ⃗D

∂t + ⃗J , (2.2)

以及高斯定律所組成的限制方程式

∇ · ⃗B = 0, (2.3)

∇ · ⃗D = ρ. (2.4)

上述方程組中， ⃗H 為磁場強度 (Magnetic field intensity)、 ⃗E 為電場強度 (Electrical field intensity)、 ⃗B、磁位移 (Magnetic displacement) 和 ⃗D 電位移 (Electrical displace- ment)。式 (2.2) 與 (2.1) 中之 ⃗J 和 ⃗J_m分別代表電流密度 (Electric current density) 和 磁流密度 (Equivalent magnetic current density)。于 (2.4) 之 ρ 代表電荷密度 (Electri- cal charge density)。在各向同性 (Isotropic) 的線性介質中需滿足如下的本構關係式

D = ϵ ⃗⃗ E, ⃗B = µ ⃗H, ⃗J = σ ⃗E, ⃗J_m = σ_mH.⃗ (2.5)

上式中，ϵ 和 µ 分別表示介電係數 (Electric permittivity) 和磁導係數 (Magnetic permeability)，σ 和 σ_m 分別為電導率 (Electric conductivity) 和磁導率 (Equivalent magnetic loss)，兩者分別代表介質的損耗和磁損耗。本文將所推導的馬克斯威爾方程之離散方法係居于顯式的格式，用來模擬電磁波於簡單介質 (Simple medium) 中的傳播行為。主要探討之材料性質係數中的導磁係數 (Magnetic permeability，µ) 以及介電係數 (Electric conductivity，ϵ) 皆視為常數，此二材料性質係數主要是在

(24)

描述材料之介電性質以及導磁特性且皆為正其中，ϵ = ϵ0· ϵr 而 µ = µ0 · µr；兩材 料性質係數決定了波的傳播速度 c (c≡ 1/√ϵµ(≈ 3.0 × 10⁸m/s))，然而，針對色散 介質而言，介質的介電常數，以及波在色散介質內之傳播速度，將與頻率相關，

我們將透過頻域下之馬克斯威爾方程組進行離散推導，並透過使用吾人針對色散介質所開發的高階顯式輔助差分方法 ADE-FDTD [19] 進行額外處理。

2.2 法拉第/安培方程組之數學特性

依上節所述之假設情況下，在無源且均勻各向同性簡單介質下之馬克斯威爾方程組可以改寫為

1

µ(∇ × ⃗E) = −∂ ⃗H

∂t , (2.6)

1

ϵ(∇ × ⃗H) = ∂ ⃗E

∂t, (2.7)

∇ · ⃗B = 0, (2.8)

∇ · ⃗D = 0. (2.9)

式 (2.8-2.9) 以及式 (2.5) 並非獨立於馬克斯威爾方程組 [20] 之外；即式 (2.8-2.9) 可視為馬克斯威爾方程組之限制條件。將式 (2.6) 及 (2.7) 兩端分別取散度算子 (divergence operator)，可得知出高斯定律係隱含於法拉第定律與安培定律之中。

在非交錯網格系統中，吾人可直接地求解雙曲線之系統方程 (2.6-2.7) 得到瞬時之電場與磁場。對於式 (2.8-2.9) 之無散度條件 (divergence-free constraint conditions) 之限制，吾人將於程式驗證中進行分析。

馬克斯威爾方程組可以利用漢彌爾頓動力系統 (Hamiltonian dynamical system) 表示之，即

∂

∂t ( H⃗

E⃗ )

=

( 0 −I I 0

) ( δH / δ ⃗H δH / δ ⃗E

)

, (2.10)

其中，漢彌爾頓方程 (Hamiltonian function，H) 可以下式表示之 [21]

H( ⃗H, ⃗E) =

∫

Ω

1 2(1

ϵ

H⃗ · ∇ × ⃗H + 1 µ

E⃗ · ∇ × ⃗E) dΩ. (2.11)

根據 Kole 等人之研究 [22]，式 (2.6-2.7) 可寫成 _∂t^∂ ψ(t) = G ⃗⃗ ψ(t)，其中 ψ(t)⃗ ≡ (⃗m(t), ⃗n(t))^T (

= µ^1/2H(t), ϵ⃗ ^1/2E(t)⃗ )T

。上述之 G 矩陣，它可以如下之反

(25)

對稱矩陣 (skew-symmetric matrix) 表式之 [22]

G =

( 0 −µ^−1/2∇ × (ϵ^−1/2) ϵ^−1/2∇ × (µ^−1/2) 0

)

. (2.12)

馬克斯威爾方程組可由 ∂ ⃗ψ(t) / ∂t = G ⃗ψ(t) 表示之，如此吾人可得該方程之通解 ψ(t) = exp⃗

(

t G ⃗ψ (t = 0)

)，其中 exp(t G) 顯示了馬克斯威爾方程組之解係隨時間 演進呈指數的變化。將向量 ⃗ψ 特徵化為∫

Ωϵ ⃗E· ⃗E + µ ⃗H· ⃗HdΩ，它與電磁場之能 量密度 (energy density，w(t)) 存在直接的理論關聯，即

w(t) =

∫

Ω

ϵ ⃗E· ⃗E + µ ⃗H· ⃗H dΩ. (2.13) 式 (2.13) 之能量密度不會隨時間而改變。對於式 (2.11) 及 (2.13) 隨時間皆不變之性質，吾人可將其作為程式驗證之用。

(26)

2.3 色散介質

常見的色散介質模型包括德拜模型 (Debye model)、洛侖茲模型 (Lorentz model) 和德魯模型 (Drude)，它們所對應的單極介電係數表達式如下:

(1) Debye 介質

ϵr(ω) = ϵ_∞+ ϵ_s− ϵ∞

1 + Iωτ + σ

Iωϵ₀ ≡ ϵ∞+ ∆ϵ

1 + Iωτ + σ

Iωϵ₀ (2.14) (2) Lorentz 介質

ϵ_r(ω) = ϵ_∞+ ∆ω²_p

ϵ²_p+ 2Iωδ_p− ω² + σ

Iωϵ₀ (2.15)

(3) Drude 介質

ϵ_r(ω) = ϵ_∞+ ω_p²

Iωγ_p− ω² + σ

Iωϵ₀ (2.16) 其中，∆ϵ = ϵs− ϵ∞。ϵs為靜態或零頻時的相對介電係數、ϵ_∞為無窮大頻率時的 相對介電係數、δ_p 為阻尼係數、ω_p 為媒質的諧振頻率 (Resonant frequency)、τ 為 鬆弛時間 (Relaxation time)、且 γ_p = 1/τ 為鬆弛時間的倒數。

考慮式 (2.14-2.16)，式 (2.2) 安培定律旋度方程均可表示為

∇ × ⃗H = ϵ₀ϵ_∞∂ ⃗E

∂t + σ ⃗E + ⃗J_d, (2.17) 其中 ⃗J_d項為介質的色散極化所引起的極化電流項 (Polarization current density)，對於不同類型的色散介質，它們分別滿足如下的不同輔助偏微分方程

J⃗_d+ τ∂ ⃗J_d

∂t = ϵ₀∆ϵ∂ ⃗E

∂t , (Debye) (2.18) ω_p²J⃗_d+ 2δ_p∂ ⃗J_d

∂t +∂²J⃗_d

∂t² = ϵ₀∆ϵω_p²∂ ⃗E

∂t , (Lorentz) (2.19) γ_p∂ ⃗J_d

∂t +∂²J⃗_d

∂t² = ϵ₀ω_p²∂ ⃗E

∂t. (Drude) (2.20) 在色散介質的 FDTD 計算方面，如前章所述，有遞推卷積 (RC) 法、輔助方程 (ADE) 法、Z 變換法、分段線性遞推卷積 (PLRC) 等方法。本文在色散介質的研究將採用二階輔助方程 (ADE) 法 [19]，其中，在實際題目應用上，將引用輔助方程法，在模擬計算時，針對 Drude model 色散介質進行特別的處理。

(27)

2.4 卷積完美匹配吸收層

在自然界中，大部分電磁波問題需要在無限域下模擬方可捕捉其真實的物理本質，礙於電腦的記憶體之限制，理論上吾人不可能於無限大空間內進行計算，

及儲存大量的數據。當初 Yee 所提出的 FDTD 方法時，計算區域的人為截斷所導致的數值回波問題即已顯現，並導致了 FDTD 在最初的十多年發展較為緩慢，直到 1981 年荷蘭科學家 G.Mur 提出了 Mur 吸收邊界條件 [6], 1994 年 Berenger 提出了高吸收效能的完美匹配層 (Perfectly Matched Layer, PML) 邊界條件 [7]，才使得 FDTD 方法煥發出新的生命力。目前，已發展出數種 PML 方法，即 Gedney 所提出的軸向完美匹配層 (uniaxial PML, UPML) [23]，Abarbanel 所提出的各向異性完美匹配層 (anisotropic PML, APML) [24] 等方法。本文所使用的吸收邊界稱為卷積完美匹配層 (convolution PML, CPML) [9]，其優點在於它對於消逝波 (evanescent wave) 具有很好的吸收能力，此外它尚能大幅度減少記憶體的占用。除了吸收效果較好之外，亦具普適性，CPML 吸收邊界完全獨立於 FDTD 計算域內的介質，

可以不做任何修改地應用到各向同性、各向異性、吸收、色散或是非線性介質的計算中。因此，CMPL 相對于 UPML 吸收邊界而言，更適用於非自由空間中目標 散射特性的研究。在此，我們給出了其張量係數 (tensor coefficient)s_w 的形式 [25]

s_w = κ_w+ σ_w

a_w+ jωϵ₀ ; w = x, y, z (2.21) 式 (2.21) 所描述的是電參數沿某方向的一個空間函數。首先，引入擴展座標 (stretched-coordinate) 轉換關係式

˜ w→

∫ w 0

s_w(w^′) dw^′ ; w = x, y, z , (2.22)

透過傅立葉轉換關係，取ℑ⁻¹ 為傅立葉逆變換算子，定義

¯

s_w =ℑ⁻¹

( 1 κ_w+ σ_w

)

= δ(t)

κ_w − σ_w ϵ₀κ²_wexp

[

− ( σ_w

ϵ₀κ_w + a_w ϵ_w

) t

]

u(t)≡ δt

κ_w + ζ_w(t) (2.23)

(28)

(其中 u(t) 和 δ(t) 分別是單位步階函數 (Unit step function) 和脈衝函數 (Delta function))，可得馬克斯威爾方程組的時域形式:

∂D_x

∂t = ( 1

k_x

∂H_z

∂y − 1 k_z

∂H_y

∂z + ζ_y ∗∂H_z

∂y − ζz∗ ∂H_y

∂z )

,

∂D_y

∂t = (1

k_z

∂H_x

∂z − 1 k_x

∂H_z

∂x + ζ_z∗ ∂H_x

∂z − ζx∗ ∂H_z

∂x )

,

∂D_z

∂t = ( 1

k_x

∂H_y

∂x − 1 k_y

∂H_x

∂y + ζ_x∗∂H_y

∂x − ζy∗ ∂H_x

∂y )

, (2.24)

∂B_x

∂t =− ( 1

k_x

∂E_z

∂y − 1 k_z

∂E_y

∂z + ζ_y∗ ∂E_z

∂y − ζz∗ ∂E_y

∂z )

,

∂B_y

∂t =− (1

k_z

∂E_x

∂z − 1 k_x

∂E_z

∂x + ζ_z∗ ∂E_x

∂z − ζx∗ ∂E_z

∂x )

,

∂B_z

∂t =− ( 1

k_x

∂E_y

∂x − 1 k_y

∂E_x

∂y + ζ_x∗ ∂E_y

∂x − ζy ∗∂E_x

∂y )

. (2.25)

在 FDTD 的模擬中，CPML 吸收邊界內電磁場的迭代計算可透過 (2.24) 和 (2.25) 的離散式來實現。計算主要的困難是離散卷積 ζw∗(∂v/∂w) 的執行，因為直接進行卷積項的計算，將會消耗大量的記憶體和計算時間。透過 Luebbers 和 Hunsberger [26] 所採用的遞推卷積 (RC) 技術，可以有效的解決此一問題。

定義式 (2.23) 中脈衝 ζ(t) 的離散響應 Z_w(m) 為

Z_w(m) =

∫ (m+1)∆t m∆t

ζ_w(τ ) dτ ≈ − σw

ϵ₀κ²_w

∫ (m+1)∆t m∆t

exp [

− ( σw

ϵ₀κ_w + aw

ϵ₀ )

τ ]

dτ (2.26)

= c_wexp [

− ( σ_w

ϵ₀κ_w +a_w ϵ₀

) m∆t

]

(2.27) 其中

c_w = σ_w σwκw+ κ²_waw

[ exp

(

− ( σ_w

ϵ0κw

+ a_w ϵ0

)

∆t )

− 1 ]

(2.28)

式 (2.24) 中 ζw 和 ∂Hv/∂w 卷積離散形式可以近似成

ψ_w,v(n) = ζ_w(t)∗ ∂

∂wH_v(t)|t=n∆t ≈

n−1

∑

m=0

Z_w(m) ∂

∂wH_v(n− m) (2.29) 由式 (2.29) 可看出，計算 n∆t 時刻的離散卷積 ψ 需經 n 次的乘法和加法，這在實

(29)

際的計算中將嚴重地影響計算速度。Luebbers 和 Hunsberger 將式 (2.29) 寫為以下的遞推形式:

ψ_w,v(n) = b_wψ_w,v(n− 1) + cw

∂

∂wH_v(n) (2.30) 其中，c_w 如式 (2.28)，b_w為

b_w = exp [

− ( σ_w

ϵ₀κ_w + a_w ϵ₀

)

∆t ]

. (2.31)

如此一來，式 (2.30) 中的 ψ_w,v 可以透過簡單的時間循環迭代來實現應用 CPML 吸 收邊界。最後，我們將 ζ_w 和 ∂H_v/∂w 的卷積離散由 ψ_w,v(n) 取而代之，並延伸至 具不同材料色散介質特性之應用，可得到如下具 CPML 吸收邊界層之馬克斯威爾方程組的時域形式:

∂E_x

∂t = 1 ϵ0ϵr

( 1 ky

∂H_z

∂y − 1 kz

∂H_y

∂z − Jd,x+ ψ_E_x,y − ψEx,z

) ,

∂E_y

∂t = 1 ϵ₀ϵ_r

( 1 k_z

∂H_x

∂z − 1 k_x

∂H_z

∂x − Jd,y+ ψEy,z − ψEy,x

) ,

∂E_z

∂t = 1 ϵ₀ϵ_r

( 1 k_x

∂H_y

∂x − 1 k_y

∂H_x

∂y − Jd,z+ ψEz,x − ψEz,y

)

, (2.32)

∂H_x

∂t =− 1 µ0µr

( 1 ky

∂E_z

∂y − 1 kz

∂E_y

∂z + ψ_H_x,y − ψHx,z

) ,

∂H_y

∂t =− 1 µ₀µ_r

( 1 k_z

∂E_x

∂z − 1 k_x

∂E_z

∂x + ψ_H_y,z− ψHy,x

) ,

∂H_z

∂t =− 1 µ₀µ_r

( 1 k_x

∂E_y

∂x − 1 k_y

∂E_x

∂y + ψHz,x − ψHz,y

)

. (2.33)

(30)

其中

ψⁿ_E_x,y = by · ψ_Eⁿ⁻¹_x,y + cy· ∂H_zⁿ

∂y , ψ_Eⁿ_x,z = bz· ψⁿ_E⁻¹_x,z + cz· ∂H_yⁿ

∂z , ψ_Eⁿ_y,x = b_x· ψ_Eⁿ⁻¹_y,x + c_x·∂H_zⁿ

∂x , ψ_Eⁿ_y,z = b_z · ψ_Eⁿ⁻¹_y,z + c_z· ∂H_xⁿ

∂z , ψ_Eⁿ_z,x = b_x· ψⁿE⁻¹z,x + c_x· ∂H_yⁿ

∂x , ψ_Eⁿ_z,y = b_y· ψⁿE⁻¹z,y + c_y· ∂H_xⁿ

∂y , ψⁿ⁺

1 2

Hx,y = b_y· ψ_Hⁿ⁻_x,y¹² + c_y· ∂Eⁿ⁺

1

z 2

∂y , ψⁿ⁺

1 2

Hx,z = b_z· ψⁿ_H⁻_x,z¹² + c_z· ∂Eⁿ⁺

1

y 2

∂z , ψⁿ⁺

1 2

Hy,x = b_x· ψHⁿ⁻y,x¹² + c_x· ∂Eⁿ⁺

1

z 2

∂x , ψⁿ⁺

1 2

Hy,z = b_z· ψⁿH⁻y,z¹² + c_z· ∂Eⁿ⁺

1

x 2

∂z , ψⁿ⁺

1 2

Hz,x = b_x· ψHⁿ⁻z,x¹² + c_x·∂Eⁿ⁺

1

y 2

∂x , ψⁿ⁺

1 2

Hz,y = b_y· ψHⁿ⁻z,y¹² + c_y· ∂Eⁿ⁺

1

x 2

∂y .

(2.34) 式 (2.34) 中所述之 b_w 及 c_w(w = x or y)，可表示成:

b_w = e

(−(

σw ϵ0kw+^aw

ϵ0

)

∆t)

; c_w = σ_w σwkw + k_w²aw

( e

(−(

σw ϵ0kw+^aw

ϵ0

)

∆t)

− 1 )

. (2.35)

(31)

其中

σw = σw,max

(d− w d

)m

, a_w = a_w,max

(d− w d

)ma

; 0≤ w ≤ d, k_w = 1 + (k_w,max− 1) ·

(d− w d

)m

. (2.36)

式 (2.36) 中的 d 代表 CPML 厚度，而 σ_max 之大小可以下式計算之

σ_w,max = 0.8(m + 1)

η₀∆√ϵ_rµ_r, (2.37) 其中 η₀ =√

µ₀/ϵ₀為真空中的波阻抗。

圖2.1為 CPML 之示意圖，利用圖2.2說明了以理想馬克斯威爾方程為核心概念，並延伸至具不同材料色散介質特性的應用。CPML 的概念圖如下。

(32)

Fig 2.1: CPML 吸收邊界之示意圖。

Fig 2.2: 求解電磁波傳遞問題的方程及計算空間之示意圖。

(33)

第三章數值方法

求解電磁波方程的數值方法主要可分為時域和頻域兩大類。頻域法包括有矩量法 (Moment method, MoM)、有限元素法 (Finite element method, FEM) 等；時域法包括有時域有限差分法 (Finite difference time domain, FDTD)、時域積分方程法 (Time domain integral equation, TDIE) 等。本論文根據差分理論 (difference theory)，

利用泰勒級數展開式 (Taylor series) 將微分方程轉換成差分方程。計算區域將劃分成有限個網格點，然後於網格點上對統御方程進行離散與計算。

3.1 非交錯網格系統下之 FDTD 離散方法

在卡氏座標下，于理想介質下考慮無源之馬克斯威爾方程式 (2.6-2.7)，可以將安培定律及法拉第定律之旋度方程展開成

∂E_x

∂t = 1 ϵ

(∂H_z

∂y − ∂H_y

∂z )

∂E_y

∂t = 1 ϵ

(∂H_x

∂z −∂H_z

∂x )

∂E_z

∂t = 1 ϵ

(∂H_y

∂x − ∂H_x

∂y )

∂H_x

∂t = 1 µ

(∂E_y

∂z − ∂E_z

∂y )

(3.1)

∂H_y

∂t = 1 µ

(∂E_z

∂x − ∂E_x

∂z )

∂H_z

∂t = 1 µ

(∂E_x

∂y − ∂E_y

∂x )

以時域有限差分方法求解三維馬克斯威爾方程組，探討手機于頭部之電磁波之比吸收率之分布

國立臺灣大學工學院工程科學及海洋工程學系 碩士論文

Department of Engineering Science and Ocean Engineering College of Engineering

National Taiwan University Master Thesis

以時域有限差分方法求解三維馬克斯威爾方程組，探 討手機于頭部之電磁波之比吸收率之分布

Calculation of Mobile Phone Induced Specific Absorption Rate Prediction in Phantom Head by a 3D FDTD code for

Maxwell’s equations

林緯皓 Wei-Hao Lin

指導教授：許文翰博士

Advisor: Tony Wen-Hann Sheu, Ph.D.

中華民國 107 年 7 月

July, 2018

誌謝

摘要

Abstract

目錄

圖目錄

表目錄

符號說明

第 一 章 序論

1.1 前言

1.2 文獻回顧

1.3 研究動機

1.4 研究目標

1.5 論文大綱

第 二 章 馬克斯威爾方程組

2.1 法拉第/安培/高斯方程組及其推導

2.2 法拉第/安培方程組之數學特性

2.3 色散介質

2.4 卷積完美匹配吸收層

第 三 章 數值方法

3.1 非交錯網格系統下之 FDTD 離散方法

國立臺灣大學工學院工程科學及海洋工程學系碩士論文

以時域有限差分方法求解三維馬克斯威爾方程組，探討手機于頭部之電磁波之比吸收率之分布

誌謝　　　

第一章序論

第二章馬克斯威爾方程組

第三章數值方法