# 流體力學 Week 13 - Internal Flow

[TOC]

# 6. Internal Flow

首先是看看 Laminar Flow 與 Turblent Flow

* Laminar Flow : smooth streamline and highly ordered motion.
* velocity fluctuation and highly disordered motion.

這倆個主要是看「雷諾數」的大小。所以先回顧一下雷諾數：

$$Re = \frac {\rho V_{avg}D_{h}}{\mu}$$

這裡要注意的點是：這裡的 $D_{h}$ 是++Hydraulic Diameter++, 定義上是：

$$D_{h} = \frac {4A_{c}}{p}$$

其中 $A_{c}$ 是截面積， $p$ 流體截面積的周長。

比如說充滿流體的圓管：

$$D_{h} = \frac {4(\frac {\pi D ^ 2}{4})}{\pi D} = D$$

但是，如果只有一半有充滿液體的話，這個會變成：

$$D_{h} = \frac {4(\frac {\pi D ^ 2}{4}\cdot \frac {1}{2})}{\frac {\pi D}{2} + D} = \frac {\pi}{\pi + 2}D$$

之前說過雷諾數是黏滯力和慣性力的比例，黏滯力會讓流體粒子之間不會亂跑，而慣性力會讓流體亂滑。所以如果雷諾數很大，那可以預期慣性力相對黏滯力很大，粒子之間容易亂滑，就容易有 turblent flow 。主要的分界是：

* Laminar Flow : Re < 2300
* Transitional Flow : 2300 ~ 4000
* Turblent Flow : 4000 +

雖然是這樣說，在很安定沒有擾動的狀況之下，實驗可以做到雷諾數 100000 還保持 laminar flow 。但這只有在實驗室控制的很小心的狀況下才做得到。

# 6.2 Entrance Region

流體剛流進管子裡面的時候，跟流體流經很長一段距離，變成 fully-developed 的時候，兩個表現出來的速度分佈樣子會不太一樣。

而流體進去管子，到變成 fully developed 中間的這段就叫做「entrance region」。

在逐漸變成 fully-developed 的過程大概是這樣的：


1. 進去之前是 uniform flow
2. 進去之後，流體就被邊界給拖累住，所以++管壁附近的流體速度變慢++。不過因為質量要守恆，所以如果周圍變慢，中間就要會變快一點點。所以靠管壁的地方就比平均慢一點點，中間一點的地方就比平均快一點點，但大致上還是平的。
3. 隨著流體持續移動，邊界的影響持續擴大，由外往內深入，最後就影響到最中間的部分。

邊界的影響是用速度分佈來看的。「邊界還沒影響到的地方」意思是「速度分佈還是平的區段」。

如果把整個過程畫出來的話，大概就是這樣：


「邊界的影響」會對稱地逐漸由管壁往管子中心深入，最後影響到中心。變成 fully-developed。這段 fully-developed 之前的這段區域就叫做「Entrance Regeon」，個逐漸發展的行為就叫做「developing flow」; 而最後面，已經完全影響到中間以後的部分，叫做「fully-developed」。

那這個入口區的長度是多少呢？這裡有兩個公式：

* Laminar Flow :
$$\frac {L_{laminar}}{D} \approx 0.05Re$$

* turblent Flow :
$$\frac {L_{turblent}}{D} \approx 1.359Re^{\frac {1}{4}}$$
如果 $Re$ 用 4000 代的話，大概會是 10 左右。所以也有用 10 直接去近似的狀況。

至於這個公式怎麼來的？課本會跟你說是某一篇 paper 上做出來的東西。

# 6.3 Head Loss : Major and Minor Lost

* Major Lost： 摩擦力造成的損失。會這樣叫是因為如果這是直的，那唯一會影響的就是摩擦力。
* Minor Lost : 其他部分造成的損失，比如說流過閥，或是管道有轉彎、縮口的部分，都會讓外部是給流體力，然後減少能量。

