一具有蠕動壁面的管道流分析

國 立 交 通 大 學

機械工程學系

碩 士 論 文

一具有蠕動壁面的管道流分析

Numerical Analysis of the flow in a channel

with a peristaltic wall

研 究 生：郭大慶

指導教授：崔燕勇 教授

一具有蠕動壁面的管道流分析

Numerical Analysis of the flow in a channel with a peristaltic wall

研 究 生：郭 大 慶 Student：Da - Ching Guo 指導教授：崔 燕 勇 Advisor：Yeng-Yung Tsui

國立交通大學 機械工程學系

碩士論文

A Thesis

Submitted to Department of Mechanical Engineering Collage of Engineering

National Chiao Tung University In Partial Fulfillment of the Requirements

For the degree of Master of Science

In

Mechanical Engineering July 2011

i 一具有蠕動壁面的管道流分析 研究生：郭大慶 指導教授：崔燕勇 博士 國立交通大學機械工程學系 摘要 本研究乃利用非結構性移動網格以有限體積法求解一個二維具有蠕動 的上管壁的管道流，並探討在上方管壁的蠕動波以不同的振幅比、波數進 行變形時與不同的雷諾數和背壓對流場的平均流量所造成的影響，同時經 由與他人的實驗數據比較可知，模擬計算的結果與實驗數據之間的差距甚 小。而由計算的結果得知不論在高背壓或低背壓的情況下，藉著提高振福 比或者波數皆可使得流量的輸出獲得改善與提高。並且在零背壓的情況 下，隨著振幅比的提昇則雷諾數對流量的影響幾乎是可以忽略的。另外在 固定振幅比的時候，雷諾數在零背壓與低波數的情況下對流量的影響是可 以忽略的，但是在高波數的時候流量則穩定的隨雷諾數增加。並且根據計 算結果可知，流量與背壓之間是呈線性遞減的關係，流量隨著背壓而減少。

ii

Numerical analysis of the flow in a channel with a peristaltic wall

student：Da-Ching Guo advisor：Dr. Yeng-Yung Tsui

Department of Mechanical Engineering National Chiao Tung University

ABSTRACT

A 2-D channel flow with a peristaltic wall was solved by using unstructured moving grid and Finite Volume method. We calculate the flow when there are several peristaltic waves on the upper wall with different wave number and amplitude ratio under different back pressure and Reynolds number. And by the comparisons of numerical results and the experiments and theory, we found the differences between the numerical results and experiments are very small. In addition, we observed that the average flow rate is increased when we increase the amplitude ratio or wave number, no matter the flow is under a situation of higher back pressure. And note that when the back pressure is zero and the wave number is low, the effects of Reynolds number to average flow rate could be ignored as the amplitude ratio increases and big enough, otherwise the average flow rate is increased with Reynolds number when the wave number is big enough. Also, according to the numerical results, the average flow rate is decreased linearly with back pressure.

iii

誌謝

感謝崔燕勇教授在我漫漫的碩士生活裡對我的研究上的一些問題的指 導，其中最重要的是這幾年裡崔老師讓我理解了做研究最重要就是嚴謹的 態度，但其實不僅僅是做研究，我想日後面對人生的各種問題時也應該以 一顆嚴謹的心來面對。另外我還要感謝胡育昌學長、林仕文學長、張祚昌 學長、陳俊佑學長、吳孝修學長在我有任何問題時都能不吝提供幫助，還 有陳信宏以及何光桓在我感到苦悶的時候替我打氣加油，雖然和你們相處 的日子不長，但我還是非常感謝你們。我還要感謝賴胤男、黃義政、林子 翔、鄭凱鴻、吳奉起、黃裕堂學弟以及陳虹汝學妹，你們驅散了沉悶的氣 息，讓整間實驗室充滿了生氣，有你們在實驗室總是不會感到無聊。最後 我要特別感謝我的家人讓我無後顧之憂的完成學業，雖然我回家的時間裡 總是在房間裡睡大覺，但是每次吃完飯要坐車回實驗室時爸總是令我有很 深的感觸，我長這麼大了確還是受到家裡跟其他人這麼多的照顧，雖然我 現在還不是個靠得住，但是我有信心讓自己成為一個可以讓人信得過靠得 住的人，謝謝你們。

iv 目錄 中文摘要………...…i 英文摘要………..…ii 誌謝……….iii 目錄……….…iv 表目錄………vii 圖目錄………...viii 符號說明………...xiv 第一章 緒論……….………1 1.1 蠕動式傳輸現象與應用……….………1 1.2 文獻回顧……….2 1.3 研究內容……….6 第二章 數學分析模型………8 2.1 物理模型………8 2.1.1 流場的基本假設條件………...8 2.2 ALE 座標系統的數學模型………..8 2.2.1 連續方程式………...9 2.2.2 動量方程式………...9 2.2.3 移動邊界………...9

v 2.3 無因次化………9 2.4 邊界條件………...10 2.4.1 出/入口邊界條件………10 2.4.2 牆邊界……….10 第三章 空間的守恆………..12 3.1 移動網格連續方程式………...12 3.2 空間守恆定理………...13 第四章 統御方程式的離散………..15 4.1 有限體積法及離散方程式………...15 4.1.1 動量方程式……….15 4.1.2 非穩態項……….16 4.1.3 對流項………..16 4.1.4 擴散項……….17 4.1.5 壓力項……….18 4.2 線性代數方程式………..18 第五章 數值方法………..19 5.1 PISO 演算法………..19 5.1.1 預測步驟……….20 5.1.2 第一次修正步驟………21

vi 5.1.3 第二次修正步驟………22 5.2 出/入口壓力邊界條件的流量計算……….23 5.2.1 出/入口壓力邊界條件的速度計算………24 5.3 PISO 演算法流程………..25 第六章 結果與討論………..26 6.1 網格測詴………26 6.2 到達周期性穩態的時間測詴………...26 6.3 不同的管長的流量測詴………...27 6.4 數值方法驗證………...27 6.5 振幅比的影響………..28 6.6 雷諾數的影響………..29 6.7 波數的影響………..30 6.8 背壓的影響………..31 第七章 結論………..32 參考文獻………34 附錄………36

vii

表目錄

表(6-1) 網格測詴…..………39 表(6-2) 週期流量表………..39 表(6-2) 不同管長的流量測詴………..39

viii 圖目錄 圖(2-1) 物理模型…...………...40 圖(3-1) 隱式法誤差區塊…...………...41 圖(3-2) 顯示法務差區塊…...………...41 圖(4-1) 擴散項離散化………..42 圖(5-1) 壓力邊界的流量計算表示………..42 圖(6-1) 平均流量隨時間變化(高波數)..………..43 圖(6-2) 出口瞬時流量隨時間變化(高波數)..………..43 圖(6-3) 數值方法驗證………..44 圖(6-4) 零背壓下，振幅比對流量的影響………...45 圖(6-5) 不同背壓下，振幅比對流量的影響………45 圖(6-6a) 背壓

 

p -0.5、振幅比

 

 0.1 流場於 0T 的壓力圖………...…46 圖(6-6b) 背壓

 

p -0.5、振幅比

 

 0.1 流場於 1/4T 的壓力圖………46 圖(6-6c) 背壓

 

p -0.5、振幅比

 

 0.1 流場於 2/4T 的壓力圖…...………….46 圖(6-6d) 背壓

 

p -0.5、振幅比

 

 0.1 流場於 3/4T 的壓力圖………46 圖(6-6e) 背壓

 

p -0.5、振幅比

 

 0.1 流場於 1T 的壓力圖………...…46 圖(6-7a) 背壓

 

p -0.5、振幅比

 

 0.1 流場於 0T 的流線圖…………...……47 圖(6-7b) 背壓

 

