中 華 大 學

中 華 大 學 碩 士 論 文

題目：數值分析震波對壁面熱傳效應之影響

系 所 別：機械與航太工程研究所碩士班 學號姓名：M09408022 陳 渝 瑋 指導教授：鄭 藏 勝 博 士

中華民國九十六年 七 月

中文摘要

本研究主要是以數值模擬對震波通過背向式階梯之流場進行分 析，探討在不同階梯高度，不同壁面溫度及不同震波馬赫數之下，震 波對壁面熱傳效應之影響。統御方程式是使用納維爾-史脫克方程式 (Navier-Stokes equation)，並於數值方法上，對黏滯項及對流項之空間 差 分 分 別 採 用 四 階 精 度 的 中 央 差 分 法 及 五 階 精 度 的 Weighted Essentially Non-Oscillatory Scheme(WENO)，而在時間積分上採用顯 式的四階 Runge-Kutta scheme，進行二維平面流場之分析。

計算結果發現對於流場之結構方面，在背向式階梯的角落部份由 於剪切層的效應會產生渦卷，當反射震波與渦卷產生交互做用會造成 複雜的流場結構。於階梯底部壁面加入溫度變化雖然對震波結構影響 較小，但影響底部壁面紐塞爾數。在相同馬赫數、階梯高度下，壁面

絕熱與壁面加至熱 時，其流場結構並無太大之改變。但對於階

梯高度之改變，會對紐塞爾數造成影響，當階梯高度越高，其紐塞爾 數會隨之增加。對瑞理數的影響參數是震波與階梯底部壁面之間的溫 度差。

K 3000

Abstract

Numerical simulations of shock wave passing through a two-dimensional backward-facing step channel are performed to investigate the effects of shock wave Mach number, channel step height, and channel bottom wall temperature variation on the flowfield structure and heat transfer. The compressible Navier-Stokes equations governing the flowfield in planar coordinates are discretized using the fifth order accuracy Weighted Essentially Non-Oscillatory (WENO) and the fourth order accuracy central difference schemes for the convective and diffusive terms, respectively. The governing equations are integrated with the fourth order accuracy Runge-Kutta scheme.

Computed results show that the interactions of the step corner induced vortex with the reflected shock waves from the walls form complicated flow structures. The effects of shock wave Mach number, channel step height, and channel bottom wall temperature variation on the heat transfer are characterized by the nondimensional parameter, Nusselt number. The variation of wall temperature is represented by Rayleigh number. For the fixed shock wave Mach number condition, the calculated local Nusselt number increases with increasing the step height and Rayleigh number. For the fixed Rayleigh number condition, the calculated Nusselt number decreases with increasing Mach number but increases with increasing step height.

誌 謝

本論文得以完成，首先我要感謝我的指導老師鄭藏勝博士，在我 就讀研究所的兩年中給我的許多指導，不論是學業上很照顧我，連生 活上老師也給我許多的幫助，老師讓我學會了如何去做好一件事情，

使我懶散的個性有了很大的改變，因此在此跟老師說一聲謝謝老師。

另外感謝蔡永培老師與黃和順老師在百忙中抽空來擔任本論文的口 試委員，並提供寶貴的意見，使本論文能更臻完善，在此深感謝意。

在論文撰寫其間很感謝本實驗室的盧韋廷學弟與陳永軒學弟，以 及鋁合金實驗室的林家宇同學與疲勞力學實驗室的柯智瑋同學，還有 CFD 實驗室的周中祺同學，在我生活上給予我很多的協助與討論。

最後要感謝我的爸爸媽媽跟我弟弟，在我論文撰寫時遇到瓶頸時 給我的鼓勵與關懷，使本論文能順利完成。

該感謝的朋友實在太多太多了，寫再多還是掛一漏萬，謝謝曾經 協助我的朋友，謝謝中華大學機航所的所有同學，謝謝你們，謝謝。

目錄

中文摘要 --- i

英文摘要 --- ii

誌 謝 --- iii

目 錄 ---iv

圖目錄 ---vi

表目錄 ---vii

符號說明 --- viii

第一章 緒論 --- 1

1-1 前言 --- 1

1-2 研究動機 --- 2

1-3 文獻回顧 --- 2

第二章 物理問題 --- 7

2-1 物理模式--- 7

2-2 基本假設--- 7

第三章 數值方法 --- 9

3-1 統御方程式--- 9

3-2 時間積分--- 12

3-3 空間差分--- 12

3-4 起始條件與邊界條件 --- 18

3-4.1 起始條件 --- 18

3-4.2 邊界條件 --- 18

3-5 局部平均紐賽爾數及紐賽爾數--- 20

3-6 瑞理數--- 20

第四章 結果與討論 --- 21

4-1 程式驗證--- 21

4-2 格點測試 --- 22

4-3 流場結構分析 --- 23

4-4 階梯高度對紐塞爾數(Nusselt number)之影響--- 25

4-5 馬赫數對紐塞爾數之影響 --- 28

第五章 結論與未來展望 --- 30

5-1 結論 --- 30

5-2 未來展望 --- 31

參考文獻 --- 32

圖目錄

圖 2-1 背向式階梯流場模擬區域示意圖 --- 36 圖 4-1 本文之模擬結果與Jiang及Takayama[1]之實驗及模擬結果比較

圖 ， M_s=1.3 (a) Jiang 及 Takayama 之 實 驗 結 果 (b) Jiang 及 Takayama之模擬結果 (c) 本文之模擬結果，無因次化時間為 2.7--- 37 圖 4-2 本文之模擬結果與Jiang及Takayama之 實驗及模擬結果比較

圖，M_s=1.3 (a) Jiang及Takayama之實驗結果 (b) Jiang及

Takayama之模擬結果 (c) 本文之模擬結果，無因次化時間為 5.35 --- 38 圖 4-3 背向式階梯底部壁面絕熱條件下不同格點數之壓力比較圖，點

一為x=3h、y=h，點二為x=2h、y=0.5h，點三為 、 --- 39

h x=3

=0 y

圖 4-4 背向式階梯底部壁面溫度 不同格點數之壓力比較圖，點

一為 、 ，點二為

K 3000 h

x=3 y=h x=2h、y=0.5h，點三為 、 --- 40

h x=3

=0 y

