網格描述與座標系統

在使用ALE（Arbitrary Lagrangian-Eulerian）法計算移動邊界的問題中，有 三種座標定義域的關係存在，分別敘述如後，或參考圖2-1 之描繪[20,21]：

1、空間定義域（Spatial Domain，Ωy）：

代表一般流力問題所陳述的定義域，在探討移動邊界問題的領域中，此 定義域會移動，而在數值方法的格點劃分上，此定義域即為節點座標，因 此又稱為「Eulerian」座標系統。

2、物質定義域（Material Domain，Ωz）：

代表流體粒子所佔有的定義域且亦會移動，而在數值方法的計算中，此 定義域即為流體本身，其座標系統稱之為「Lagrangian」座標系統。

若是存在某物質定義域座標點zi，則其與空間定義域之座標點yi的關係 可由函數f1轉換如下：

( )

_{i t}

( )

_i

i f z t f z

y = ₁ , = ₁ （2-1a）

3、參考定義域（Referential Domain，Ωx）：

代表固定不動的定義域，因此在數值方法中為絕對的座標系統，通常稱 之為「Referential」座標系統。

若是存在某座參考定義域標點x，則其與空間定義域之座標i yi的關係可

由函數f2轉換如下：

1、流體速度（Material Velocity，u）

zi

3、對流速度（Convective Velocity，cˆ）

zi

1、Referential 座標與 Lagrangian 座標之轉換 根據微分定理可知

則

2、Eulerian 座標與 Lagrangian 座標之轉換 同理，根據微分定理可知

3、Referential 座標與 Eulerian 座標之轉換 同理，根據微分定理可知

將上式代入式（2-9）即可獲得 Referential 座標與 Eulerian 座標之轉換 關係式

2.2 ALE 中的統御方程式

根據轉換公式（2-10），將 Eulerian 座標之統御方程式（2-11）～（2-13）轉 換為ALE 形式之統御方程式，即可得到：

y 方向

Hughes [20,21]等對於 ALE 的運動學理論有詳細的探討，並定義了三種座標 定義域：在 ALE 法中，參考座標系統的移動速度是可以任意給定的，因此計算

f

1

f

₂

f

₃

第三章 物理模式

本研究將利用研究振動對過濾效率的影響，物理模型如圖 3-1，其中 CDIJ 為過濾材，本研究中濾材假設為多孔性介質中的纖維狀假設，BCJK 與 DEHI 為 韌性體，可使過濾材適當地移動，利用前章所述的ALE 法計算震動濾材的流場 計算，並利用粒子受到流場的影響下，計算粒子在流場中位置與速度，在利用分 子碰撞理論及單一纖維收集效率，計算粒子在濾材中的收集機率分佈與濾材的收 集效率。

3.1 流場統御方程式

為了簡化分析，本研究做了以下的假設：

（1） 工作流場為水平流場，不考慮重力方向。

（2） 工作流場為二維不可壓縮流場。

（3） 工作流體各項性質皆為常數。

（4） 流體與物體之界面滿足無滑動條件（no-slip condition），亦即移動面之流 體速度等於管道壁面之運動速度。

（5） 多孔性介質為纖維狀不可變形，且與流體不起化學作用。

（6） 多孔性介質具有等向性。

（7） 流場為二為穩定流場，不因粒子的加入而有所改變。

基於以上的假設，整個流場可以區分為多孔性介質內部之流場（簡稱內流 場），多孔性介質以外的流場（簡稱外流場）。

3.1.1 內流場之統御方程式 連續方程式：

=0

=0

=0

3.2 粒子運動之統御方程式

在本研究中，必須了解粒子在流場中的位置與速度，才能計算通過濾材時，

粒子在不同位置被捕捉的機率。由於微小粒子的運動極為複雜，為了便於簡化分 析，下列基本假設視為成立：

(1) 粒子為球體剛體，群體中個體運動仍視為單一球體剛體運動，運動中 相互間不發生碰撞。

(2) 流體的物理性質為定值。

(3) 不考慮溫度、濃度的影響。

(4) 不考慮粒子撞擊壁面後所產生之影響。

基於上述假設，粒子的運動方程式[24][25]可用下列各式表示粒子運動方程 式：

X 方向：

s CD ρf A vr

(

uf-us

)

dt

m du = ⋅ ⋅ ⋅ ⋅ 2

1

s (3-9a) Y 方向：

s D f r

(

f s

)

s C ρ A v v -v