p -0.5、振幅比

 

 0.1 流場於 1/4T 的流線圖……….47 圖(6-7c) 背壓

 

p -0.5、振幅比

 

 0.1 流場於 2/4T 的流線圖….………47

ix 圖(6-7d) 背壓

 

p -0.5、振幅比

 

 0.1 流場於 3/4T 的流線圖…....…………47 圖(6-7e) 背壓

 

p -0.5、振幅比

 

 0.1 流場於 1T 的流線圖…..….…...…….47 圖(6-8a) 背壓

 

p -0.5、振幅比

 

 1.0 流場於 0T 的壓力圖………...…..….48 圖(6-8b) 背壓

 

p -0.5、振幅比

 

 1.0 流場於 1/4T 的壓力圖..……….48 圖(6-8c) 背壓

 

p -0.5、振幅比

 

 1.0 流場於 2/4T 的壓力圖……..……….48 圖(6-8d) 背壓

 

p -0.5、振幅比

 

 1.0 流場於 3/4T 的壓力圖……...………48 圖(6-8e) 背壓

 

p -0.5、振幅比

 

 1.0 流場於 1T 的壓力圖..………...…….48 圖(6-9a) 背壓

 

p -0.5、振幅比

 

 1.0 流場於 0T 的流線圖..………...…….49 圖(6-9b) 背壓

 

p -0.5、振幅比

 

 1.0 流場於 1/4T 的流線圖..……….49 圖(6-9c) 背壓

 

p -0.5、振幅比

 

 1.0 流場於 2/4T 的流線圖..……….49 圖(6-9d) 背壓

 

p -0.5、振幅比

 

 1.0 流場於 3/4T 的流線圖..……….49 圖(6-9e) 背壓

 

p -0.5、振幅比

 

 1.0 流場於 1T 的流線圖..………...…….49 圖(6-10a) 背壓

 

p 1.5、振幅比

 

 0.1 流場於 0T 的壓力圖…………...…...50 圖(6-10b) 背壓

 

p 1.5、振幅比

 

 0.1 流場於 1/4T 的壓力圖………...50 圖(6-10c) 背壓

 

p 1.5、振幅比

 

 0.1 流場於 2/4T 的壓力圖………...50 圖(6-10d) 背壓

 

p 1.5、振幅比

 

 0.1 流場於 3/4T 的壓力圖………...50 圖(6-10e) 背壓

 

p 1.5、振幅比

 

 0.1 流場於 1T 的壓力圖…………...…...50 圖(6-11a) 背壓

 

p 1.5、振幅比

 

 0.1 流場於 0T 的流線圖…………...…...51 圖(6-11b) 背壓

 

p 1.5、振幅比

 

 0.1 流場於 1/4T 的流線圖………...51

x 圖(6-11c) 背壓

 

p 1.5、振幅比

 

 0.1 流場於 2/4T 的流線圖………...51 圖(6-11d) 背壓

 

p 1.5、振幅比

 

 0.1 流場於 3/4T 的流線圖………...51 圖(6-11e) 背壓

 

p 1.5、振幅比

 

 0.1 流場於 1T 的流線圖………...……...51 圖(6-12a) 背壓

 

p 1.5、振幅比

 

 1.0 流場於 0T 的壓力圖…………...……52 圖(6-12b) 背壓

 

p 1.5、振幅比

 

 1.0 流場於 0/4T 的壓力圖………52 圖(6-12c) 背壓

 

p 1.5、振幅比

 

 1.0 流場於 2/4T 的壓力圖………52 圖(6-12d) 背壓

 

p 1.5、振幅比

 

 1.0 流場於 3/4T 的壓力圖………52 圖(6-12e) 背壓

 

p 1.5、振幅比

 

 1.0 流場於 1T 的壓力圖…………...……52 圖(6-13a) 背壓

 

p 1.5、振幅比

 

 1.0 流場於 0T 的流線圖………...…53 圖(6-13b) 背壓

 

p 1.5、振幅比

 

 1.0 流場於 1/4T 的流線圖………53 圖(6-13c) 背壓

 

p 1.5、振幅比

 

 1.0 流場於 2/4T 的流線圖………53 圖(6-13d) 背壓

 

p 1.5、振幅比

 

 1.0 流場於 3/4T 的流線圖………53 圖(6-13e) 背壓

 

p 1.5、振幅比

 

 1.0 流場於 1T 的流線圖………...…53 圖(6-14) 零背壓下，雷諾數對流量的影響………54 圖(6-15) 不同背壓下，雷諾數對流量的影響………54 圖(6-16) 零背壓下，波數對流量的影響………55 圖(6-17a) 波數

 

 0.03125、雷諾數

 

Re 20 流場於 0T 的壓力圖………...…56 圖(6-17b) 波數

 

 0.03125、雷諾數

 

Re 20 流場於 1/4T 的壓力圖…………56 圖(6-17c) 波數

 

 0.03125、雷諾數

 

Re 20 流場於 2/4T 的壓力圖…………56

xi 圖(6-17d) 波數

 

 0.03125、雷諾數

 

Re 20 流場於 3/4T 的壓力圖…………56 圖(6-17e) 波數

 

 0.03125、雷諾數

 

Re 20 流場於 1T 的壓力圖………...…56 圖(6-18a) 波數

 

 0.03125、雷諾數

 

Re 20 流場於 0T 的流線圖……...……57 圖(6-18b) 波數

 

 0.03125、雷諾數

 

Re 20 流場於 1/4T 的流線圖…………57 圖(6-18c) 波數

 

 0.03125、雷諾數

 

Re 20 流場於 2/4T 的流線圖…………57 圖(6-18d) 波數

 

 0.03125、雷諾數

 

Re 20 流場於 3/4T 的流線圖…………57 圖(6-18e) 波數

 

 0.03125、雷諾數

 

Re 20 流場於 1T 的流線圖………...…57 圖(6-19a) 波數

 

 0.28125、雷諾數

 

Re 20 流場於 0T 的壓力圖………...…58 圖(6-19b) 波數

 

 0.28125、雷諾數

 

Re 20 流場於 1/4T 的壓力圖…………58 圖(6-19c) 波數

 

 0.28125、雷諾數

 

Re 20 流場於 2/4T 的壓力圖…………58 圖(6-19c) 波數

 

 0.28125、雷諾數

 

Re 20 流場於 3/4T 的壓力圖…………58 圖(6-19c) 波數

 

 0.28125、雷諾數

 

Re 20 流場於 1T 的壓力圖…………...58 圖(6-20a) 波數

 

 0.28125、雷諾數

 

Re 20 流場於 0T 的流線圖………...…59 圖(6-20b) 波數

 

 0.28125、雷諾數

 

Re 20 流場於 1/4T 的流線圖…………59 圖(6-20c) 波數

 

 0.28125、雷諾數

 

Re 20 流場於 2/4T 的流線圖…………59 圖(6-20c) 波數

 

 0.28125、雷諾數

 

Re 20 流場於 3/4T 的流線圖…………59 圖(6-20c) 波數

 

 0.28125、雷諾數

 

Re 20 流場於 1T 的流線圖………...…59 圖(6-21) 不同背壓下，波數對流量的影響………60 圖(6-22a) 背壓

 

p 1.5、波數

 

 0.03125 流場於 0T 的壓力圖…………...61

xii 圖(6-22b) 背壓

 

p 1.5、波數

 

 0.03125 流場於 1/4T 的壓力圖…………...61 圖(6-22c) 背壓

 

p 1.5、波數

 

 0.03125 流場於 2/4T 的壓力圖…………...61 圖(6-22d) 背壓

 

p 1.5、波數

 

 0.03125 流場於 3/4T 的壓力圖…………...61 圖(6-22e) 背壓

 