圖 4-5 震波馬赫數 ，階梯高度 ，背向式階梯底部壁面為絕熱之 密度等位圖 --- 41

0 .

2 2h

圖 4-6 震波馬赫數 ，階梯高度 ，背向式階梯底部壁面溫度 之密度等位圖 --- 42

0 .

2 2h 3000K

圖 4-7 不同階梯高度不同溫度的階梯底部壁面瞬間紐塞爾數對時間 之分布比較圖，(a)階梯高度 (b)階梯高度 (c) 階梯高 度 (d)階梯高度 (e)階梯高度 --- 43

h 5 .

0 1.0h

h 0 .

2 3.0h 4.0h

圖 4-8 圖 4-8 不同階梯高度時，階梯底部壁面之局部平均紐賽爾數之 比較圖，(a)階梯高度 (b)階梯高度 (c) 階梯高度 (d)階梯高度 (e)階梯高度 --- 45

h 5 .

0 1.0h 2.0h

h 0 .

3 4.0h

圖 4-9 階梯底部壁面溫度 3000K，不同階梯高度之無因次化 近壁面三點溫度比較圖，其中

5 .

=2 Ms

K

T₀ =300 --- 47 圖 4-10 不同震波馬赫數之紐塞爾數對瑞理數比較圖 (a) 階梯高度

(b) 階梯高度 (c) 階梯高度 (d)階梯高度

(e)階梯高度 --- 48 h

5 .

0 1.0h 2.0h 3.0h

h 0 . 4

表目錄

表 4-1 絕熱壁面不同格點數壓力誤差比較 --- 35 表 4-2 加熱壁面不同格點數壓力誤差比較 --- 35

符號說明

h ：背向式階梯入口高度

k ：熱傳係數

Ms ：震波馬赫數

Nu ：紐塞爾數(Nusselt Number) P ：氣體壓力

Pr ：普蘭特數(Prandtl Number)

Ra ：瑞理數(Rayleigh Number) s ：背向式階梯階梯高度

S ：黏滯係數之 Sutherland 常數

S1 ：熱傳係數之 Sutherland 常數 T0 ：參考溫度

T1 ：震波進入背向式階梯之溫度

Tw ：背向式階梯底部壁面溫度

u ：

x

方向之速度分量

v ：

y

方向之速度分量 ρ ：流體密度

希臘符號

γ ：比熱值，γ =1.4 μ ：黏滯係數

第一章 緒論

1-1

前言

由於近年來電腦的發展日新月異，數值方法運算的效率也越來越 高，使得在數值模擬的部份有越來越多人去研究，使用的數值方法也 有許多。在現今的許多研究中，會以數值模擬之結果與實驗數據做比 較，發現數值模擬之結果與實驗結果有一定程度的近似，所以在實驗 設備不足或是實驗量測不易時，很適合使用數值模擬來預測流場物理 現象。

由於震波在數學上屬於一種不連續的問題，所以在選擇數值方法 時，需考慮到此方法對震波位置的捕捉能力，又需要考慮到其精確 度，因此格點數量不可過於稀疏，但過多的格點數會使電腦運算時間 加長，所以數值方法的效率亦是一個選擇數值方法的重要考量。此 外，當面對流場邊界問題時，亦有許多不同的考量，如震波管的入口 與出口邊界條件的設定，以及壁面是否有熱傳等相關問題，所以在出 口邊界條件設定上，必須考慮起始震波能順利傳出流場，不會有預期 外的震波反射干擾流場而影響答案的正確性。

1-2

研究動機

突 張 式 ( sudden-expansion ) 流 場 ， 亦 稱 背 向 式 階 梯 ( backward-facing step )流場，由於其簡單的幾何構造同時具有多種不 同的流場型態，不僅在學術上成為廣泛探討的主題，在工業上也被廣 泛應用在燃燒室的駐焰方面。

當震波通過一背向式階梯流場時，震波會因為流場幾何形狀的改 變而產生震波的膨脹，並於背向式階梯角落部份因為剪切層之效應造 成渦卷產生，當震波接觸壁面時會造成震波反射，與渦卷產生交互作 用的複雜流場結構。

過去對背向式階梯流場的研究主要著重在迴流區長度及迴流區 對熱傳效應影響的探討，或是震波與渦卷交互做用的研究，甚少對具 有震波與熱傳效應之流場做深入探討，因此本研究將利用數值方法，

模擬不同震波馬赫數，不同背向式階梯高度，不同壁面溫度之下，震 波對流場熱傳效應之影響。

1-3

文獻回顧

如何有效的模擬出震波通過背向式階梯，及加入熱傳對此一流場 會產生何種影響，因此整理出一些相關之文獻，概略分成以下兩個部 份：第一個部份為有關震波模擬之文獻，第二部份為有關背向式階梯

3 .

震波的模擬部份，Jiang 等人[1]以M_s =1 與M_s =1.5之震波，模 擬此震波由震波管進入膨脹管，由於體積驟變導致壓力變化引起震波 膨脹，且因為剪切層效應造成渦卷，當震波接觸壁面時會造成震波反 射，針對不同膨脹比率，探討其反射震波與渦卷交互作用和震波與剪 切層交互作用所產生的流場結構。Jiang 等人[2]以M_s =1.5及

之三維震波流場經由小體積的管路進入大體積的空腔，觀察到在空腔 成直角部份的剖面角落處有一特殊的震波反射系統，此一現象包含四 種震波：起始震波、兩個反射震波及馬赫圓盤(Mach disk)，此一震波 系統可以代表在三維流場中的馬赫反射。Jiang 等人[3]以 與 之震波，模擬此震波由震波管進入膨脹管後再拋射一固體方 塊入此流場內，觀察流場內震波反射、震波聚焦、震波與接觸表面之 交互作用等現象。Liang 等人[4]模擬震波由高壓地區進入低壓地區，

由於壓力改變引起震波膨脹產生渦卷，對渦卷各部份現象做名稱的定 義，並針對不同馬赫數探討在渦卷中心的壓力變化，以及不同馬赫數 時反射震波的位置。Teng 等人[5]以

0 .

=2 Ms

0 .

=2 Ms

4 .

=2 Ms

0 .

=2

Ms 到M_s =5.0之震波由一震 波管進入膨脹管，使用 Dispersion Controlled Dissipation (DCD) scheme 的數值方法來模擬，此篇文獻主要針對震波聚焦過程與流場結構做探 討，此震波聚焦過程包含震波繞射、聚焦以及震波反射。