這裡要注意的是，Minor Lost 不一定會是 Minor, 也有可能 Minor Lost 比 Major Lost 還要多。

# 6.4 Major Lost of an incompressible fully developed laminar flow

## Recall : Poisuelle Flow

這裡是用 Poiseulle Flow 為模型。Poiseulle Flow 的速度分佈：

$$\vec{u{}} = u(r)\hat{z{}}$$

(如果已經知道他是一條拋物線的話，那這段可以跳過去)。

然後塞進 Momentum Equation 裡面：

$$0 = -\frac {dP(x)}{dz} + \frac {\mu}{r}\frac {d}{dr}(\frac {rdu(r)}{dr})$$

因為一個是 z 的函數，一個是 r 的函數，所以：

$$\frac {dP(z)}{dz} = \frac {\mu}{r}\frac {d}{dr}(\frac {rdu(r)}{dr}) = const.$$

首先看：

$$\frac {-dP(z)}{dz} = const.$$

所以壓力分佈是線性的，所以

$$\frac {-dP(z)}{dz} = const. = \frac {\Delta P}{L}$$

$\Delta P = P_{upstream} - P_{downstream} = 壓力差$

接著看右手邊：

$$\frac {\mu}{r}\frac {d}{dr}(\frac {rdu(r)}{dr}) = const = \frac {\Delta P}{L}$$

邊界條件是 $u(0)$ 有限，$u(r) = 0$。

解出來就是：

$$u(r) = \frac {(\Delta P / L)R^2}{4 \mu}(1 - \frac {r^2}{R^2})$$

(其實就是之前解的那個)。

然後就可以開始算一堆有的沒的東西了～～

## P 的力量

因為 Poissuelle Flow 一句話版本就是「壓力差推動流體移動」，所以就會忍不住想猜所有東西都是壓力的函數。實際上也差不多是。

1. 流速：

$$u_{max} = 2u_{avg} = \frac {\frac {\Delta P}{L}R^2}{4\mu} = \frac {(\Delta P/L)D^2}{16}$$

這裡可以注意的是平均流速是最大流速的一半，所以可以在積分時偷懶一下。

2. Volume Flow Rate：

$$Q = \frac {1}{\pi R^2}\int_{0}^{R}u(r)2 \pi r dr = \pi R^2 = \frac {(\Delta P / L)\pi R ^ 4}{8 \mu} = \pi R^2 = \frac {(\Delta P / L)\pi D ^ 4}{128 \mu}$$

3. Shear Stress 也可以算，利用：

$$\tau_{zr}(r) = \mu\frac {du(r)}{dr} = \frac {(-\frac {\Delta P}{L})r}{2}$$

因為牆壁受到的 shear stress 跟邊界方向相反，所以：

$$\tau_{rz}(r) = \frac {(\frac {\Delta P}{L})r}{2} = \frac {(\frac {\Delta P}{L})D}{4}$$

<br>

算到這裡就可以發現，上面 3 個東西都是跟壓力或壓力梯度有關的函數。所以整理一下：

$$\frac {\Delta P}{L} = \frac {32\mu u_{avg}}{D^2} = \frac {4\tau_{w}}{D} = \frac {128 \mu Q}{\pi D ^ 4}$$

只要有壓力梯度就可以算出所有的東西。壓力差跟平均流速都有的話：

$$\Delta P = (\frac {64\mu}{\rho u_{avg}D})\frac {L}{D}\frac {\rho u_{avg}^2}{2} = \frac {8\tau_{w}}{\rho u_{avg}^2}\frac {L}{D}\frac {\rho u_{avg}^2}{2} = \frac {128 \mu Q L}{\pi D ^ 4}$$

然後你可能會想沒為什麼不把上面那幾坨東西化簡一下？再仔細看一下可以發現沒化簡的部分最後面黏的項都是動壓。這就在暗示接下來要做的事了。

定義 Darcy Friction Factor :

$$f = \frac {8\tau_{wall}}{\rho u_{avg}^2}$$