p 1.5、波數

 

 0.03125 流場於 1T 的壓力圖………...…...61 圖(6-23a) 背壓

 

p 1.5、波數

 

 0.03125 流場於 0T 的流線圖……...……...62 圖(6-23b) 背壓

 

p 1.5、波數

 

 0.03125 流場於 1/4T 的流線圖…………...62 圖(6-23c) 背壓

 

p 1.5、波數

 

 0.03125 流場於 2/4T 的流線圖…………...62 圖(6-23d) 背壓

 

p 1.5、波數

 

 0.03125 流場於 3/4T 的流線圖…………...62 圖(6-23e) 背壓

 

p 1.5、波數

 

 0.03125 流場於 1T 的流線圖………...…...62 圖(6-24a) 背壓

 

p 1.5、波數

 

 0.28125 流場於 0T 的壓力圖………...…...63 圖(6-24b) 背壓

 

p 1.5、波數

 

 0.28125 流場於 1/4T 的壓力圖…………...63 圖(6-24c) 背壓

 

p 1.5、波數

 

 0.28125 流場於 2/4T 的壓力圖…………...63 圖(6-24c) 背壓

 

p 1.5、波數

 

 0.28125 流場於 3/4T 的壓力圖…………...63 圖(6-24c) 背壓

 

p 1.5、波數

 

 0.28125 流場於 1T 的壓力圖…………...63 圖(6-25a) 背壓

 

p 1.5、波數

 

 0.28125 流場於 0T 的流線圖……...……...64 圖(6-25b) 背壓

 

p 1.5、波數

 

 0.28125 流場於 1/4T 的流線圖…………...64 圖(6-25c) 背壓

 

p 1.5、波數

 

 0.28125 流場於 2/4T 的流線圖…………...64 圖(6-25c) 背壓

 

p 1.5、波數

 

 0.28125 流場於 3/4T 的流線圖…………...64 圖(6-25c) 背壓

 

p 1.5、波數

 

 0.28125 流場於 1T 的流線圖…………...64

xiii

圖(6-26) 不同波數下，背壓對流量的影響………65 圖(6-27) 不同振幅下，背壓對流量的影響………65

xiv

符號說明

v ：流體速度 g v ：網格速度 c v ：對流速度 wall v ：壁面速度  ：流體密度  ：流體黏滯係數 p ：流場壓力 t ：時間 1 y ：上壁面邊界移動方程式 h ：一半的入口高 A ：振幅 f ：頻率  ：波長 c _：波速 * v ：無因次流體速度 * g v ：無因次網格速度 * 1 y ：無因次上壁面邊界移動方程式 * p ：無因次壓力

xv * t ：無因次時間  ：振幅比  ：波數 Re ：雷諾數 T ：週期

1 第一章 緒論 1.1 蠕動式輸送現象(peristaltic transport)與應用： 蠕動式傳輸現象(Peristaltic Transport)，所謂蠕動式傳輸現象即在充滿流 體的管道內於管壁上產生周期性的、連續的波動來達成流體輸送的目的， 即管內的流動為由連續的波動在壁面所產生的位移造成。藉由此種傳遞波 動來完成流體的輸送的現象可以在大部份的生物系統上發現，比如說食道 的吞嚥食物、尿道以及腸道的蠕動。 蠕動式傳輸現象除了應用在流體輸送與幫浦上外，在 Yi et al. [12]中亦 指出其對精密實驗中的晶圓加工過程中的微流道化學反應具有促進的效 果。由於化學反應大都為兩種以上不同的詴劑進行混合與相互作用，而在 以往的加工過程中此部分大都只依賴擴散效應來完成，因此希望藉由蠕動 的方式所產生的對流效應加速過程以節省混合反應所需要的時間。 由於多數生物系統皆是藉由蠕動的方式來進行流體的輸送，故在醫學研 究上較為感興趣的即討論血液於血管內的流動所產生的剪應力與壓力等機 械性質對管壁的影響。又近來傾向應用在動脈硬化上，醫療裝置與血流的 相互作用上的設計與計算。在處理動脈硬化的外科手術中，藉由安置特殊 的支架使得因脂肪囤積而造成流道窄小的部位得到擴張。在這方面的計算 又分為假設血管為小量變形的剛性壁面與大量變形之壁面兩種，其中大量 變形之血管壁面在計算上較剛性壁面為複雜，因為此類計算通常都牽涉到

2 了流固耦合的問題。故在 Figuero et al. [14]中提到，巨觀的來說，要排除流 固耦合的因素的話則壁面的運動就不可太過劇烈，並且經由實驗觀察可得 知此並不會影響到速度場的計算。 1.2 文獻回顧： 在文獻回顧裡，主要為：

(1) Shaprio et al. [1]的以波座標(wave frame)及長波長與低雷諾數假設，對簡 化過的 N-S 方程式進行求解。

(2) Fung and Yih[8]的低振幅假設限制，將解建立在以振幅與入口高比值為 參數展開的冪級數近似解。

(3) Yi et al. [12]的封閉管道內的蠕動流場，其採用的求解方式是以 Fung and Yih [8]的低振幅限制與有限元素法的兩種分析方法。

(4) 其他研究生物心臟、血管與腸道內血液或流體的流動情況，以套裝軟體 或者有限元素法來求解 N-S 方程式。

Shaprio al.的波座標(wave frame)與低雷諾數假設：

Shaprio et al. [1]在蠕動流場乃建立在極大的波長與口比與忽略慣性力 (亦即低雷諾術)的假設之下求解 N-S 方程式，並依波動座標(wave frame)可 以將流場視為穩態來進行分析。其並於文中指出，在以波動座標(wave frame) 對進行分析所得到的結果可以將流場視為一純粹由壓力推動之 Poiseuille 流 動。

3

Misra et al. [2]在分析小血管中的二維蠕動流場中，以 Casson 模型與牛 頓流體的數學模型將流場分為核心與周邊兩部分，並在極大的波長與口比 與忽略慣性力的假設之下，以波座標(wave frame)來分析穩態的流場，並發 現當降低外層流體黏滯性時可促進整理流量的提升。

Mishra and Rao [3]在極長的波長(與入口高相比)與低雷諾數的假設下分 析二維不可壓縮蠕動黏性流場。同樣的亦使用波座標(wave frame)為主要參 考座標將整個流場視為穩態來分析，並探討在上下管壁以不同的振幅與相 角、頻率的傳遞波對壓力梯度、回流與流量等現象的影響。結果顯示回流 與陷流(trapping)層只有在管道截面具有差異時存在，且當管道越來越趨近 對稱管道的時候回流與陷流(trapping)區有增大的趨勢，並在對稱管道時達 極大值。

Hayat and Ali [4]在波座標(wave frame)之下以極長的波長(與入口高相 比)與低雷諾數的假設，並經過 order-of-magnitude 分析下簡化流場的統御方 程式，形成一 closed form 求得積分形式的解析解，並借助數值軟體 MATLAB 求解積分方程式型態的解析解。其中黏滯係數在此被假設為一依變空間變 數的函數並以數值解進行分析。結果顯示當流體黏滯性增高則管道內之壓 力梯度隨之增高。

Joel Jimenez-Lozano[5]以 Shaprio et al. [1]的理論討論二維、上下對稱蠕 動的牛頓流體流場，並求解流線函數，所得到的解為建立在長波長、低諾

4 數的近似解，並探討由對稱蠕動所產生的回流(back flow)以及陷流(trapping) 和增流(augmented flow)等現象，其中特別討論以上三種流動彼此之間的轉 換與穩定性問題。 Weinberg S. L et al.[6] 在實驗中驗證[1]的理論，在當雷諾數小於 1 的時 候成立，流動可視為無慣性力的(inertia-free)，並討論在不同的背壓下對平 均流量所造成的影響，同時藉由流體的可視化技術觀察由蠕動管壁造成的 回流(reflux)與陷流(trapping)。