背向式階梯流場模擬之部份，Hsrtnett 及 Minkowycz[6]以k−ε 紊

流模式模擬一個二維背向式階梯流場，於階梯底部以及流場上方平板 皆加入一固定溫度，觀察在不同的司徒哈數(Strouhal number)之下對 紐塞爾數(Nusselt number)及表面摩擦係數(skin friction coefficient)之 分佈有何影響。由此一文獻可瞭解到當 Strouhal number 越大時其紐 塞爾數會隨之增加，而表面摩擦係數是當 Strouhal number 為 0.25 時 呈現最大值 Strouhal number 繼續增加其表面摩擦係數則開始下降。

Sun 及 Takayama[7]使用尤拉方程式(Euler equation)與納維爾-史脫克 方程式(Navier-Stokes equation)以及k−ε紊流模式做數值模擬比較，在 格點數的選擇也選用兩種不同的格點數來分析，文中發現使用尤拉方 程式及較精細之格點會產生小渦卷，此小渦卷的產生是由於使用尤拉 方程式運算時所造成的不穩定現象，而二次渦卷的產生是因為黏滯效 應的影響，因此使用尤拉方程式無法觀察到二次渦卷的產生，在實驗 中可以觀察到二次渦卷但無法觀察到小渦卷。Abu-Hijleh[8]模擬二維 背向式階梯於底部設定一等向的多孔隙熱源加熱，並改變不同的熱源 高度以及長度，藉以觀測在不同條件下壁面紐塞爾數有何不同之改 變，由此一模擬中瞭解到當熱源高度增加時，其紐塞爾數亦隨之增 大。Chang 及 Tsay[9]模擬在背向式階梯內，自然對流之流場於靠近階 梯的垂直壁面為一高溫壁面，其相對面為低溫壁面，針對不同的水平 壁面長度，以及絕熱之垂直壁面高度，探討其紐塞爾數之變化，最後

歸納出對此流場影響之主要因素為瑞理數(Rayleigh number)、普蘭特 數(Prandtl number)、封閉區域幾何形狀大小。Avancha 及 Pletcher[10]

利用大尺度渦流法(Large eddy simulation)來模擬背向式階梯流場，並 且於背向式階梯底部加一熱傳量，探討對階梯底部通過不同熱傳量 時，不同截面位置之溫度的分佈情形。Tsay 等人[11]模擬背向式階梯 流道內，於階梯垂直壁面與階梯底部加上一固定熱傳量，並於流道上 方平板加入一阻礙物，藉此觀察阻礙物之大小以及阻礙物位置，對於 紐塞爾數之影響。瞭解到當阻礙物越靠近階梯部份時，其等溫線的分 佈較鬆散，阻礙物越離開階梯部份等溫線分佈越趨於集中，阻礙物面 積增大時在階梯底部平面之紐塞爾數會增大，而阻礙物寬度的變化對 熱傳遞影響甚小。Chen 等人[12]模擬流體通過背向式階梯流道時，於 階梯底部加入一熱傳量，觀察不同之階梯高度對於迴流區之大小有何 影響，在不同階梯高度時流場之 x 方向速度與 y 方向速度分佈有何不 同，最後再比較不同階梯高度時，對於階梯底部溫度之分佈有何種改 變，由此研究歸納出當階梯高度增大時迴流區將增大，橫向速度的峰 值將會變小，溫度最大值會增大。

由以上文獻回顧吾人瞭解到過去對於震波流場的研究多數是以 假 設 壁 面 為 絕 熱 來 模 擬 震 波 流 場 ， 並 使 用 尤 拉 方 程 式 (Euler equation)，即未將黏滯項考慮在內。而背向式階梯流場則是對於壁面

僅固定單一溫度做模擬，因此本文嘗試以納維爾-史脫克方程式 (Navier-Stokes equation)來模擬震波通過背向式階梯流場時，探討在不 同階梯高度，不同壁面溫度及不同震波馬赫數之下，震波對壁面熱傳 效應之影響。

第二章 物理問題

2-1

物理模式

本文係以納維爾-史脫克方程式(Navier-Stokes equation)模擬二維 背向式階梯流場，當超音速震波通過一背向式階梯時，在階梯處因流 場區域增大導致壓力減小，而引起震波的膨脹，在階梯角落的地方由 於剪切層的作用而產生渦卷。背向式階梯流場之模擬區域示意圖如圖 2-1 所示。其中 h 為入口高度，s 為階梯高度，s 的變化為 ，

為入口區域至階梯的長度為 ， 為階梯底部至出口的長度為 ，

並 於 階 梯 底 部 壁 面 部 份 加 以 不 同 之 溫 度 變 化 ， ， 其 中

。

h h~4 5 .

0 L₁

h

2 L₂ 8h

0 0 ~10

2T T

K T₀ =300

2-2

基本假設

統 御 方 程 式 主 要 是 以 納 維 爾 - 史 脫 克 方 程 式 (Navier-Stokes equation)為依據，對流場的特性，作了下列幾點的假設：

(1) 流體為理想氣體。

(2) 為一黏性流體。

(3) 為一可壓縮流場。