這個東西（看上面的算式）可以知道在 Poisseulle Flow 會等於：

$$f_{laminar} = \frac {64}{Re}$$

至於其他的狀況，就自己去查表拔～～

那這個東西訂出來有什麼用處呢？如果已經有 $f$ 的話，那麼壓力差就可以這樣算：

$$\Delta P = f\frac {L}{D}\frac {\rho u_{avg}^2}{2}$$

只要給定管壁的動壓，給定管壁的尺寸，加上查出相應的 $f$ ，就可以知道壓力減少多少。而且還是只要給平均流速就算得出來，實驗上平均流速也很好量。真是太親民了。

然後可以更進一步定一一個頭出來：

$$H_{laminar, poisseulle, cylinder} = \frac {\Delta P}{\rho g} = f\frac {L}{D}\frac {u_{avg}^2}{2g}$$

這個東西只要速度頭知道，$f$ 查表查到，就可以算出來了。而最有趣的是維持這個流要給的功率：

$$\dot{W} = Q\Delta P = \dot{m}\rho g \Delta H$$

就可以把 head loss 當成高度差去算能量問題。

# Major Loses for Turblent Pipe Flow

* turblent : Highly disordered flow

可能會有一個主要的震盪，不過不同方向會有 fluctuation :+1:


通常的做法就是用一個平均與他的 fluctuation 的相加。不過這個東西很複雜，所以這裡就暫時先不討論。

那在館內流當中， turblent flow 的結構會是什麼樣子呢？這個東西可以分成幾層：

1. 壁面附近：叫 viscous sublayer 。很亂~~所以沒人知道怎麼算他~~。
2. 接著是以外的區域。這裡平均速度的分佈大致滿足 power law :
$$\frac {\bar{u}}{U} = (1 - \frac {r}{R})^\frac {1}{n}$$

一般來說，雷諾數越大, n 越大，圖形的鼻部會越緩。

因為 turblent flow 實在是太亂了，關係很難直接用推的。所以這裡用 dimensional analysis 取而代之。假定壓力差會跟管子的尺寸( $D$ )、流體流經的距離( $L$ )、表面粗糙度( $\epsilon$ )、流體平均速度($\bar{u}_{avg}$)還有流體性質($\rho$、$\mu$)有關，也就是：

$$\Delta P = \Delta P(D, L, \epsilon, \bar{u}_{avg}, \rho, \mu)$$

然後開始無因次化。不過因為無因次化的過程有點無聊，所以就省略掉。無因次化完的結果就是：

$$\frac {\Delta P}{\frac {1}{2}\rho \bar{u}_{avg}^2} = f(\frac {L}{D}, \frac {\epsilon}{D}, Re)$$

因為實驗做出來大致上來說，「$\Delta P$ 正比於 $L$」, 所以可以把這個關係更近一步寫成：

$$\frac {\Delta P}{\frac {1}{2}\rho \bar{u}_{avg}^2} =f(\frac {\epsilon}{D}, Re)\frac {L}{D}$$

因此：

$$\Delta P = f(\frac {\epsilon}{D}, Re)\frac {L}{D}\frac {1}{2}\rho \bar{u}_{avg}^2$$

然後就發現這東西跟剛剛訂出 Darcy Friction Factor 時，形式根本是一樣的！只是剛剛是說：如果知道是 Laminar Flow ，那 $f$ 可以直接用 $\frac {64}{Re}$ 代掉, 而 Turblent Flow 是...恩...不知道是什麼東西。

那這個 $f$ 要怎麼找？這要使用另一個~~沒有為什麼~~從別人 paper 引用來的公式。

## Colebrook Equation

大致上就是在跟你說那個 $f$ 跟 $Re$, $\epsilon$, $D$ 的關係。

$$\frac {1}{\sqrt{f}} -2.0log_{_{{10}}}(\frac {\epsilon / D}{3.7} + \frac {2.51}{Re\sqrt {f}})$$