M.V. Subba Reddya et al. [7] 研究在對稱管道中，由上下具有不停振幅 大小、波速以及相位差的蠕動波所造成的流動在管道的軸向速度以及對流 場的流量所造成的影響，其中所使用的流體為 Ostwald–deWaele power-law model 的流體。

Fung al.的低振幅假設限制：

Fung and Yih [8]相較於 Shaprio et al. [1]的研究，引入流線函數的定義求 解 N-S 方程式，且可以達到較高的雷諾數，不過在振幅上就沒有那麼自由 了，其與管道入口的比值必須非常的小。但是在雷諾數趨近於零以及相同 振幅之下的情況下，其結果皆與 Shaprio et al.[]的文章相符合。 Yin C. C. et al.[9]用自製實驗設備模擬只具有一側蠕動的流場(下側)，研 究不同的背壓、深度的管道在不同的振幅比下的流量，並將實驗結果與[8] 的理論解進行比較得到結果可知，[8]的理論在高振幅比的時候準確度極

5 低，誤差率甚至達到 10%，如果要在高振幅比得到更準確的理論解，則需 要更高階的數值解。 S. Nadeem * et al.[10] 藉 由 [1] 的 理 論 假 設 ， 利 用 小 擾 動 展 開 法 (perturbation expansion) 的 方 式 求 解 流 線 函 數 ， 其 中 [10] 的 流 體 乃 採 用 Williamson 為模型的非牛頓牛體來進行計算，並發現流體在高 Williamson 模型參數的時候，壓力降對流量的曲線呈非線性變化，且流體性質與一般 牛頓流體相似。

Mokhtar A. Abd Elnaby et al.[11]在研究中使用了具有彈性的壁面，並且 在流體與固體的界面上，同時考慮了流體與固體彼此之間的作用力對蠕動 壁面的變形所造成的影響，另外在流體力學的計算上，利用了擾動展開 (perturbation expansion)來求解流線方程式，Mokhtar A. Abd Elnaby et al. [11] 在研究中考慮了壁面的阻尼、彈性係數與張力等對流場的軸向速度等其他 方面的影響。 Yi et al.的封閉管道內的蠕動流場： 在 Yi et al.[12]的文章裡主要分析兩端封閉之蠕動流場，主要是將速度 與壓力場對振幅比展冪級數代入求解，並配合有限元素法分析流場。其流 場容器外壁依預先設定好的運動方式作動，使得容器內的流體進行攪拌。 文章指出此類流場適合應用在低雷諾數下的微型攪拌與促進化學反應。 其他生物方面相關的研究：

6

Melicher V. et al.[13]研究一維度的血管流動，文中提到了在大型動脈中 的血液流動可以用一維的流動模型來求解，主要的數值方法乃利用了有限 元素法的方法(explicit Taylor–Galerkin scheme)來求解。

Figuero C. A. et al. [14]在研究三維度的大型管動脈的血液流動中，考慮 了流體力學與固體力學之間的作用力，並實際求解了流力與固力兩者之間 的流固耦合的問題，並深入的探討在不同的壁面的邊界條件下所產生的現 象。 Heys J.J. et al. [15]在求解流固耦合的問題中，由於牽涉到了許多非線性 的方程式，因此導致了在計算上的不穩定，但是可以藉由多重網格，以及 有限元素法的方法(first –order system least squares)簡稱 FOSLS 法來獲得數 值計算上的最佳化的穩定性，因此 Heys J.J. et al. [15]利用 FOSLS 法來求解 三維的大型動脈裡的血液流動，同時並考慮進了血液與動脈壁面的流固耦 合問題。

1.3 研究內容：

從文獻回故可發現，不論是 Shprio et al. [1]的利用低雷諾數(Inertia free) 假設與長波長來對 N-S 方程式進行簡化求解，或者是 Fung and Yih [8]的將 解建立在以振幅與入口高比值為參數展開的冪級數近似解。這些方法不是

以特殊的條件來對 N-S 方程式進行簡化求解，就是以展級數的方式對振幅

7

et al. [12]的文章中有提出以 Fung and Yih [8]的假設和方法以及有限元素法 來直接對求解 N-S 方程式兩種方式。其中 Yi et al. [12]雖以有限元素的方式 可以直接求解 N-S 方程式，不過所解流場是一個兩端封閉的容器。 綜觀以上，並未有真正求解開放管道流的 N-S 方程式之問題，因此在 本研究裡，將嘗詴以有限體積法的方式直接求解開放管道流場的 N-S 方程 式。後面將討論在不同的管壁波動振幅與頻率下，管道內的壓力變化與整 體淨流量的表現。 本研究將以計算流體力學的方式來對蠕動式管道流場進行分析。流場外 壁已依預先設定好的運動模式來對流場進行擾動，造成管道內的壓力變化 達成流體輸送的目的。數值方法將採用非結構性網格，並以 ALE 座標系統 建立統御方程式。接著透過有限體積法對統御方程式離散，並以空間守恆 定理的滿足來達成質量守恆的要求，然後利用 PISO 演算法求解 N-S 方程 式。接著並討論各個參數之間的關係與對流量所造成的影響。 本研究主要的參數為背壓、振幅比、波數、雷諾數，以上四種參數各自 都會對流量產生影響，因此後面的分析將各別討論在固定其中一種參數後 改變另外三者的情況下會對流場造成何種影響。最後藉由本研究希望能得 在背壓極大的情況下，如何藉由調整振幅比或者波數來獲得較佳的流量的 結果。

8 第二章 數學分析模型 本章將介紹研究中所使用的物理模型與數學模型。在 2.1 節中將敘述流 場的物裡模型與基本假設。接著 2.2 節中，將描述在由前節流場的基本假 設為前提之下所成立的統御方程式，其中振動邊界的運動方程式也將於本 節中詳細介紹。2.3 節主要為流場的無因次化分析。最後 2.4 節則為流場的 邊界條件的描述與設定。 2.1 物理模型 流場如圖(2.1)所示為二維管道流，上端管壁將依照一時變函數的設置隨 時間產生一個連續不斷由右至左傳遞的波動。流場的進/出口位於管道的兩 端，藉由上端管壁的蠕動來驅動流體進出管道進而達到輸送的目的。其中 波動方程式的波長、波速分別為c、。 2.1.1 流場的基本假設條件： (1) 牛頓流體。 (2) 二維非穩態、不可壓縮。 (3) 層流。 (4) 忽略重力項的影響。 2.2 ALE 座標系統的數學模型 由 於 網格 將隨 時間 產 生變 形與 移動 的 關係 (移 動 網格 )， 故 在 使用 Eulerian 座標系統的時候對流場進行無法完盡的描述。因此必須採用 ALE

9 座標系統，而所謂的 ALE(Arbitrary Lagrangian-Eulerian)座標系統即當網格 移動速度與流體相同時，則可視為 Lagragian 的座標系統，若網格為固定不 動的時候，則可視為 Eulerian 的座標系統，下面將介紹在 ALE 座標系統描 述下的統御方程式。 2.2.1 連續方程式： 0 v   (2.1) 2.2.2 動量方程式：



g



 

v v v v p v t  _{ }  _ _{   } _    (2.2) 其中為流體密度，、p則為流體黏滯係數與壓力、vg為數值計算中網格 移動的速度。其中注意到上面的(2.1)式與 Eulerian 座標系統的連續方程式 相同，並不包含移動網格項，其理由可見下一章，在考慮空間守恆得到證 明。 2.2.3 移動邊界： 本流場於上端壁面的運動方程式可看做整數個隨著時間由左端入口連 續不斷的向右端出口以一固定波速進行傳遞的弦波，其數學式如下。

 