dt

m dv = ⋅ ⋅ ⋅ ⋅ 2

1 (3-9b)

其中

² 4d^s A=π

8 3 4 _{s s}³

s

πd m = ρ

vr =

(

uf −up

) (

²+ vf −vp

)

²

為了數值計算的方便，導入下列無因次變數，其定義如下：

ν0

U_s =u^s ， ν0

V_s = v^s ，

ν0

U_f = u^f ，

ν0

V_f = v^f ， H τ=tv⁰

經無因次化的粒子運動方程式為：

X 方向：

r

(

f s

)

摩擦係數(Drag force Coefficience) C_D與相對速度、粒子大小有關，可由下 式表示： 率(single-collector efficiency)，計算接近過濾材的懸浮粒子被捕捉的機率，另一 部份計算懸浮粒子有多少機率會接近過濾材，本研究將兩者結合，計算懸浮粒子

因此研究中，將利用Zhu[34]所提出的擴散項、攔截項與慣性單一纖維收集 效率，作為本研究中的所須考慮的收集效率項。

擴散項(Diffusion)單一纖維收集效率

2/3

攔截項(Interception)單一纖維收集效率

2

慣性項(Inertia)單一纖維收集效率

( )

2

由於擴散項與慣性項不可能獨立出現，因此研究中利用Hinds[35]所提出擴 散項與攔截項單一纖維收集效率的修正項

3.3.2 懸浮粒子被收集效率

本節將利用分子碰撞理論[27]計算懸浮粒子有多少機率會接近過濾材，並利 用3.3.1 節所計算的結果，計算懸浮粒子在行進過程中，所被過濾材捕捉到的機 率，即是由單一纖維收集效率轉換成單一懸浮粒子被收集的機率。

在分子碰撞理論中，定義在距離ξ到ξ+dξ區間內，發生粒子發生碰撞的機 率為βd ，其中ξ β 為ξ到ξ+dξ 區間內的碰撞機率密度函數(probability density function)，基於此定義將可知，從原點行進到ξ+dξ 距離時，分子不被碰撞的機 率(F )可寫為

(

ξ dξ

)

F

( )(

ξ βd^ξ

)

F + = 1− (3-15) 對左式做Taylor 展開

( ) ( )

dξ F

( )

ξ β F

( )

ξ dξ dξ

ξ ξ dF

F + = − (3-16) 簡化後

( )

_{β F(ξ)}

dξ ξ

dF =− (3-17)

此解必須滿足初始條件F

( )

0 =1，因此可得

( )

ξ e ^βξ

F = ⁻ (3-18) 定義 f

( )

ξ 為 0 到ξ區間不發生碰撞，但在ξ時間點上發生碰撞的機率，因 此計算行進0 到ξ區間不發生碰撞，但在ξ到ξ+dξ 區間發生碰撞的機率，可寫 為

( )

^ξ d^ξ F

( )

^{ξ β}d^ξ

f = 因此可得

β βξ