以不同馬赫數的震波進入一背向式階梯流場，其震波之關係式如 下所示：

) 1 1(

1 2 ²

0

1 +

+ +

= M_s

P P

γ

γ (2-1)

1 ) 1 2 (

1

0

1 − +

= +

P M_s P

γ

γ (2-2)

⎟⎟

⎟⎟

⎠

⎞

⎜⎜

⎜⎜

⎝

⎛

− + +

− + +

=

0 1 0 1

0 1

0 1

1 1 1

1 1

P P P P

P P T T

γ γ γ γ

(2-3)

0 1 0 1

0 1

1 1

) 1( 1 1

P P P P

− + +

− + +

= γ γ

γ γ ρ

ρ (2-4)

上式下標定義如下：

0：表示震波前方之氣體狀態 1：表示震波後方之氣體狀態

第三章 數值方法

3-1

統御方程式

本研究使用二維納維爾-史脫克方程式(Navier-Stokes equation)為 统御方程式，其守恆方程式如下所示：

v v

v H

y F x H E y F x E t

Q r r r

r r r r

α

α +

∂ +∂

∂

= ∂

∂ + +∂

∂ +∂

∂

∂ (3-1)

當α =0為一個二維平面流場，當α =1為一個二維軸對稱流場。

式中，Qr

為守恆變數定義為：

⎥⎥

⎥⎥

⎦

⎤

⎢⎢

⎢⎢

⎣

⎡

= e

v Q u

ρ ρ ρ

r (3-2a)

上式中，ρ為密度，u、v分別為x與y方向速度分量，e為單位體積之 能量，可由理想氣體方程式得來。Er Fr

、 為非黏性通量，Hr

為非黏性 通量之源項定義如下：

(

^{2 2}

)

2 1

1 u v

e p + +

= − ρ

γ (3-2b)

( )

⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢

⎢⎢

⎢

⎣

⎡

+

= +

v P e

uv P u u

E ρ

ρ ρ r 2

(3-2c)

( )

⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢

⎢⎢

⎢

⎣

⎡

+

= +

v P e

P v uv u

F ₂

ρ ρ ρ

r (3-2d)

( )

⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢

⎢⎢

⎢

⎣

⎡

+

=

=

v p e v uv v

H y

ρ t

ρ ρ ρ

2

r 1

(3-2e)

v

v F

Er r

、 為黏性通量，Hr_v

則為黏性通量之源項定義如下：

⎥⎥

⎥⎥

⎥

⎦

⎤

⎢⎢

⎢⎢

⎢

⎣

⎡

=

x xy xx

Ev

β τ τ 0

r (3-3a)

⎥⎥

⎥⎥

⎥

⎦

⎤

⎢⎢

⎢⎢

⎢

⎣

⎡

=

y yy xy

Fv

β τ τ 0

r (3-3b)

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎦

⎤

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎣

⎡

⎟⎟⎠

⎜⎜ ⎞

⎝

⎛

∂

− ∂

⎟⎟⎠

⎜⎜ ⎞

⎝

⎛

∂

− ∂

−

⎟⎟⎠

⎜⎜ ⎞

⎝

⎛

∂

− ∂

⎟⎟⎠

⎜⎜ ⎞

⎝

− ⎛

−

⎟⎟⎠

⎜⎜ ⎞

⎝

⎛

∂

− ∂

=

y uv y x

y v y y

y v

y v y y

y v y v y x

H y

y yy xy

v

μ μ

μ β

μ μ

τ τ

μ τ

θθ

3 2 3

2 3

2

3 2 3

2 3 2 0

1

2 2

r (3-3c)

其中τ 為剪應力及β為熱通量

(

x y

)

x

xx λ u v μu

τ = + +2 (3-4a)

(

^{u v}

)

^μv

λ

τ = + +2 (3-4b)

(

x y

)

xy =μ u +v

τ (3-4c)

⎥⎦

⎢ ⎤

⎣

⎡ ⎟⎟−

⎠

⎜⎜ ⎞

⎝

⎛

∂ +∂

∂

= ∂

y v y v x

u 2

λ

τ_θθ (3-4d)

xy xx

x u v

x

k T τ τ

β + +

∂

− ∂

= (3-4e)

yy xy

y u v

y

k T τ τ

β + +

∂

− ∂

= (3-4f)

其中λ μ 3

−2

= ，μ為黏滯係數，k為熱傳係數。

黏滯係數以及熱傳係數與溫度之關係式如下：

S T

S T T

T

+

⎟⎟ +

⎠

⎜⎜ ⎞

⎝

≈⎛ ⁰

5 . 1

0

μ0

μ (3-5)

1 1 0 5 . 1

0

0 T S

S T T

T k

k

+

⎟⎟ +

⎠

⎜⎜ ⎞

⎝

≈⎛ (3-6)

其中

μ0:氣體之溫度為 T0之黏滯係數 k0: 氣體之溫度為 T0之熱傳係數

S:黏滯係數之 Sutherland 常數，S = 111k

S1:熱傳係數之 Sutherland 常數，S = 194 k

3-2 時間積分

本研究中選用二維平面流場，因此α =0，對統御方程式( 3-1 )之 時間積分，使用四階 Runge-Kutta scheme：

y F x E y F x E t

Q _{v v}

∂ +∂

∂

=∂

∂ +∂

∂ +∂

∂

∂r r r r r

^{Qr =Qrn + ΔtL}

( )

^Qrn

2

) 1

1 (

( )

( )¹

) 2 (

2 1 tLQ Q

Qr rⁿ r

Δ +

=

( )

( )²

) 3 (

2 1 tLQ Q

Qr rⁿ r

Δ +

=

( ) ( ) ( )

(

^{1 2 3}

) ( )

^{( )3}

1

6 2 1

3

1 Q Q Q Q tLQ

Qrⁿ rⁿ r r r r

Δ + + +

+

−

=

+

( )

^Q

dt L Q dr r

=

y F x F E

y F E

x E Q

L ^{v v}

j i j i j

i j

i ∂

+ ∂

∂ +∂

⎥⎦

⎢ ⎤

⎣

⎡ −

− Δ

⎥⎦

⎢ ⎤

⎣

⎡ −

−Δ

= + − + −

r r r

r r

r r

2 , 1 2 , 1 2,

, 1 2 1

1 ) 1

( (3-7)

3-3 空間差分

對於納維爾-史脫克方程式中對流項( convective terms )之空間差 分的近似值計算，本文是採用 WENO scheme ( Weighted Essentially Non-oscillatory scheme )，源自 1987 年由 Engquist 等人[13]提出 ENO scheme 在數值不連續處容易出現失真，所以 Liu 等人[14]與 Jiang 及 Shu[15]在 1994 及 1995 年分別對於 ENO 的差分方式進行平滑化 ( Smoothness )，其方法為在對流項Er(Q)

及Fr(Q)

進行空間差分時，對 於計算差分時所需的格點皆分配一個小小的重量( weight )ω 於其

上，以進行平滑化，並藉以由補差的方式使得精密度能夠達到( 2r-1 ) 階，本文計算時選則 r = 3，所以在空間上的精度可以到達 5 階。

為了求得方程式( 3-7 )中的

j i

E

2, +1

r 的高階近似，因此 Jiang 及 Shu [15]

定義：

∫

−^+ΔΔ

( )

= Δ ²

2

) 1 (

x x

x xh d

Q x

E^r ξ ζ

( ) [ ( ) ( ) ]

x x x x h

x h Q

E Δ

−Δ Δ −

= + 2 2

r (3-8)

接 著 對 全 流 場 進 行 通 量 分 離 ( flux splitting ) 本 文 係 採 用 Lax-Friedrichs ( LF ) scheme 其^Er(Q)之通量分離式如下：

) ( ) ( )