1 2 , 2 cos y x t h A  x ft     (2. 3) c f   為振動頻率(c為波速)，h為入口高度的一半長度，A為波動振幅，見 圖(2-1)。 2.3 無因次化： 取移動方程式中的波速(c)以及管道的一半高度(h)為特徵速度、特徵長

10 度進行無因次化，其中無因次化參數關係為： * x x h  ， * y y h  ， * v v c  ， * g g v v c  ， * 2 p p c   ， * ct t h  ， 1 * h    。 統御方程式： * * ( )v 0   (2.4)





 

* * * * * * * * * * * 1 Re g v v v v p v t  _{        }      (2.5) 其中無因次參數雷諾數(Re)定義為Re  c   。 移動邊界：









* * * * * 1 , 2 cos 2 y x t    x t _(2.6) A h   為無因次的振幅比，波數( )定義為  h  ，週期

 

T 則為 1  。 2.4 邊界條件 2.4.1 出/入口邊界條件 壓力邊界條件 給定出/入口壓力，在此定壓下討論整體進出/口的淨流量與流場內部壓 力變化，而進/出口速度的則由對流邊界的條件計算，數學式如下。 pconst (壓力邊界條件) (2.7) 0 c v v v t x  _  _   ，vc為對流速度。 (對流邊界條件) (2.8) 在第五章將介紹出/入口在對流邊界下的速度計算與對流邊界條件的離散 方法。 2.4.2 牆邊界：

11

採用無滑移邊界條件(no-slip)，即流體速度與壁面的運動速度相等，數 學式如下。

wall

12 第三章 空間的守恆 3.1 移動網格連續方程式 引用參考文獻[9]，在說明由移動網格造成的空間守恆問題之前先考慮 最初不加設任何條件的連續方程式：

 

0 div v t   _{ }  (3.1) 利用有限體積法將(3.1)式進行體積分可得：

 

0 d div v d t      _{   } _   





(3.2) (3.2)式經由高斯散度定理與 Leibniz 的微分法則上式可表式為：



g



0 d v v d s t     _{   } 





(3.3) 上式即為對一移動網格的積分型態連續方程式，積分中的v d sg 項是由於所 考慮的控容其體積非固定，其容積表面會隨時間改變位置，將本式進行離 散後可表示如下：

   

_

_

0 n o g f f v v s t       _{ }  _  _  



(3.4) 由於移動網格的緣故，所以連續方程式(3.4)式裡於對流項速度採用相對網 格移動的相對速度來計算。但不論(3.4)式採用隱式法(implicit scheme)或者 顯式法(explicit scheme)進行差分，其結果則如圖(3.1)與圖(3.2)所示，將產 生多餘或減少的質量使得計算上無法滿足質量守恆定律。這些並不滿足質 量守恆的值可視為人為的質量，或說質量源。 為了解決由於這些人為非保守的質量源，可以縮短在時間上的分割，如

13 此一來即可減少質量源的產生。不過縮短時距的作法無法完全的避免這方 面問題的產生，因此在這引入空間守恆的觀念。 3.2 空間守恆定理 空間守恆定理，其概念可由(3.5)式來表達： 0 g d d v d s dt



  



  (3.5) 上式所表示的即網格在隨時間變形時，其體積變化率將等同於控制體積上 的控制表面以其移動速度於單位時間內所掃過的體積。其離散化如下：

 

0 n o g f f v s t      



(3.6) 從上式可知，若在滿足空間守恆的情況下，網格速度所造成的體積流率可 由網格於前後時步下的體積差率所求得：

 

n o g f f f f v s t t    _{ }    





   (3.7) 單一表面積於單位時間內掃過的體積可表示為：

 

g f f v s t     _{ } _ (3.8) 亦即(3.4)式的對流項不論是以縮短時距、顯式法或者隱式法都勢將產生的 人為質量流率的問題將消失，如此一來滿足了質量守恆的條件。 則將(3.7)以及(3.8)式代入(3.4)式可獲得：





0 f f v s  



(3.9) 此即 (v)0的離散式，而在不可壓縮的條件下可進一步獲得： 0 v   (3.10)

14 上式即一般不可壓縮流的連續方程式，這也是在上一章中的(2.1)式，由此 式發現在滿足空間守恆定理的條件之下推導出的連續方程式，其移動網格 所產生的速度不復存在，在移動網格的問題中，藉由空間守恆的達成來替 代依網格移動速度來計算所造成的人為質量源。此結果最後導出了在同時 滿足空間守恆及質量守恆的條件下，連續方程式幾與一般的固定網格、不 可壓縮流動的問題相同，因此在接下來第五章的 PISO 演算法中，網格速度 將不會對壓力修正方程式造成任何的影響與改變，不過動量方程式中對流 項計算仍然必須採用相對速度獲得。

15 第四章 統御方程式的離散 在本章第一節的有限體法分析中，統御方程式被重新表示為積分方程的 型態，如此一來方程式較能表現出在流量的計算上的物理意義，這就是有 限體積法的優異之處，且利用本法可方便計算出通過小單元體積的面上各 種物理通量。以下將進行有限體積法的統御方程式的離散，並之後第五章 節中進行流場求解。 4.1 有限體積法及離散方程式 透過有限體積法對統御方程式進行體積分。本法的觀點主要是從一小單 元的控制體積上的各種物理通量的計算上出發，在獲得新的統御方程式後 再以各種不同的方法進行差分與離散。其中由於網格速度的關係，導致在 單一控制體積元素的流量的變化，這對接下統御方程式的離散與數值分析 中極為重要。 4.1.1 動量方程式 對動量方程式進行體積分得到下式：



g



d vd v v v ds v ds pd dt



  



  



   



 (4.1) (非穩態項) (對流項) (擴散項) (壓力項) ：控制體積積分範圍。  ：控制表面積積分範圍。 以下將仔細依動量方程式中的非穩態項、對流項、擴散項、壓力項來進行

16 離散。 4.1.2 非穩態項 在非穩態採用一階近似離散，如下所示：



n o



p p d vd v v dt t        



(4.2) 其中、t分別代表體積變量、時間變量，vnp、 o p v 則為新舊速度。 4.1.3 對流項





r g _f f f v v v ds m v   





c f f F 



(4.3) 上式中s代表面法向量、mr_f 為面上的質量流率、F_fc乃對流通量、下標 f 為控制體積上的表面。其中 r f m 定義如下： r f m  _f



g



_f f v v s   (4.4) 注意到上式為了保持空間守恆，故由網格速度

 

vg 所產生的相對流量將以第 三章中的體積變化率 (3.8)式來計算，可將其表式程如下形式： r f m  _f f m t     _{ } _ (4.5) 上式中的mf 代表不考慮網格移動時流過面上的絕對流率，其計算將於下一 章介紹。 對流項的差分上採用一階上風法與中央差分法的線性組合來表示：

 

n





o c UD CD UD f f f f F  F _ F F _ (4.6) 上式中 為介於 0~1 之間的因子，負責調控上風、中央差分法的比例， 為 1中央差分法、 為0則為上風法。而FfCD、 UD f F 則分別代表上風差分與中央

17 差分法的通量計算，其中等號右邊第二項為舊的時間下已知的計算結果， 故將其移至源項處理。上風與中央差分法表示如下：









max , 0 min , 0 UD p c f f f F  m v  m v (4.7)



1



CD p c f f p p F m _  f v  f v _ (4.8) 下標nb代表與位置p鄰近的位置， f_p為權函數。 4.1.4 擴散項 (4.2)式中的擴散項離散如下：





_f f v ds v s  







(4.9) 其中sf  d



sf d