這是一個隱函數的關係，所以不能直接把 $f$ 帶出來。不過可以用計算機迭代來找出 f 。

有一個比較近似的顯示關係是：

$$\frac {1}{\sqrt{f}} \approx -1.8log_{_{{10}}}(\frac {6.9}{Re} + (\frac {\epsilon / D}{3.7})^{1.1})$$

如果表面很光滑，也就是 $\epsilon = 0$，這個就可以做一些化簡。變成：

$$\frac {1}{\sqrt{f}} = 2.0log_{_{{10}}}(Re\sqrt{f})$$

這個叫做 Prantl Equation。因為是 Prantl 做出來的。

近似的顯式表示就是：

$$\frac {1}{\sqrt{f}}\approx 1.8 log_{_{{10}}}(\frac {2.51}{Re\sqrt{f}})$$

另外一個特例就是 $Re \to \infty$，那麼就換另外一項變成 0。所以就變成：

$$\frac {1}{\sqrt{f}} = 1.8log_{_{{10}}}(\frac {\epsilon / D}{3.7})$$

這個叫做 Von Karbon equation。

>> 蔡：這樣看起來 Prantl 跟 Kaman 都在偷懶...不過實際上是 Prantl 跟 Kaman 先做出一部分的結果，最後 Colebrook 把他們兩個整合再一起的...所以這個啟示就是「適時」的「統合」前人的結果...

上面那一堆醜的要命的公式真的是很難帶。所以有一個叫 Moody 的人做了更懶的事：


這個東西就叫做 Moody Chart

這個圖的概念是這樣：
1. 每一條曲線代表一個 $\frac {\epsilon}{D}$ 值
2. 每一個曲線經過一個界限之後，都會變平的。把那些準備轉平的點通通找出來，連成一條線。就是上面那個切開的虛線
3. 左邊那一坨意義不明的東西是 Transition Flow 的意思。

講了這麼多，這些東西到底要怎麼用呢？通常處理問題的方法是這樣的：

1. 從 $Q$ 來找 $\bar{u}_{avg}$
2. 算雷諾數，看是哪一種流。
3. 找出對應的東西代：
* 如果是 Laminar, 就帶 $\frac {64}{R}$
* 如果是 turblent flow, 就手動迭代出 $f$。不然就查表(如果有表可以用的話)。
4. 把對應的 f 帶進 $\Delta P = \frac {L}{D}\frac {\rho \bar{u}_{avg}^2}{2}$

所以可以發現關鍵都在查表、代公式、迭代。

迭代的話，先用顯示猜一個值，然後一直輪流迭代。大概的演算法像這樣：

    #define ERROR 10e-5 // 通常是有效位數 4 位以下
    double pre = explicit(Re(u, L, D, mu), roughness, D, L);
    double f = implicit(Re(u, L, D, mu), roughness, D, L, init);
    while(abs(pre - f) > ERROR) {
        double tmp = Re(u, L, D, mu);
        if(tmp < 4000){
            printf("exception occur.\n");
            break;
        }
        f = implicit(Re(u, L, D, mu), roughness, D, L, f);
    }

這個東西的基本原理其實就是 fixed-point method。

Major Lost 主要都是在教怎麼迭代。這...有點尷尬啊。