ξ)= e⁻ (

f (3-19) 因此利用(3-19)計算懸浮粒子發生碰撞為止所行進的平均距離< ξ >

∫ ( )

∞

=

⋅

>=

<

0

1 dξ β ξ f ξ

ξ (3-20)

如圖 3-3 碰撞體積示意圖可知，其中 d 為懸浮粒子的有效碰撞直徑，懸浮

粒子在行進dt過程中，被捕捉的總粒子數，亦可視為所受到的碰撞次數，而碰 撞次數等於碰撞體積乘與單位體積纖維數(n)與總單一纖維收集效率(η )，下式_T 表示

T 2

12 η

πd v dt×

n v_r (3-21) 因此在行進過程中，分子的碰撞頻率ν 為

T 2 12 η π

ν =n d vv_r (3-22) 由碰撞頻率ν 與平均行進距離 ξ 的關係式為 ξ ν_s

v = 懸浮粒子平均行進 速度，可得

β η π

ξ ν ¹

n ₁₂² _T

s =

=

vr

d (3-23) 在本研究中

2 2

s s

s u v

v = + ，

( )

12 2

s

p d

d d +

=

(

s p

) (

² s p

)

²

r u u v v

v = − + − (3-24)

將(3-23)帶入(3-19)可得懸浮粒子在行進過程被收集的機率密度函數 ^f

( )

^ξ

表 3-1 不同

^R_er

值下，

C₀

、

C₁

、

C₂

之值

R er C ₀ C₁ C₂

<0.1 0 24 0 0.1~1 3.69 22.73 0.0903

1~10 1.222 29.1667 -3.889

10~10² 0.6167 46.5 -116.67

10²~10³ 0.3644 98.33 -2778 10³~5×10³ 0.357 148.62 -4.75×10⁴

5×10³~10⁴ 0.46 -490.546 -5.79×10⁵

3-1 物理模式

圖 3-2 單一纖維收集效率定義圖

圖 3-3 直徑為 d

₁₂

速度 v

_r

、單位時間 dt 碰撞體積示意圖

d12

vr dt

第四章 數值方法

本研究的數值方法將分為兩部份，第一部份將計算速度場，採用葛拉金有限 元素法（Galerkin finite element method）。所有元素均為八節點二次等參元素，並 引入處罰函數（penalty function）[28]處理壓力項和連續方程式，使所需計算的 變數僅為Ｕ、Ｖ，對於時間項則採用後向差分隱式法（backward different implicit method）。此外，應用牛頓拉斐遜（Newton-Raphson）迭代法[29]處理動量方程

( )

( )

則（4-10a）與（4-10b）式變為：

( )

為 了 計 算 處 理 上 的 方 便，可 以 把 原 本 非 線 性 積 分 方 程 式 加 以 線 性 化 ，

1

∫

(

1

) (

1

)

2

=1 ²

2 −ξ −η

N (4-22b)

(

1

)(

1

)(

1

)

4

= 1

3 +ξ −η ξ−η−

N (4-22c)

(

¹

) (

¹

)

2

=1 ²

4 +ξ −η

N (4-22d)

(

1

)(

1

)(

1

)

4

= 1

5 +ξ +η ξ+η−

N (4-22e)

(

1

) (

1

)

2

=1 ²

6 −ξ +η

N (4-22f)

(

−ξ

)(

+η

)(

+ξ−η

)

− 1 1 1

4

= 1 N7

(4-22g)

(

¹

) (

¹

)

2

= 1 ²

8 −ξ −η

N (4-22h)

如 此 即 可 獲 得 單 一 元 素 之 矩 陣 方 程 式

[ ]

^{( )}

[ ]

^{( )}

[ ]

^{( )}

(

^{C e + K e +λ L e}

) ^{{ }}

^{q ( )e =}

^{{ }}

^{f ( )e (4-23）}

其 中

{ }

^{q ( )e =}

[

^U₁^{, U}₂ ^{, .. . , U}₈ ^,V₁ ^,V₂^{, . . . ,V}₈

]

^t

[ ]

C ^{( )e} 表 非 線 性 迭 代 U 與 V 所 組 成 之 矩 陣 ，

[ ]

K ^{( )e} 表 純 由 幾 何 形 狀 函 數 與 時 間 項 所 組 成 之 矩 陣 ，

[ ]

L^{( )e} 表 帶 有 處 罰 函 數

λ

項 所 組 成 之 矩 陣 ，

{ }

f ^{( )e} 表右半部已知向量之矩陣。

若是將計算區域內的所有元素結合成單一矩陣方程式，則上式變為

[ ] [ ] [ ]

(

^{C + K +λ L}

){ } { }

^{q = f} (4-24）

至於各矩陣的詳細內容可以參考附錄A 所示。

由於（4-24）式為一組非常龐大的聯立方程組，為減小計算所需的記憶體空 間，本研究採用鋒面法（frontal method）[30-32]並配合高斯喬登消去法來求解。

另外，程式的收斂條件則定為

{ } { }

(

^{q m+1− q m}

) ^{{ }}

^{q m+1 < 10−3}

MAX (4-24）

而為避免數值計算中所造成的發散現象，求解高雷諾數的流場時，是利用所求得 的低雷諾數之流場作為初始值，再用此穩態流場之計算結果，當作求解暫態流場 之初始值。

數值計算流程圖如圖4-2 所示。詳細的步驟如下：

1. 在穩態流場下，經由網格測試決定最佳的計算網格分佈和元素數量。

2. 求得穩定流場下之流場，作為暫態解的初始值

3. 計算各網格點的網格速度，並檢查初始條件與邊界條件有無錯誤。

4. 求得所需的各項參數（矩陣、形狀函數等）。

5. 反覆解聯立方程組求速度場，直到每一格點的速度及溫度值滿足下列收 斂條件。

3 1

1

10⁻

+

+ − <

m m m

φ φ

φ ，式中

φ

表示U 、_f V 。 （4-26） _f

7. 繼續下一個時間的運算，直到達到預設的時間。

4.2 粒子運動軌跡與捕捉機率計算

利用上述方法求得流場分佈，代入式(3-10a)、(3-10b)則可計算粒子運動軌 跡。用數值方法聯立解(3-10a)、(3-10b)時，本文採用有限差分法(finite difference)，

以反覆計算求得粒子軌跡各位置速度U 、_s V 的收斂解。 _s

(

f X DXY DY s^mX DXY DY

)

r s s f D Y

X s m

DY Y DX X

s V d U U

d C H

U

U _,⁺¹_{+ , +} = _{, ,} + ⋅ _{, + , +} − _{, + , +} 4

3

τ

ρ

ρ

(4-27a）

(

f X DX Y DY s^mX DXY DY

)

r s s f D Y

X s m

DY Y DX X

s V d V V

d C H

V

V_,⁺₊¹ _{, +} = _{, ,} + ⋅ _{, + , +} − _{, + , +} 4

3

τ

ρ ρ

（4-27b）

其中

2

(

^{X DX,Y DY}

)

^{P X Y}

(

^P

(

^{X Y}

) )

^F_T

P + + = ( , )+ 1− , × 。 6. 判斷是否碰壁或離開過濾材。

7. 如果無將既序計算下一點速度位置，並重複計算下一點捕捉情形，

反之則停止。

總計算懸浮粒子運動軌跡與捕捉機率計算的流程圖如圖4-3

4-1 元素節點排列方式示意圖

圖 4-2 程式架構圖

開始

建網格

計算形狀函數與 Jacobian 矩陣

網格是否正確

解速度場

速度場是否收斂

輸出資料

移動網格點

是否需要網格重建 是

計算新網格點之 速度值

否

是否繼續計算

否 停止

結束 初始條件

網格點位置 邊界條件

是

是

否

是 否

是否為多孔性介質

使用外流場方程式 使用內流場方程式

否 是

4-3 粒子運動軌跡與捕捉機率計算的流程圖

開始

計算懸浮粒子所 在位置的流場速

計算位置是否碰壁 計算速度與位置

是否繼續計算

否 結束

初始條件

是

是

是 否

是否為過濾材

計算懸浮粒 子被收集效 否

第五章 結果與討論

本研究為了使過濾器設計人員更了解過濾器的原理，因此在本章第一節將 利用以上章節所推導出的公式，進行說明過濾各參數對過濾效率的影響。由於在 過濾器的設計中過濾器所產生的壓力差對設計風機時會有所影響，尤其是壓差過 大時，會使得風機設計的過於大，佔去過多廠房空間，因此壓差越小過濾效率越 高是過濾器的主要目的，因此本章第二節也會說明各參數對壓力的影響。在本章 第三節中將說明增加振動對過濾器的影響，並與Kim[20]的實驗值做比較，進行 討論。

5.1 過濾各參數對過濾效率的影響

對各參數的設計於表5.1 列出，研究中將對進口速度、孔隙率、纖維直徑、

懸浮粒子粒徑與過濾材厚度五項參數，作為影響對過濾效率的主要因素作為討 論，其中表5-1 所使用的參數值，為了與 Kim[20]所使用的參數值所計算出的結 果相呼應，因此使用與Kim[20]相同的參數值。

由圖 5.1 進口速度與收集效率關係圖中，可看出速度對總過濾效率(Total efficiency、solid line)、擴散項(Diffusion efficiency、dash-dash line)、攔截項 (Interception efficiency、dash-dot line)、慣性項(Inertial efficiency、dot line)、擴散 與攔截修正項(DR Interaction efficiency、gray line)的關係，其中圖中的各作用項 是利用各項”單獨存在”作用力時所造成的收集效率，總過濾效率以實線表示，擴

由圖 5.1 進口速度與收集效率關係圖中，可看出速度對總過濾效率(Total efficiency、solid line)、擴散項(Diffusion efficiency、dash-dash line)、攔截項 (Interception efficiency、dash-dot line)、慣性項(Inertial efficiency、dot line)、擴散 與攔截修正項(DR Interaction efficiency、gray line)的關係，其中圖中的各作用項 是利用各項”單獨存在”作用力時所造成的收集效率，總過濾效率以實線表示，擴

在文檔中 利用振動增加濾材之捕捉效率之理論研究 (頁 23-0)