，d為方向由主網格中心至鄰近網格中心的向量，其大小 d 可作為評估擴散作用的一個因子，並且對數值解的穩定性有著極為重要 的影響，其數學式表示如下： 2 f pnb pnb f s d s    (4.10) 其中pnb為由主網格中心點與鄰近網格中心點連線所構成之向量，如圖(4.1) 所示。則(4.10)式可以整理成如下所示：





D f f f f v s F  





(4.11)





f 2



n n

  

o





D nb p f _f f f pnb f s F v s v v v s d s            (4.12) (4.12)式中 D f F 為擴散通量。等號右邊第二項乃於前一個時間下的計算結 果，需移動至源項做處理。

18 4.1.5 壓力項 在壓力項裡面的壓力梯度皆為已知值，則離散如下：

 

pd p        



(4.13) 由於壓力梯度計算乃採取前一個時間下的計算值，故將本項移至源項處理。 4.2 線性代數方程式 將上述離散完畢的各項代回(4.2)式的動量方程式中可得到一線性代數 的方程式，且注意到離散的過程中，代數方程式裡具有特定的幾項因為採 用上一次時步(time step)下的計算數值給定，故這些項將是構成源項的重要 成員，則所要求解的代數方程式如下： p nb p nb f A v 



A v Q (4.14)





2 max , 0 f p f pnb f f s A m t s       _{ }        



(4.15)





2 min , 0 f nb f pnb f s A m s     (4.16) 其中





_{CD UD}



o

 

o







o

 

p f f f f f Q F F v s d v p t        _  _        



(4.17) (4.17)式即所謂的源項，由動量方程式離散後的已知項組成。其中由(4.14) 式所構成的矩陣即為數值解中所要求解的。

19 第五章 數值方法 由空間守恆定理排除了人為質量源的問題後，在第三章節中所離散的代 數方程式可以用矩陣的線性代數的聯立系統來表示。也就是說只需要求解 此矩陣所構成的線性代數方程式即可。而矩陣的求解方式有直接求解法與 數值上疊代的方式兩種，其中直接求解的方式在運算上相當的消耗電腦資 源與運算時間，故在此將以疊代的方式取代之。又現存主要的疊代方式與 數值方法有多種，這裡我們所要求解的系統乃聯立了動量方程式與連續方 程式，且在動量方程式中還有壓力項與速度項的耦合關係，本章將以數值 方法中的 PISO 演算法來處理壓力與速度項的耦合關係。當然了最後還是 必須要滿足連續方程式的質量守衡，下面將簡略介紹本問題所使用的 PISO 演算法與求解過程。 5.1 PISO 演算法 概念上，PISO 演算法的發展乃是為了求解非穩態流場問題。在處理壓 力及速度的耦合上是不需如其他演算法進行疊代的，因此速度上更加快 速。PISO 演算法與 SIMPLE 演算法相較之下其所求得的解乃暫態之並滿足 質量守衡，不過只需要將 PISO 於第二個修正步驟之下所得到的速度場、壓 力場重新代回設定為初始猜值繼續進行 PISO 的演算，如此一來反覆數次的 結果就會趨近於穩態問題的解了，因此 PISO 演算法在穩態問題上同樣適 用。結構上 PISO 演算法包含了一個預測步驟、兩個修正步驟，通常在第一

20 個修正步驟就可以滿足質量守衡的問題了，以下將詳細介紹 PISO 演算法的 運算步驟。 5.1.1 預測步驟 先利用現有壓力p由(4.13)式求出速度場v*：

 

* * p c p c p f A v 



A v   S p  (5.1) 不過由於此時的速度場 * v 並未滿足質量守恆定律，故需要在下列修正步驟 中使其滿足質量守衡。 將預測步驟(5.1)式移項整理可得：

 

* * p p _p p v H p A     (5.2) 其中 * 1 * p _c c f p H A v S A    _  _ 



，對每一控制體積要求質量守恆須先估算面上的流 速及流量，此面上流速可以下式求得：

 

 

* f f _{f f} p _f p _f v v p p A A  _  _    _ _  _ _   _{ }  _{ }   (5.3) 其中上標” ”表示面上的值是由其相鄰兩控制體積的中心點的值內插 求得。 其中pf 以及vf 為主格點與鄰近格點的內插值：



1



f p c p p p w p w p       (5.4)



1



f _{p c p p} v w v  w v (5.5) 上式中w_p為權重因子。

21 其中面上壓力梯度的係數值則是取平均值： 1 2 p _f p _p p _nb A A A   _ _ _             _    _   _    _ (5.6) 利用上述面上的速度預測步驟，可得面上的流率如下： * * f f f f m  v s





2 f f f c f f c p f c p _f f s v s p p p A s     _     _ _{ }    _   (5.7) 5.1.2 第一次修正步驟 本次修正所產生之速度場以符號 ** v 來表示之，且相對應之壓力場以符號 p* 來表示，其向量方程式如下：

 

** * _* p p _p p v H p A     (5.8) 將上式扣除預測步驟方程式(5.2)可得：

 

' _' p _p p v p A    (5.9) 同樣的，控制體積表面上的速度v'f 亦可表式為：

 

' _' f _f p _f v p A _ _ _    (5.10) 上式中， ' ** * p p p v v v 、p'p  p*ppp。 則第一次修正的質量流率可表示如下： ** * ' f f f m m m (5.11) 其中 _' ' f f f f m  v s 。 修正後的質量流率如下：

22









2 ** * f ' ' ' f f f f nb p f f pnb p _f f p _f s m m p p p s d A s A    _ _   _ _   _ _         (5.12) 修正後的面上質量流率滿足質量守恆定理： ** 0 f f m 



(5.13) 因此將(5.13)式代入上式可整理得： ' ' ' p p c c f a p 



a p Q (5.14) 2 f p f c f p f s a A s   _  _ _  



(5.15) 2 f nb f c p f s a A s   _  _ _   (5.16)





* ' ' _{f f f} f f f p _f Q m p s d A      _ _     





(5.17) 由 方 程 式 (5.15) 式 獲 得 第 一 次 的 壓 力 修 正 因 子 的 壓 力 場 ' p (pressure correction)。 5.1.3 第二次修正步驟 此步驟之速度v***滿足下列式子：





*** ** _** p p _p p _f v H p A _  _ _    (5.18) 將(5.19)式與(5.9)式相減可得：

 

" ' _" p p _p p _f v H p A _  _ _    (5.19) 其中， " *** ** p p p v v v 、p"p  p**p p*p。

23 流量修正表示如下：





_{ }

_

_

' 2 " " " " c f c f _{c p} f f f f f f f _f c p _f f p _f p A v s p p m p s d s A s A A        _ _     _{ }  _ _       _{ }      



(5.20) 第二次的流量修正如下： *** ** " f f f m m m (5.21) 上式滿足質量守恆定理，故可得： ** " 0 f f f f m  m 





(5.22) 由於 ** f f m



早已滿足了質量守恆，故 ** 0 f f m 



，因此上式可化為： " 0 f f m 



(5.23) 則由(5.21)式代入上式可整理為： " " " " " p p c c f a p 



a p Q (5.24) 2 " f p f c f p s a A s    



(5.25) 2 " f c f c p s a A s     (5.26) (5.24)式即第二次壓力修正方程式，則第二次修正後之壓力場 ** p 即求，同時 速度亦可由(5.18)式求得，其中等號右邊第二項 " Q (source term)表示如下：

 





' " 1 " c f f f c f f p f f p Q A v s p s d A A          _ _           







(5.28) 5.2 出/入口壓力邊界條件的流量計算 給定在鄰近邊界上的網格中心壓力等同於邊界上的壓力，接著在修正步

24 驟中令壓力修正p'為零，即令(5.14)式中源項與鄰近網格壓力項前面的係數 為零，如此一來邊界網格中心的壓力恆等於邊界上的壓力，接著利用所求 得之速度從連續方程式中求出出入口之流量，如圖(5.1)，推導如下： ' ' ' p p c c f a p 