順帶一提，如果你計算機是用~~工院神器~~ Casio 991 的話，他有內建[解方程式](https://georgegarside.com/blog/casio-calculator-tutorials/fx-991es-plus/equation-solve/)的功能。可以減輕迭代按到手抽筋的症狀。

# 6.7 Minor Losses

他可能比 majot lost 還要大。

Minor Lost 就有點尷尬，因為狀況很多種，所以不一定導的出來。所以作法就是定義一個常數然後把所有因素都推給他。這個倒霉鬼叫做「Lost Coefficient」：
$$k_{L} = \frac {H_{L}}{\frac {u_{avg}^2}{2g}} \Rightarrow H_{L} = k_{L}\frac {u_{avg}^2}{2g}$$

或者是用等效長度。這種想法是「轉彎造成的 head lost, 看起來就像長度增長了多少」。所以就是：

$$H_{L} = k_{L}\frac {u_{avg}^2}{2g} = f\frac {L_{eq}}{D}\frac {u_{avg}^2}{2g}\Rightarrow L_{eq} = \frac {D}{f}k_{L}$$

詳細可以看表 8-4


不過這種做法僅限於「沒有分流」的狀況。

這裡有一些注意事項：有導角 vs. 沒導角的 $k_{L}$ 會差很多，理由是因為這樣：


沒有導角的狀況就會產生兩個 Recirculating bubble, 流體流進去就會一直空轉跑不出來(而且還會一直白白消耗能量)，這個效果就好像管壁變小了，這個效應叫做「Vena Contracta」，中文叫「頸縮現象」。

另外一個可能讓 head 減少很多的狀況較「sudden expansion」，像這樣：


這個也會造成~~送頭~~頭損失。這裡的 $k_{L}$ 是：

$$k_{L} = \alpha (1 - \frac {A_{small}}{A_{large}})^2$$

其中 $\alpha$ 叫做「Kinetic Energy Correction Factor」，這個跟之前的 $\beta$ 很像，只是變成能量版：

$$\alpha = \frac {\dot{kE_{tot}}}{\dot{kE_{avg}}} = \frac {\frac {1}{2}\rho \int u^3 dA}{\frac {1}{2}\rho u_{avg}^3A} = \frac {1}{A}\int(\frac {u}{u_{avg}})^3 dA$$

這個東西的值：

$$\alpha = \begin{cases}
2 & for\ fully-developed\ laminar\ flow \\
1.05 & for\ fully-developed\ turblent\ flow
\end{cases}$$

# 6.8 total head loss and pressure drop

這裡要介紹怎麼結合剛剛學到的東西做一堆 combo 技能。

首先是「如果有很多個部分都會造成頭損失，要怎麼算？」，這裡的作法其實滿直覺的，就是全部把她們加起來：

$$H_{lost} = H_{major} + H_{minor} = \sum f_{i}\frac {L_{i}}{D_{i}}\frac {u_{avg,i}^2}{2g} + \sum k_{L,j}\frac {u_{avg,j}^2}{2g}$$

然後注意一點：這個頭損失是考慮到「整個截面」之下的平均頭損失。不過印象中如果伯努力方程要長這樣的話：

$$H_{tot} = \frac {p}{\rho g} + \frac {|\vec{u{}}|^2}{2g} + z$$

必須要「沿同一條流線」，也就是這個頭是單一流線的 case。所以這個也要平均起來才行。所以這時候頭要用平均的頭「average elevation head」：

$$H_{tot,avg} = \frac {1}{\dot{m}}\int H_{tot}(\rho u)dA = \frac {1}{\dot{m}}\int (\frac {p}{\rho g} + \frac {|\vec{u{}}|^2}{2g} + z)(\rho u)dA$$

(就是仿照平均速度、平均動能的方式下去平均)。

然後想當然就是做一堆不乾淨的積分，然後用：

$$\alpha = \frac {1}{A}\int (\frac {u}{u_{avg}})^3dA$$

把動能那一項的積分代掉。所以這時候的頭就要用下面這個東西了：

$$H_{tot, avg} = \frac {p}{\rho g} + \alpha \frac {u_{avg}^2}{2g} + z_{avg}$$

所以這時候，就可以寫成：

$$H_{upstream} = H_{downsream} + H_{lost}$$

也就是：

$$\frac {p}{\rho g} + \alpha \frac {u_{avg}^2}{2g} + z_{avg}\bigg{|}_{upstream} = \frac {p}{\rho g} + \alpha \frac {u_{avg}^2}{2g} + z_{avg}\bigg{|}_{downstream} + H_{lost}$$

反正後面就是自己把他帶開就對了。