(Q E Q E Q

Er = r⁺ + r⁻

) ) ( 2( ) 1

(Q E Q Q

Er^± = r ±α

(3-9) 其中 maxEr'(Q)

α =

也 就 是 必 須 再 差 分 以 得 到 左 邊 邊 界及 右 邊 邊 界 的 數 值 通 量 ( numerical flux)。

) ( )

( )

(

2, , 1

2 , 1

2

1 Q E Q E Q

E

j i j

i j

i

− + +

+

+ = r + r

r (3-10)

在計算 ( )

2,

1 Q

E

j i

+ +

r 之正負數值通量時，採用下式以得到其近似值：

∑

⁻

= + + − + +

+

+ = ¹ ⋅ ⋅

0

, ,

1 , ,

2

1 ( ,..., )

r

k

j k i s j

r k i s r k s j k

i q l E l E

E^r ω (3-11)

∑

⁻

= − + −+ +

−

+ = ¹ ⋅ ⋅

0

, ,

1 , ,

2

1 ( ,..., )

r

k

j k i s j r k i s r k s j k

i q l E l E

E^r ω (3-12)

) ,...,

( _{1, 1,}

,s k s i r j s i r j

k⁺ =ω l ⋅E_{− +} l ⋅E_{+ −}

ω (3-13)

) ,...,

( _{2, ,}

,s k s i r j s i r j

k⁻ =ω l ⋅E_{− +} l ⋅E₊

ω

(

k=^0,1,...,r−¹

)

(3-14) 所以 ^±

+2 i 1

Er

可寫為：

∑

=

± +

±

+ = ^m ⋅

s

j s i j

i E r

E

1 ,

2 , 1

2 1

r r

(3-15)

上述的l_s及r_s為特徵矩陣。

而ω_{k ,}^±_s為一個對於計算差分時所需的格點分配之一個小重量，其 定義如下：

±

±

±

± ±

+

= +

2 1 0

, α α α

ω_ks α^k

(

k =0,1,2

)

(

0

)

²

0

1

+ +

= + ε IS α

(

1

)

²

1

6

+ +

= + ε IS α

(

2

)

²

2

3

+ +

= + ε IS α

對於上式的ε其存在的意義在於為了避免分母為零導致程式發 散，所以必須假設一極小量值，Jing 及 Shu[14]建議ε =1×10⁻⁶。

( ) (

2, 1, ,

)

²

2 , , 1 ,

2

0 4 3

4 2 1

12

13 _{+ +}

+ − + −

+− +−

+ = E_{i j} − E_{i j} +E_{i j} + E_{i j} − E_{i j} + E_{i j}

IS r r r r r r

( ) (

1, 1,

)

²

2 , 1 ,

, 1

0 4

2 1 12

13 ₊

+ + + −

+ + +−

+ = E_{i j} − E_{i j} +E_{i j} + E_{i j} −E_{i j}

IS r r r r r

( ) (

, 1, 2,

)

²

2 , 2 ,

1 ,

0 3 4 3

4 2 1

12

13 ₊

+ + + +

++ ++

+

+ = E_{i j} − E_{i j} +E_{i j} + E_{i j} − E_{i j} + E_{i j}

IS r r r r r r

(

0

)

²

0

3

−

−

= + ε IS α

(

1

)

²

1

6

−

−

= + ε IS α

(

2

)

²

2

1

−

−

= + ε IS α

( ) (

1, , 1,

)

²

2 , 1 ,

, 1

0 4 3

4 2 1

12

13 ₋

− +

−−

−+

−

−−

− = E_{i j} − E_{i j} +E_{i j} + E_{i j} − E_ij + E_{i j}

IS r r r r r r

( ) (

, ,

)

²

2 , 2 ,

1 ,

1 4

2 1 12

13 _{− − −}

− +

− +

− = E_{i j} − E_{i j} +E_{i j} + E_ij −E_{i j}

IS r r r r r

( ) (

1, 2, 3,

)

²

2 , 3 ,

2 ,

1

0 3 4 3

4 2 1

12

13 ₋

− +

− +

− +

− +

− +

− = E_i+ _j − E_{i j} +E_{i j} + E_{i j} − E_{i j} + E_{i j}

IS r r r r r r

所以( 3-8 )式可改寫為：

( ) ( )

( )

6 5 2

6

2 5

6

11 7

2

, 2 ,

1 ,

2

, 1 ,

, 1 1

, ,

1 ,

2 0 , 2 1

+ + +

+ +

+

+ + +

+

− + +

+

− +

− + +

+

− + +

+ +

+ − +

= −

j i j i j i

j i j i j i j

i j i j i j

i

E E

E Q

E E

E Q E

E E

E Q

r r

r

r r

r r

r r r

(3-16)

( ) ( )

( )

6 5 2

6

2 5

6

11 7

2

, 2 ,

1 ,

2

, 1 ,

, 1 1

, ,

1 ,

2 0 , 2 1

+ + +

+ +

+

+ + +

+

− + +

+

− +

− + +

+

− + +

+ +

+ − +

= −

j i j i j i

j i j i j i j

i j i j i j

i

E E

E Q

E E

E Q E

E E

E Q

r r

r

r r

r r

r r r

(3-17)

將式( 3-16)及( 3-17)代入式( 3-10 )中便可得

2 +1 i

Er

的正負通量 值，進而得到Ev_i,j

的值。至於Fr(Q)

之通量分離步驟與方法和Er(Q) 類 似，因此不再重覆敘述。

對於納維爾-史脫克方程式中黏滯項( viscous terms )之空間差 分的近似值計算，本文是採四階中央差分( forth-order centered difference scheme)[20]，其方程式定義如下：

( ) ( )

[

2 , 8 , 8 ( , ) ( 2 , )

]

12

1 u x h y u x h y u x h y u x h y h

dx

du = − − − + + − +

( ) ( )

[

^{, 2 8 , 8 ( , ) ( , 2 )}

]

12

1 u x y h u x y h u x y h u x y h h

dy

du = − − − + + − +

) , ( 30 ) , ( 16 ) , 2 ( 12 [

1

2 2

2 2

y x u y h x u y h x h u

dx u

d = − − + − −

+16u(x+h,y)−u(x+2h,y)]

) , ( 30 ) , ( 16 ) 2 , ( 12 [

1

2 2

2 2

y x u h y x u h y x h u

dy u

d = − − + − −

+16u(x,y+h)−u(x,y+2h)]

) , ( 10 ) , ( 10 ) , ( 10 24 [

1

2 2

h y h x u h y h x u h y h x h u

dxdy u

d = − − − − + − + −

+10u(x+h,y+h)−u(x−2h,y−h)+u(x−2h,y+h) +u(x+2h,y−h)−u(x+2h,y+h)−u(x−h,y−2h) +u(x−h,y+2h)+u(x+h,y−2h)+u(x+h,y+2h)]

)]