a p Q 令a_c0、Q'0，則求得p'0，p_p  p_b 其中pb為邊界壓力。 接著再利用邊界的壓力值求解出的速度場計算流量。假設在二維的問題 中所使用的網格具有四個面，則連續方程式： 1 2 3 4 0 f f m m m m m 



(5.28) 其中假設m4為出/入口流量的話，則可表示為：





4 1 2 3 m   m m m (5.29) (5.29)式即為出/入口流量。 5.2.1 出/入口壓力邊界條件的速度計算 壓力邊界出/口端的速度則是透過對流邊界條件(2.8)式來計算： 0 c v v v t x  _  _   離散後的(2.9)式：



 



1 0 n o c n n b b b c v v v v v t   x     (5.31) 上標n、o與下標b、c分別代表新、舊時間的計算值與邊界和其相鄰的網 格，則方程式可以重新整理為：

25 1 o n n b c b v Crv v Cr    (5.32) 其中 vc t Cr x    為 Courant 數。 5.3 PISO 演算法流程 1. 給定初始值。 2. 計算新的網格，更新網格資料。 3. 由上一步的計算值求解動量方程式。 4. 第一次的修正步驟，求解壓力修正方程式得到p'壓力場，並得到滿足質 量守衡的速度場與壓力場。 5. 第二次的修正步驟，由第一次的修正值求解第二次壓力修正方程式得到 " p 壓力場，並對壓力及速度再做一次修正。 6. 重複第三到第五步的步驟，並進行演算直到滿足設定的時間。

26 第六章 結果與討論 本章首先將進行網格與管長對流量的測詴，確認管長與網格在接下來的 實驗中不會對結果造成太大的影響後，就可以對本數值方法的結果和[8]、 [9]的實驗數據以及理論值做比較，同時並驗證本數值方法的準確度。 在完成了網格與管長的測詴，並驗證了數值方法的準確度以後將針對各參 數之間的關係以及這些不同的參數設置將會對流量所造成的影響進行討 論，主要的參數為：波數

 

 、振幅比

 

 、雷諾數

 

Re 、背壓

 

p 6.1 網格測詴 在管道高

 

2h 為2、長度為32的網格測詴中使用了三種不同密度的網格 10 80 、20 160 、40 640 來進行測詴，其對應的 x 0.4、0.2、0.05， y 0.2、 0.1、0.05，其中時間間距

 

t 為103、波數

 

 為0.0907、振幅比

 

 為0.3、 管道左右兩端不設置壓差來測詴雷諾數1、雷諾數2的結果，測詴結果如表 (6-1)所示，此三種網格測詴的結果差距微小，因此為了節省計算上的時間 以及電腦資源的花費，在接下來的計算採用的網格皆相較本測詴的10 80 的 網格密度( x 0.4、 y 0.2)來計算。 6.2 到達週期性穩態的時間測詴 在背壓為零的情況下，測詴振幅比

 

 為1.0、雷諾數

 

Re 20、波數

 

 為 0.3125，出口平均流量與瞬時流量結果如圖(6-1)、(6-2)所示。從圖(6-1)、(6-2) 可知，在高振幅比、高雷諾數與高波數的情況，管道出口端的平均流量與

27 瞬時流量約莫在第四到第六週期內達到穩定表現，因此接下來在高波數的 計算皆採用第六週期的平均流量，因此接下來所有平均流量的計算皆為第 六週期的為準，其數據如表(6-2)。 6.3 不同管長的流量測詴 在完成網格測詴後，測詴四種不同管長L8、16、32、64、波數

 

 為 0.0907、振幅比

 

 為0.3在零壓差下不同雷諾數的流量結果，結果如表 (6-2)。由表(6-2)可知，四種不同長度的管道的平均流量彼此差距微小，不 過短管的流量較長管的多，不過差距是非常微小的，因此為了節省計算資 源與時間採用L32的管長來進行接下來的所有的計算。 6.4 數值方法驗證 使用10 80 的網格密度( x 0.4、 y 0.2)、管長L8、波數

 

 為0.0907、 振幅比

 

 為0.3、0.41，管道背壓為零的情況下，測詴不同的雷諾數的流量， 並與[9]的結果進行比較，結果如圖(6-3)所示。圖(6-3)包含了[8]、[9]的實驗 數據以及理論值，在[8]的理論解是根據小擾動展開法(perturbation expansion) 所求得，於圖(6-3)中以曲線(theory)表示，另外曲線(theory modified)乃由於 [9]考慮了在實驗中問題是三維的，因此在管道的深度方向上會造成流量損 失，故為了修正此問題，[9]在討論裡依管道深度制定阻力係數(resistance factor)，並用來修正[8]的理論解。在這裡主要的比較對象為[9]的實驗數據， 由圖(6-5)可知道本數值方法所得到的結果與實驗數據十分逼近，不論在高

28 振幅比還是低振幅比的情況皆如此，故可知本數值方法在準確度上是十分 可靠的，但是注意到本數值方法的計算結果為二維的，因此與實驗數據難 免存在著差距，此差距有可能是因為實驗設備的管道在深度上兩側管壁以 及上下壁面所造成的，故導致實驗結果的平均流量低於數值計算的結果。 6.5 振幅比的影響 (a)改變雷諾數： 固定波數

 

 為0.0907，管道背壓為零，測詴不同的振幅比

 

 在三種雷 諾數0.5、10、20下所得到的平均流量，結果如圖(6-4)所示，雷諾數的增減 對平均流量所造成的影響不大，僅只於振幅比微小的時候



 0.1



才可明顯 看出差距，不過隨著振幅比的增加，不同雷諾數的曲線開始接近，最後看 起來就像疊合在一起，由此可發現在高振幅比的時候雷諾數對流量的影響 微小，並可知平均流量整體趨勢是隨著振幅比而增加。 (b)改變背壓： 固定波數

 

 為0.0907，設置背壓在0.5、0、1.5的情況下，測詴不同的 振幅比

 

 在高雷諾數



Re20



下下所得到的平均流量，結果如圖(6-5)所 示，負值的背壓



0.5



表示入口壓力大於出口壓力，此種情況下壓力將幫助 流體的流動成為助力，這在小振幅比的時候可以明顯發現，管道下方的流 體在壓力的助益下毫無阻礙的往出口流動。而在高振幅比的時候，管道中 的流體因為管壁大幅度的變形，迅速提升管內流體於軸向流動的速度，使

29 得週期內流量獲得提升。反之在正值的背壓的情況下，因為壓力在流場裡 成為一種阻力，在流動的方向上阻止流體的流動，管道下方的流體開始往 入口流動，造成了週期內平均流量上的減少。不過在振幅比提升的情況下， 此情況會獲得改善。因此由圖(6-5)可知，負值的背壓所得到的平均流量大 於其他正值的背壓，而隨著背壓的增加，平均流量越來越少。另外當振幅 比越來越高的時候，由圖(6-5)可發現這三種不同背壓的曲線開始緩緩逼近 彼此，推測當振幅比



 2



，即蠕動壁面的振幅完全將管道封閉的情況下， 則流場的平均流量完全不受雷諾數以及背壓的影響了。 6.6 雷諾數的影響 (a) 改變波數： 固定振幅比

 

 為0.3、管道背壓為零，測詴與比較兩種波數

 

 0.0907、 0.3125在不同雷諾數下的結果，如圖(6-14)所示。從圖(6-14)可見在低雷諾數 的情況下，高波數的流量明顯高於低波數，波數越高則壁面振動頻率越高， 因此對低雷諾數流動來說可輸送的流量越高。 (b) 改變背壓： 固定波數

 

 為0.0907、振幅比

 

 為0.3，管道背壓為0.5、0.5、1.0、1.5 的情況下，測詴不同的雷諾數下的平均流量，結果如圖(6-15)所示，在負值 的背壓助力下，平均流量隨著雷諾數而增加，雷諾數越高則可看作黏滯性 越弱，因此導致流場的平均流量越大，反之當管道的背壓為正值的時候，

