本章將介紹及說明本研究中所使用的方法以及其步驟。第一節為 模型的介紹。將經簡化過後的板式熱交換器切割成若干個格點,以及 模型相關的何參數等。第二節為方程式的建立。參考 Bassiouny and Martin[9] 的 U 型熱交換器分流理論,並於分流區與匯流區在動量方 程式的部分做修正,並建立修正後的方程式。第三節為解題的流程。
將方程式分成流場及溫度場兩大部分,先解流場,再解溫度場,最後 將兩者的結果反覆互相迭代,以求得板式熱交換器內部真正的流動情 形。第四節為牛頓法(Newton method))的介紹。由於根據第二節所建 立的方程式,為一非線性聯立方程組,因此在此使用牛頓法來求解此 方程組。同時,針對牛頓法,介紹了牛頓法在應用上的缺點,並針對 缺點引入牛頓下降法來改良,以增加程式的收斂範圍。第五節為程式 流程。此節先分別介紹了流場與溫度場的程式流程,並於最後說明流 場及溫度場整合後的程式流程以及收斂條件等。
綜合上述,本研究在數值方面的計算過程為,首先根據使用者所 給定的板片數,包含熱通道與冷通道的數量,以及板片的幾何形狀等,
配合 Zahid [11]於 2003 所提及的板式熱交換器熱傳與壓降相關經驗式,
可以分別解出冷側及熱側流場,再將冷側及熱側流場進行熱交換,緊 接著重新更新冷側及熱側流場的速度、壓力及其它相關熱力性質,並
再次代入溫度場中,重覆步驟直到流場的速度及壓力等變數與上一次 計算完溫度場後相同時,則程式結束。
3-1
模型介紹首先根據 2-2 節的假設,流體為連續流體且流動型式為一維,忽略重 力效應所造成的影響,且已知入口壓力及溫度並與外界絕熱。接著根 據圖 1.13U 型板熱交換流場分佈,可建立流場模型。
由於冷熱流體並不互相接觸,且行走於相異通道,為了避免混肴,
在此將冷側與熱側通道的模型配置圖分開來討論。
冷側(工作流體為水)模型如圖 3.1 所示,其中左下方為入口;由於 為 U 型配置,故出口設於左上方。在模型上方及下方分別為匯流區 以及分流區;而在中間則為板片通道,此部分佔了流場約 80%以上的 體積,而根據 2-2 假設,絕大部分的壓降以及所有的熱傳現象,也都 發生在板片通道之間。
熱側(工作流體為二氧化碳)模型如圖 3.2 所示。
由於流動方向為逆向,故熱側模型的入口在左上角而出口在左下 角,上下及下方分別為分流區及匯流區,其餘的配置則與冷測模型相 同。
圖 3.1 冷側數值模型配置圖
圖 3.2 熱側數值模型配置圖
3-2
方程式建立本文方程式的建立為以 Bassiouny and Martin[9] U 型板式熱交換 器分流理論為基礎。Bassiouny and Martin[9] 的所解出的理論解雖然 可以符合大部分流體在板式熱交換器的分流情形,二氧化碳在實際應 用上仍有一些較不方便的地方:
(i) [9]在分流區(匯流區亦同)動量方程式的部分,設定平均速度 比𝛽為定值(如式 3-1),由於固定𝛽代表分流區及匯流區的損 失係數(loss coefficient)為定值,然而在真實情況下,損失係 數會因為流量不同而有所改變,因此並不能假設為定值。
方程式推導如下:
PA − �P +𝑑𝑃𝑑𝑑∙ ∆Z� A = ρA �W +𝑑𝑑𝑑𝑑 ∙ ∆Z�2− 𝜕𝐴𝑊2+ 𝜕𝐴𝑢𝑐β𝑊 (式 3-1) 其中
β = 𝑑𝑑𝐶 = 𝑢𝑐𝐴𝑐𝑑𝐴 𝑢𝑐 =𝛽𝑑𝐴𝐴
𝑐
若將(式 3-1)加入離散的概念,整理後可以得到
𝑃𝑖A − 𝑃𝑖+1𝐴 = 𝜕𝐴𝑊𝑖+12− 𝜕𝐴𝑊𝑖2+ 𝜕𝐴𝑐𝑢𝑐β𝑊𝑖 (式 3-2) 移項可得
𝑃𝑖A − 𝑃𝑖+1𝐴 = 𝜕𝐴𝑊𝑖+12− 𝜕𝐴(1 − 𝛽2)𝑊𝑖2 (式 3-3)
因此β2項相當於分流區損失係數(loss coefficient),在文獻中作者 將β假設為定值,但實際應用上,假設為定值並不恰當,可能導致無 法反應出二氧化碳板式熱交換器內部真實流動情形。
(ii) 流動的性質如密度、黏滯係數等在文獻中皆假設為常數。
但當二氧化碳位於超臨界區時,性質在溫度及壓力變化非 常劇烈,假設為常數則會容易造成相當大的誤差。
以冷側為例,針對以上兩點,建立並修正後的方程式及控制體積 如下:
圖 3.3 分流區控制體積示意圖
圖 3.4 匯流區控制體積示意圖
圖 3.5 流道控制體積示意圖
方程式如下:
除此之外,根據[9],本研究還考慮轉入(出)流道的能量損失(Turning loss)如下:
(i) ,1
A ( 1 )
,1 2,10
i i
in ch c t i i
P A − P + + C ρ Au =
(ii) ,
( 1 )
, 2,0
i N i
ch c out to i N i N
P A − P A + + C ρ Au =
另外在熱場的部分,由於從入口到流道之間的距離非常短,因此 假設熱傳僅發生在板片通道的地方,至於方程式的建立,則參考[10]
的模型,所建立的控制體積如下圖所示:
圖 3.6 冷側熱場控制體積示意圖
圖 3.7 熱側熱場控制體積示意圖
方程式如下:
3-3
解題流程3-4
牛頓法的介紹牛頓法[12],[13](Newton method)一般又稱為切線法,為非線性方 程式相當重要的一種方法,許多解非線性方程式所使用的迭代法,皆 以牛頓法為基礎而發展而來,以下將對牛頓法有約略的介紹。
假設有一非線性方程式f(𝑥) = 0, 若存在𝑥∗使得f(𝑥∗) = 0,則稱𝑥∗ 為f(𝑥) = 0的一解。
今假設𝑥0為f(𝑥)的一近似解,將f(𝑥)在𝑥0處做泰勒展開(Taylor expension),可得
𝑓(𝑥) = f(𝑥0) + (𝑥 − 𝑥0)𝑓′(𝑥0) +(𝑥 − 𝑥0)2
2! 𝑓′′(𝑥0) … 若忽略高次項,則上式可改寫為
𝑓(𝑥) = f(𝑥0) + (𝑥 − 𝑥0)𝑓′(𝑥0) 移項可得
x = 𝑥0− 𝑓(𝑥0) 𝑓′(𝑥0) 此為牛頓法的迭代式。
3-4-1 牛頓-拉弗森法(Newton-Raphson method)
其中等號右側
𝑓1 = 𝑓1�𝑥1(𝑘), 𝑥2(𝑘), … 𝑥𝑖(𝑘)�, 𝑓2 = 𝑓2�𝑥1(𝑘), 𝑥2(𝑘), … 𝑥𝑖(𝑘)� … 等。
等號左側𝑓1𝑥1 =𝜕𝑥𝜕𝑓1
1, 𝑓1𝑥2 = 𝜕𝑥𝜕𝑓1
2… … 𝑓𝑖𝑥𝑗 =𝜕𝑥𝜕𝑓𝑖
𝑗…。此一矩陣又稱為 Jacobian 矩陣。
而𝑥𝑖𝑘表示𝑥𝑖在第 k 次迭代時的猜值,而∆𝑥𝑖𝑘即為修正值,即
𝑥1(𝑘+1) = 𝑥1𝑘 + ∆𝑥1𝑘 , 𝑥2(𝑘+1) = 𝑥2𝑘 + ∆𝑥2𝑘 , 𝑥3(𝑘+1) = 𝑥3𝑘 + ∆𝑥3𝑘… … 而在矩陣中間偏微分的部分本研究則採用中央插分法來表示:
(𝑓1)𝑥1 = 𝑓1(𝑥1+ 𝛿, 𝑥2, 𝑥3, … ) − 𝑓1(𝑥1− 𝛿, 𝑥2, 𝑥3, … ) 𝛿
(𝑓1)𝑥2 = 𝑓1(𝑥1, 𝑥2+ 𝛿, 𝑥3, … ) − 𝑓1(𝑥1, 𝑥2−𝛿, 𝑥3, … )
⋮ 𝛿
(𝑓𝑖)𝑥𝑗 = 𝑓𝑖�𝑥1, 𝑥2, 𝑥3, … 𝑥𝑗 + 𝛿, … � − 𝑓𝑖�𝑥1, 𝑥2−𝛿, 𝑥3, … 𝑥𝑗 + 𝛿, … � 𝛿
由於牛頓法為局部收斂,當初始猜值在根附近時,收斂速度相當 快速,而當猜值與根的距離相對較遠時,則容易會有無法收斂的情形 發生,而實際上,並無法保證初始猜值的位置,這使得牛頓法應用在 解較複雜的非線性方程組時,有相當大的困難。
3-4-2 牛頓下降法
為了解決牛頓法過於依賴猜值所造成的問題,牛頓下降法便是針 對未能給定較佳的初始猜值時,所對牛頓法的一種修正。此法乃是藉 由增加牛頓法的收斂範圍來降低牛頓法對猜值的依賴性,此法在應用 上非常的重要。
若將迭代式
𝑥�(𝑘+1) = 𝑥�(𝑘)− �∇𝑓̃(𝑘)�−1𝑓̃�𝑥�(𝑘)� 修正為
𝑥�(𝑘+1) = 𝑥�(𝑘)− 𝑤�∇𝑓̃(𝑘)�−1𝑓̃�𝑥�(𝑘)� , 0 < 𝑤 ≤ 1
其中 w 為一參數,選擇適當的 w 滿足
�𝑓̃�𝑥�(𝑘+1)�� < �𝑓̃�𝑥�(𝑘)��
一般而言,w 依次為 1 ,1
2 ,1
22,…,此稱為牛頓下山法。
其中
𝑓̃ = (𝑓1, 𝑓2, … 𝑓𝑖) 𝑥� = (𝑥1, 𝑥2, … 𝑥𝑖)
3-5
程式流程在程式開始之前,首先必須把水及二氧化碳的熱力性質繪出 [4](如圖 1.9-圖 1.12)。接著找出可以描述此一方程式曲線並建構在程 式中,以便程式在迭代時更新熱力性質使用。另外,根據 3-2 節所述,
程式開始之前,也必須將分流區及匯流區不同截面積及流量之下所對 應的損失係數[14] (如圖 3.8、圖 3.9)建立於程式中。
圖 3.8 分流區不同流速下主流方向所對應的損耗係數[14]
圖 3.9 匯流區不同流速下主流方向所對應的損耗係數[14]
流程圖(流場):
圖 3.10 流場程式流程圖
② 輸入板片內速度 及壓力(猜值)
⑥利用牛頓下降法計 算出新的速度及壓力
值代回③中
⑦當⑤符合時,表示連續及動量 方程式均收斂,猜值即為解。
程式結束
① 程式開始
④由輸入的板 片數量建立方 程式(f1~fn)
⑤∑ |𝑓𝑖𝑖=1 𝑖|≤收斂標準? NO
③由速度及壓力計算內 部流體的性質
YES
流程圖(溫度場):
圖 3.11 熱場程式流程圖
② 輸入板片內熱通 道及冷通道流速
⑤利用牛頓法求解溫 度場
⑧程式結束
① 程式開始
④由熱場方程式建 立矩陣
③由速度計算 Nusselt Number 等熱力性質,輸
入溫度猜值
⑥判斷是否滿足能量方程式
⑦利用牛頓下降法計 算出新的溫度猜值代
回③中
②輸入板片內速 度、壓力及初始溫 度,並得到流體性質
③計算水的流場 ④計算二氧化碳的流場
⑤計算熱交換器內之溫度場
⑥判斷流場速度與壓力與上 次計算,其差值是否收斂標準
⑧程式結束
① 程式開始
流程圖(流場及熱場整合)
圖 3.12 整合後程式流程圖
NO
⑦更新熱場溫度並代入②中 重覆此動作直至⑥符合收斂
標準
YES