, 2 ( ) , ( 8 ) , ( 8 ) , 2 ( 12 [

1 x h y x h j x h j x h j

h dx

dμ = μ − − μ − + μ + −μ +

)]

, 2 ( ) , ( 8 ) , ( 8 ) , 2 ( 12 [

1 k x h y k x h j k x h j k x h j h

dx

dk = − − − + + − +

對 y 方向之差分類似於 x 方向之差分，因此不再重覆敘述。

因此統御方程式中 dx

E dr_v

及 dy F dr_v

可寫成下式：

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥

⎦

⎤

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢

⎣

⎡

+ +

+ +

+

+ +

+

+

− +

+

−

=

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎦

⎤

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎣

⎡

=

2 2 2

2 2

2

2 2

2 2

) (

) (

)]

3( 2 2

[ )]

3( 2 2

[ 0 0

dx T kd dx dT dx dk dx vd dx

dv dx ud dx

du

dxdy u d dx

v d dy

du dx dv dx d

dxdy v d dx

u d dx

u d dy

dv dx du dx

du dx d

dx d

dx d

dx d

dx E d

xy xy

xx xx

x xy xx

v

τ τ τ τ

μ μ μ μ

β τ τ r

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥

⎦

⎤

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢

⎣

⎡

+ +

+ +

+

+

− +

+

−

+ +

+

=

⎥⎥

⎥⎥

⎥⎥

⎥⎥

⎥

⎦

⎤

⎢⎢

⎢⎢

⎢⎢

⎢⎢

⎢

⎣

⎡

=

2 2

2 2 2

2 2 2 2 2

)]

3( 2 2

[ )]