30 雷諾數越高則流體越容易受到壓力在流動的方向上的影響，導致流場因為 流體在流動方向上的減速而流量減少。另外注意到當背壓為零的時候，流 量幾乎不受雷諾數影響，此結果與[1]的理論解相符合，本研究的流場在依 照[1]的理論所推導出的流量方程式可見附錄 理論流量方程式。 6.7 波數的影響 (a)改變雷諾數： 固定振幅比

 

 為0.3，管道背壓為零，測詴不同的雷諾數

 

Re 0.5、10、 20 所得到的平均流量，結果如圖(6-16)所示，在波數小於 0.1 的時候，流量 曲線彼此差距不大，雷諾數此時對流量所造成的影響不大，然而隨著波數 的增加，雷諾數對流量的影響越來越大，並且雷諾數越高者其週期內的平 均流量越高，此乃由於高波數的管壁蠕動速度快，因此造成管內壓力變化 大，並加速了流體造成流量上的提升。 (b)改變背壓： 固定振幅比

 

 為0.3、雷諾數

 

Re 為 10，測詴背壓在-0.5、0.0、1.5 的 情況下，波數對流量的影響，結果如圖(6-21)所示，背壓為負值的時候壓力 會推動流體的流動，使得流量明顯大增，此乃由於流場內流體在壓力的助 益下，再加上管壁的蠕動使得流體獲得加速，進而提升流量。而高波數的 流場則因為壁面的高振動頻率使得流量高於低波數的流場。不過當背壓為 正值的時候，低波數的流場由於流體速度慢，即使管壁有在蠕動，但是因

31 為出口背壓的關係，流體在靠近出口的時候大量減速，使得整體流量減少， 並且注意到高波數的流場於管道下方的流體是朝入口方向流動，此現象將 造成平均流量的減少，雖然流量還是會隨著波數而增加，但是仍然會少於 負值背壓的流量。 6.8 背壓的影響 (a)改變波數： 固定振幅比

 

 為0.3、雷諾數

 

Re 10，測詴在不同的波數

 

 為0.0907、 0.1814，背壓對平均流量的影響，結果如圖(6-26)所示，當波數增加，則平 均流量的整體獲得提升，然而相較於波數較低的結果，波數的增加並未明 顯減緩流量的下降速率，兩條流量曲線彼此的斜率相當接近。 (a)改變振幅比： 固定波數

 

 為0.0907，測詴在不同的振幅比

 

 0.3、0.5下，背壓對平 均流量的影響，結果如圖(6-27)所示，當振幅比提高則流量曲線斜率會微微 的改變，乍看之下就像整條曲線往上平移，流量整體趨勢是隨著振幅比提 高而越大，但是隨著背壓的增加而減少。

32 第七章 結論 1. 本研究所採用的數值計算方法幾乎並不存在高雷諾數與高振幅比等因 素限制，並且計算出來的結果接近實驗數據。 2. 流量隨振幅比呈非線性增加，振幅比越大則流量也越大，同時在零背壓 與的情況下，不同的雷諾數流場在相同的振幅比越高時，彼此的差距越 來越小，流量曲線幾乎重疊再一起，且低雷諾數的流量曲線在平均流量 的成長速率上較高雷諾數為略快，因此雷諾數的影響對流量的影響亦可 忽略。 3. 在設置背壓後的平均流量曲線，在振幅比不大的時候可以明顯看出彼此 的差異，然而隨著振幅比的提升，不同背壓下的流量曲線還是會彼此互 相靠近。 4. 在零背壓的情況，流量幾乎不受雷諾數影響，此結果與[1]的理論流量 方程式相符合，並且在低雷諾數的情況下，流量會隨著波數增加。然而 當背壓改變後，整體流量則明顯受雷諾數所影響，並且因為背壓的提升 快速減少。 5. 流量隨著波數穩定的增加，並且在高雷諾數的流場，隨著波數越高流量 越高，此時雷諾數對流場的影響逐漸增強。 6. 平均流量整體表現隨著背壓呈線性遞減，並且當提高波數時，整體流量 亦線性的增加，但是流量隨背壓下降的速率並未明顯緩和。

33

7. 改變振幅比可以大幅線性提升流量，可是流量整體表現還是會隨著背壓 的提升遞減，同時流量的下降速率幾乎並未獲得改善，情況類似提升波 數。

34

參考文獻

[1] Shaprio A. H., et al., ”Peristaltic pumping with long wavelength at low Reynolds number ”, Journal of Fluid Mechanics , v37, pp.799-825,1969.

[2] Misra J. C., Pandey S. K., “Peristaltic Transport of Blood in Small Vessels： Study of a Mathematical Model”, Computers and Mathematics with Application, v43, pp.1183-1193, 2002.

[3] Mishra M., Rao A. R., ”Peristaltic transport of a Newtonian fluid in an asymmetric channel”, Z. angew. Math. Phys., v54, pp.532-550, 2003. [4] Hayat T., Ali N., ”Effect of variable viscosity on the peristaltic transport of a Newtonian fluid in an asymmetric channel”, Applied Mathematical Modelling, v32, pp.761-774, 2008.

[5] Lozano J. J., et al., ”Streamline topologies of two-dimensional Peristaltic flow and their bifurcations”, Chemical Engineering and P r o c e s s i n g : P r o c e s s i n t e n s i f i c a t i o n , v 1 6 , p p . 1 2 , 2 0 0 9 . [6] Weinberg S. L., et al., “An experimental study of peristaltic Pumping”,

Fluid Mech, v49, pp.461-497, 1971.

[7] Subba Reddya M. V., et Al., “Peristaltic motion of a power-lawfluid in an asymmetric channel”, International Journal of Non-Linear Mechanics, v42, pp.1153 – 1161, 2007.

[8] Fung Y. C., Yih C. S., “Peristaltic Transport”, Journal of Applied Mechanics, v35, pp.669-679,1968.

[9] Yin C. C., Fung Y. C., “Comparison of theory and experiment in peristaltic transport”, Fluid Mechanics, v47, pp.93-112, 1971.

[10] Nadeem S., Akram S., “Peristaltic flow of a Williamson fluid in anasymmetric channel”,CommunNonlinear Sci Numer Simulat,

v15, pp.47-50, 2009.

[11] Mokhtar A., et al., “A new model for study the effect of wall properties on peristaltic transport of a viscous fluid”, Communications in Nonlinear S c ien c e a nd Nume r i c a l Si mu l a tion , v13 , p p .752– 76 2 , 2008 . [12] Mingqiang Yi, et al., “Peristaltically induced motion in a closed cavity with two vibrating wall”, Physics of Fluids, v14, pp.184-197, 2004.

35

[13] Melicher V., Gajdosik V., “A numerical solution of a one-dimensional blood flow model-moving grid approach”, Journal of Computational and Applied Mathematics, v215, pp.512-520, 2008.

[14] Figuero C. A., et al., ”A coupled momentum method for modeling blood flow in three-dimensional deformable arteries”, Computational Methods of Applied Mechanical Engineering, v195, pp.5685-5706, 2006. [15] Heys J.J., et al., ”First-order system least-square(FOSLS) For

modeling blood flow”, Medical Engineering & Physics, v28, pp.495-503, 2006.

[16] Ferziger J. H., Peric M., Computational Methods for Fluid Dynamics. New York : Springer-Verlag, 1996.