3( 2 2

[

) (

) (

0 0

dy T kd dy dT dy dk dy vd dy

dv dy ud dy

du

dy v d dxdy

u d dy

v d dy

dv dx du dy

dv dy d

dy u d dxdy

v d dy

du dx dv dy d

dy d

dy d

dy d

dy F d

yy yy

xy xy

y yy xy

v

τ τ τ τ

μ μ μ μ

β τ τ r

3-4 起始條件與邊界條件

背向式階梯流場之示意圖如圖 2-1 所示。

3-4.1 起始條件：

進入背向式階梯的氣體條件為不同馬赫數換算出來的值，其換 算關係式可由方程式(2-1)~(2-4)得到。

背向式階梯內的氣體皆假設為靜止之空氣體，壓力為 1atm，密 度為 1.177kg/ m³，速度為 0 m/s。

3-4.2 邊界條件：

在入口處之邊界：其值為起始條件所設之狀態，並保持一固定 值。

壁面邊界採用無滑動邊界 ( non-slip boundry condition )，其壁面 上速度為零，除了L₂壁面以不同溫度變化假設為不同熱傳量之外，其 餘壁面皆假設為絕熱狀態。

流場出口邊界：在此是使用黎曼不變量邊界條件(Riemann invariant boundary condition)：此邊界法是假設一二維流場在邊界時可 視為近似一維自由流，因此根據特徵線理論，沿特徵線 R⁺及 R^-各有 一個黎曼不變量，其可表示為

2 + −

=

+ c

V R

1 2

− −

=

−

r V c

R _n

其中c為音速，V_n =V⋅nr，亦即表示垂直於流場邊界的速度分量，以向 外為正，R⁺及R⁻相對應的特徵速度為：

c V_n −

1 = λ

c V_n +

1 = λ

此時已有兩個方程式，但是為了求解四個變數值，所以還需要另外兩 個方程式，此時可以選擇熵(entropy)的關係式 _⎟⎟

⎠

⎜⎜ ⎞

⎝

= ⎛ ρ

S ln p 及切線速度(V_t) 為另外兩個關係式，來求解四個特徵值。

當超音速流出時(outflow)，λ₁及λ₂皆大於0，所以所有邊界值皆 由流場內部插分而得，當次音速流出時，λ₁<0，λ₂ >0，此時R⁻由自 由流條件而來，R⁺、V_t及熵(S)則由流場內部給定，如此將可求得各 變數之邊界值。

3-5 局部紐賽爾數及平均紐賽爾數(Local Nusselt number and Average Nusselt number)

為了瞭解不同加熱溫度對壁面熱傳效應之影響，本文使用局部紐 賽爾數(Local Nusselt number)以及平均紐賽爾數(Average Nusselt

number)來觀測震波對流場熱傳效應之影響，其定義如下：

0 2

) (

w

y

L

T T T x h

Nu ∂

∂

−

= −

(3-18)

∫

= ²

2 0

0 1 ^L ( )

dx x L Nu

Nu (3-19)

3-6 瑞理數(Rayleigh Number)

為了瞭解流場中浮力與黏滯力之相對關係，因此利用瑞理數來

分析，其定義如下：

αμ β

ρg H³(T T₀)

Ra ^w −

= (3-20)

=1

第四章 結果與討論

本研究主要是探討震波通過一背向式階梯會產生何種震波結 構，以及當壁面部份加上熱傳效應時，其熱傳效應對震波傳遞時有何 影響。由於數值模擬是否能正確預測流場結構，取決於使用之數值方 法的精度以及邊界條件的正確性，因此先就所撰寫之程式作程式驗證 與基本格點測試，以確定所撰寫之程式無誤，以及選用之格點能準確 的模擬出流場之結構。

4-1 程式驗證

模擬流場結構前為了確定所撰寫之程式是否正確，吾人將(3-1) 式之統御方程式選用α 之軸對稱型式，並忽略黏滯效應(亦即尤拉 方程式)做數值模擬，使用之格點數為665(x)×200(y)，以震波馬赫數 1.3 及震波管的管徑比為R/r=2.0之密度等位圖，於無因次化時間為 2.7 以及 5.35 時之結果與 Takayama 和 Jiang[1]之實驗與模擬結果做比 對，其結果如圖 4-1 及 4-2 所示。圖 4-1(a)與 4-2(a)為 Takayama 和 Jiang 利用全相干涉所拍攝出來的實驗結果。圖 4-1 為震波由震波管近入膨 脹管的初期等密度圖，其無因次化時間為 2.7，於圖中可見在階梯部 分由震波膨脹及剪切層做用所產生的渦卷，與震波接觸壁面後所造成 的反射產生交互作用，圖 4-2 為兩反射震波產生交互做用之後的密度

等位圖，其無因次化時間為 5.35，於圖中可見剪切層在靠近背向式階 梯附近，因為不穩定的擾動現象造成在主要渦卷後方產生小渦卷，以 及震波分別接觸到兩邊壁面造成的反射震波產生的交互作用，由以上 兩圖中可瞭解到本文所作的模擬與 Takayama 等人之模擬非常近似，

因此可證實程式撰寫無誤。

4-2 格點測試

流場中的格點數量對於計算結果影響很大，較多的格點數雖然 可以詳細的描述流場中各區域之物理量，但卻會增加運算的時間，過 少的格點則可能因為格點間的相互間距過大，而導致無法準確的描述 流場內各物理量的值，造成計算上的誤差。因此，基於精確度與時效 性的考量，在進行數值模擬計算之前，須先做格點測試。

為了研究二維背向式階梯流場之結構，所以選用(3-1)式統御方 程式中α =0之二維平面流型式，以階梯高度為 ，階梯底部至出口處 為L

h

2，分別選定三個不同位置，以及階梯底部壁面絕熱邊界條件與階

梯底部壁面加熱 條件，觀察其壓力波的分佈情形做格點測試。

其結果如圖 4-3 及圖 4-4 所示，各點誤差值如表 4-1 至 4-2 所示。點

一為 、

K 3000

h

x=3 y=h，點二為x=2h、y =0.5h，點三為x=3h、 ，其 中四種不同格點分別為 Grid 1：

=0 y

) ( 80 ) (

400 x × y ；Δx_min =Δy_min = 0250. h (2)

Grid 2：800(x)×160(y)；Δx_min=Δy_min =0.0125h (3) Grid 3： ； (4) Grid 4：

) ( 200 ) (

1000 x × y

h y

x_min =Δ _min=0.01

Δ 1333(x)×266(y)；Δx_min=Δy_min=0.0075h， 由表 4-1 及 4-2 結果顯示，階梯底部壁面絕熱邊界條件與階梯底部壁 面加熱3000K條件下，格點數800×160與1000×200及 與

於點一、點二、點三之誤差值分別為 、 、

與 、 、 ，由此可瞭解到以格點數 來模擬

會較為精確，由於格點數1333

1000×200

1333×266 0.108% 1.899% 0.054%

% 191 .

1 2.784% 0.112% 1333×266

×266每分鐘約運算約 3 個步階( time step )，若改用格點數 則每分鐘計算量提昇至 5-6 個，因此本

文最後選定以， 之格點數來做為本文數值模擬之格點

數。

1000×200

) ( 200 ) (

1000 x × y

4-3 流場結構分析

本文將就震波流場做一些介紹，並且比較壁面加熱是否對震波

結構產生影響，因此選定震波馬赫速 ，階梯高度 ，階梯底部壁

面分別為絕熱條件與壁面加熱10 之條件做比較計算結果分別為圖 4-5 所示。首先對流場一些現象做名稱的定義，如圖 4-5(a)-(j)所示，

起始震波(Precursor SW)、主要之渦卷(Primary Vortex Ring)、接觸表 面(Contact Surface)、剪切層(Shear Layer)、反射震波(Reflected SW)、

規則反射現象(Regular Reflection)、切斷震波現象(Split SW)、馬赫反 h

2 0

. 2

T0

射現象(Mach Reflection)、二次震波(Secondary SW)、第二渦卷

(Secondary VR)。以上所產生之各種現象在 Jiang 及 Takayama[1]在二 維軸對稱之震波管的數值模擬中亦可觀察到，吾人之模擬與二維軸對 稱模擬最大的不同點，在於吾人的模擬結果中不會出現馬赫圓盤 (Mach Disk)，也不會有震波的聚焦現象。

對於流場結構的分析本文將分成以下三個部份來討論，震波與 渦卷之交互做用(Shock Wave/Vortex interaction)、震波與剪切層之交互 做用(Shock Wave/Shear Layer interaction)、由規則反射震波轉換成馬 赫反射震波。

震波由於通過一背向式階梯，因體積膨脹導致壓力改變，引起 震波的膨脹，使得震波在背向式階梯的角落部份由於剪切層的影響產 生渦卷如圖 4-5(a)，當震波不斷向前傳遞會接觸到底部的壁面，進而 產生反射震波如圖 4-5(b)，當反射震波回傳接觸到渦卷如圖 4-5(c)，

會造成反射波與渦卷之交互做用如圖 4-5(e)，反射震波會在渦卷的部 份造成一切斷震波，震波會切開並通過渦卷然後再重新連，結此一現 象可由圖 4-5(g)觀察到。

震波與剪切層之交互做用由圖 4-5(e) ~ 圖 4-5(j)可觀察到，此一 現象可以在震波向剪切層擴展時，觀察到剪切層被從震波管角落被顯 現出來。推論原因是由於震波向剪切層擴展時，放射微粒在震波之後

速度增加，因此剪切層快速移動到震波中，此結果表示剪切層分離與 黏滯性無關。其他尚有一些震波與剪切層之交互現象如圖 4-5(e)，例 如一個渦卷在剪切層部分引起震波，並且此震波從剪切層分離了渦卷 如圖 4-5(j)，此一現象經常在自由噴流內觀察到，當流體有強烈的擴 張就會有產生的可能。

當一震波由規則反射轉變成馬赫反射由圖 4-5(e) ~ 4-5(g)可觀察 到，考慮一個反射震波在上方分裂，某種反射波的型式可能被觀察 到，這是由於反射波與分裂的部分其夾角為固定不變，在本文中可觀 察到起始震波(Precursor Shock Wave)的曲率越傳遞變的越來越小，所 以此一角度在壁面與震波間的變化將會引起反射震波的轉折。

當改變階梯高度時，對流場結構主要之影響有馬赫反射的易觀 察與否，當階梯高度越增加越多，反射波會較延遲接觸到壁面，其馬 赫反射的現象會越不明顯。另一影響則為當階梯高度越小時，反射波 會不斷在震波管內撞擊壁面，以及與管內之渦卷不斷交互做用，使其 結構更加複雜。

由圖 4-5 與 4-6 中可發現到，在階梯底部變化壁面溫度對震波結 構並不會有太大的影響，僅會在邊界部份產生熱邊界層。

0

4-4 階梯高度對紐塞爾數(Nusselt Number)之影響

在對流場的研究中，通常會將統御方程式無因次化，而為了減 少變數的數目，也會將一些變數結合在一起成為新的無因次化參數。

對流熱傳係數通常會被無因次化為紐塞爾數，其定義如方程式(3-18) 式所示。

本研究針對固定震波馬赫數(Ms=2.5)，改變五種不同的階梯高 度作探討，藉此用於觀察不同之階梯高度對於紐塞爾數分佈有何種影 響。如圖 4-7 及 4-8 所示發現，當階梯高度增加，其紐塞爾數會增大，

則對流就越有效率，其原因是由於特徵長度為階梯之高度，所以特徵

長度增加因此紐塞爾數會隨之增大，由此可知當階梯高度由 到

時，紐塞爾數會隨之增大。若改變特徵長度為階梯底部至出口處 之長(L

h 5 . 0

h 0 . 4

2)，則對紐塞爾數之影響參數變為溫度差，當溫度差越小時其 紐塞爾數會增大，反之則紐塞爾數會減小。

圖 4-7 為階梯高度由 到 時，不同溫度的階梯底部壁面瞬 間紐塞爾數對時間之分佈比較圖，由圖 4-7(a)與(b)中可發現，當階梯 高度為 以及 時，震波馬赫度為 ，溫度為 ，在 x 方向之加 熱壁面，會產生一紐塞爾數為負值之現象，此時不再是壁面將熱傳至 震波，而是由震波把熱傳向壁面。此一現象的發生是由於壁面之溫度 較震波之溫度低，如圖 4-9(a)與(b)所示，計算區域內近壁面之三點溫

h 5 .

0 4.0h

h 5 .

0 1.0h 2.5 2T

T0

0

0 w

度無因次化後高於壁面之溫度無因次化，因此其溫度梯度( )為

正值，導致紐塞爾數為負值，所以產生吸熱之現象，當階梯高度不斷 增加，其壁面之無因次化溫度將會高於震波之無因次化溫度，如圖 4-9(c)所示，就不再有溫度梯度為正之現象產生，所以紐塞爾數將不 會再有負值之現象。

Y T ∂

∂ /

由圖 4-7(b)到圖 4-7(e)中觀察到，當階梯高度增加時紐塞爾數會

有振盪的現象，尤其在階梯高度 與 時，震波馬赫數為 ，階

梯底部壁面溫度為 、3 、 ，有一明顯之紐塞爾數驟增現象，以

壁面溫度 最為顯著，其造成之原因是由於二次震波接觸到階梯底

部壁面時，二次震波之溫度較階梯底部壁面溫度低出許多，因此二次 震波與階梯底部壁面溫度之溫度梯度增大，造成紐塞爾數會驟增。

h 0 .

3 4.0h 2.5

2T0 4T₀ 2T

壁面溫度對紐塞爾數之變化，由圖 4-7 與 4-8 中可觀察到當壁面 溫度增加時紐塞爾數亦會增加，因為對紐塞爾數影響的參數除了特徵 長度之外，還有溫度差與溫度梯度會對其造成影響，因此當階梯底部 溫度增高會使得溫度差(T −T )增大，當溫度差增大時紐塞爾數會減 小，但是階梯底部壁面與震波間之溫度梯度較溫度差增加幅度更大，

因此使得紐塞爾數增大。

4-5 馬赫數對紐塞爾數之影響

欲探討不同階梯底部壁面加熱溫度，不同階梯高度，不同馬赫 數，對紐塞爾數之分佈有何不同之影響，所以選定三種不同馬赫數之

震波進入背向式階梯流場，並選用六種不同溫度( )

對階梯底部壁面做加熱，其中

0 0 0 0 0

0 3 4 6 8 10

2T、T、T、T、T、 T

K

T₀ =300 。溫度與瑞理數之關係為：當 溫度為2 ₀、 、4 、₀ 6 、₀ 及 時，瑞理數分別為 ⁷、 ⁷、

、 、 及 。

T 3T₀ T T 8T₀ 10T₀ 2.8×10 5.6×10

107

4 .

8 × 1.4×10⁸ 1.9×10⁸ 2.5×10⁸

由瑞理數定義(3-20)式會發現到，對瑞理數會有影響之主要參數 為溫度差以及特徵長度，於本研究中計算瑞理數所選用之特徵長度為 背向式階梯底部至流場出口之長度(L2)，在此對瑞理數影響之參數僅 有背向式階底部梯壁面與參考溫度的溫度差，因此在階梯底部溫度相 同情況下，其瑞理數亦相同，當階梯底部溫度增高，則瑞理數會增高。

圖 4-10 為階梯高度 到 時不同馬赫數，紐塞爾數對瑞里 數之關係，由圖中可明顯看出，當瑞里數較高時，紐塞爾數會較大。

進入背向式階梯區域內震波在固定瑞里數之情況下，震波馬赫數為 時之震波紐塞爾數會較小，其原因是由於瑞理數固定則階梯 底部壁面溫度為定值，則影響紐塞爾數之參數為溫度梯度，當震波馬

赫數為 ，震波與階梯底部壁面間之溫度差較小，所以其溫度

梯度較小，因此紐塞爾數會較小。由圖 4-10(d)與圖 4-10(e)中觀察到 h

5 .

0 4.0h

5 .

=2 Ms

5 .

=2 M_s