7.2. 以特徵變數通量限制函數法求解
7.2.1. 非黏性流場之驗證測試
7.2.1.5. 斜震波流場解析－流經三角柱之高速流場
當一超音速自由流流經一向上折角
θ
時會產生如圖7-167 所示的斜震波 (摘自 Modern Compressible Flow [111]一書)，而此震波的物理型態及其斜震 波角度(β
)則會因自由流之馬赫數(M1)及壁面折角(deflection angle：θ
)而有所 差異，當還處於接觸斜震波(Attached Shock)的情況下(θ
＜θ
max)會滿足所謂的
θ
-β
-M 關係式，如下所示：2 2
1 2 1
sin 1 tan 2 cot
( cos 2 ) 2 M
M
θ β β
γ β
⎡ − ⎤
= ⎢⎣ + + ⎥⎦ (7.18)
其中
γ
=1.4。若已知一馬赫數M1，則存在有一個最大壁面折角
θ
max，如果物理幾何 外形折角大於最大壁面折角(即θ
＞θ
max)，則不會有直線的斜震波解，而此 情況下斜震波會產生彎曲變形且脫離壁面，即所謂的脫離震波(Detached Shock)，以下吾人就針對此一物理現象進行超音速流流經一正三角柱之流場 探討，以驗測本文發展之解子具有震波捕捉之能力。三角柱之迎流角為60 度，故其物理壁面折角(θ
)為 30 度，對於一固定折角時超音速自由流之馬赫 數(M1)由高遞減時，斜震波角度(β
)會由小遞增，直到形成接觸斜震波 (Attached Shock)的最小值 Mmin時，其所對映的θ
＝θ
max，當馬赫數持續減少 到小於Mmin時則斜震波會脫離三角柱而產生曲形之脫離震波(Detached Shock)，馬赫數愈小，其與三角柱前緣之距離則愈遠。由θ
-β
-M 關係曲線圖 [111]壁面折角(θ
)為 30 度時，可得知自由流之馬赫數(M1)在介於 2.4~2.6 之 間會產生折角恰為30 度之極大值θ
max，由Wikipedia 網站資料庫查得之θ
-β
-M 關係曲線圖如圖7-168 所示，可以較明確看出要使θ
max=30 時 M1 必須大於 2.5，另由式(7.18)可得到當自由流之馬赫數(M1)=2.5192 時θ
max =30 度，震波 角β
= 64.8 度。本項測試之計算網格使用7202 個外場三邊形網格及近場 3600 個四邊 形網格所組合而成的混合無結構性網格(如圖 7-169 所示)，採用 MD-1 特徵 變數通量Van Albada 限制函數法(dm=0.8)進行測試，測試的自由流馬赫數 (M1)由 2.6 開始進行並持續下降至 2.4(如圖 7-170~圖 7-175 所示)，由計算結 果顯示可以發現當馬赫數M1＝2.53 時(如圖 7-172 所示)斜震波已到達接觸 震波之極限，即三角柱之前緣頂點之馬赫數恰等於1，此與理論值 2.5192 相當接近，誤差只有在0.4%以內。M1＝2.52 時斜震波已開始脫離三角柱前 緣形成脫離型震波(如圖 7-173 所示)，馬赫數遞減時，其脫離現象愈明顯，
如M1＝2.4 時(如圖 7-175 所示)，測試結果完全滿足
θ
-β
-M 關係式，另在 M1＝2.53 時三角柱附近之流線及速度向量分佈如圖 7-176 所示，背部有形 成不規則的低速渦流現象。由測試結果顯示本流場解子除具有高速接觸斜 震波及脫離之曲線震波解析能力外，同時亦可處理極低速之背部渦流現 象，充份證明具有全速流之計算能力。7.2.2. 黏性流之驗證測試
7.2.2.1. 低速空穴流流場
本測試條件與前節原始變數求解法測試情況相同，以特徵變數通量限 制法(包括 MD1 及 MD2 兩種方法)，採用 SUPERBEE 及 Van Albada 通量限 制子(dm=1)進行四邊形網格及三邊形網格的計算，測試結果(a)流線圖；(b) 流場向量圖；(c)在中心線 x=0.5,y=0.5 之 v,u 速度分佈與 Chia 等人[115]之計 算結果比較圖；(d)計算殘值收歛情形：使用四邊形網格在 MD1 及使用 SUPERBEE 通量限制子之計算結果如圖 7-177 所示；Van Albada 通量限制 子之計算結果如圖 7-178 所示。在 MD2 及使用 SUPERBEE 通量限制子之 計算結果如圖7-179 所示；Van Albada 通量限制子之計算結果如圖 7-180 所 示。另外使用三邊形網格在 MD1 與 SUPERBEE 通量限制子之計算結果如 圖7-181 所示；在 MD2 與 Van Albada 通量限制子之計算結果如圖 7-182 所 示。沿著垂直中心線上之 u-速度及水平中心線上之 v-速度值分佈量化的比 較結果，其與Ghia 等人[115] 之計算結果相當一致且無數值振盪情形。
在本測試案例的計算結果顯示，不論是原始變數或特徵變數通量限制 法均能執行低速流場計算，唯在特徵變數限制法之 MD2 特徵變數求解方式 所得到的計算結果有不易收歛之現象，故在低速流時 MD1 特徵變數求解方 式優於MD2。
7.2.2.2. 流經圓柱之低速流外流場
以混合網格分別採用 Van Albada 及 SUPERBEE 限制函數來進行特徵 變數通量限制法之測試，參數設定為CFL = 1000、梯度修正次數 kgcor = 1(即 無修正)、壓力修正次數 kpcor = 3、鬆弛因子ηu=ηv=ηT=0.5；ηp=0.1，各種 條件下計算結果之(a)流線圖、(b)下游中心線之 u-速度分佈、(c)收歛殘值及 (d)圓柱表面之壓力分佈：圖 7-183 及圖 7-184 為特徵變數 MD1 使用 Van Albada 通量限制法在 Re=20 及 Re=40 之計算結果；圖 7-185 及圖 7-186 為 特徵變數MD2 使用 Van Albada 通量限制法在 Re=20 及 Re=40 之計算結果；
圖7-187 及圖 7-188 為特徵變數 MD1 使用 SUPERBEE 通量限制法在 Re=20 及Re=40 之計算結果；圖 7-189 及圖 7-190 為特徵變數 MD2 使用 SUPERBEE 通量限制法在Re=20 及 Re=40 之計算結果。
由計算結果與前節原始變數法比較顯示，原始變數法的收歛殘值在 10-3 即持平，而特徵變數法則可以再下降至 10-5以下；SUPERBEE 限制子所得 之渦流泡長度較Van Albada 限制子所得之渦流泡長度長。採用相同的限制 子在原始變數或特徵變數所得之結果差異不大。
表7-1 所示為本文計算結果與學者 Tritton(exp) [116] 、Takami 與 Keller [35] 、Nieuwstadt 與 Keller [117]、Coutanceau 與 Bouard [118]、Fornberg [119]
及Braza[120]等人所計算及實驗結果的量化數據比較情形，其中包括黏性阻 力係數CD、泡狀渦流終點位置 XE及流體分離角度 θsep，本文計算結果與參 考文獻之結果相當接近。
7.2.2.3. 雙喉部噴嘴內流場
計算網格同樣使用 315×40 四邊形無結構性網格，參數設定如下：
(1)特徵變數 MD1 限制法：CFL=1000；鬆弛因子ηu=ηv=ηp=0.5；ηT=0.8；
dm=2；壓力修正次數 kpcor=3；梯度修正次數 kgcor=2(2 表示修正 1 次)。
(2)特徵變數 MD2 限制法：CFL=1000；鬆弛因子ηu=ηv=ηp=0.5；ηT=0.8；
dm=2；壓力修正次數 kpcor=3；梯度修正次數 kgcor=1(1 表示無修正)。
各雷諾數的計算結果，其馬赫數、壓力及溫度等值圖，Reo=100 如圖 7-191 所示；Reo=400 如圖 7-192 所示；Reo=1600 如圖 7-193 所示。沿著壁 面之壓力分佈、摩擦係數(skin friction coefficient)分佈、沿著對稱中心面之 壓力與馬赫數分佈如圖7-194 至圖 7-197 所示，與原始變數通量限制函數法 之計算結果差異並不大且其與GAMM Workshop [113]的數據比較亦相當吻 合，顯示本解子對於黏性次音速、穿音速及超音速且與邊界層流混雜之複 雜流場同樣具有良好的解析能力。另數值收歛情形則在 MD1 時如圖 7-198 及MD2 時如圖 7-199 所示。
由原始變數及兩種特徵變數法進行 SUPERBEE 與 Van Albada 兩種限 制函數所計算之三種雷諾數結果，可以獲得流體分離點 Xs 及再接觸點 Xr 之位置如表7-2 所示，與 GAMM Workshop [113]之研究人員所得到的結果 相當一致。
7.2.2.4. NACA 0012 翼型外流場
與原始變數法同樣採用下述兩種條件來進行黏性流流經 NACA 0012 翼型之計算：(i) M∞=0.8, α=10°, Re∞=500；(ii) M∞=0.50, α=0°, Re∞=5000。
在第一個案例時，使用SUPERBEE 限制函數來進行測試，在特徵變數 MD1 計算時 dm=0.7；特徵變數 MD2 時 dm=1，計算結果如圖 7-200 及圖 7-201 所示，在翼型上表面同樣可以解析一個很大的渦流；其表面壓力 Cp 分佈情形亦與 Jawarhar 等人[75]用密度基底法之高階限制法之比較亦相 近；由收歛情形比較結果顯示特徵變數MD2 有較佳的收歛性。
第二種案例時，特徵變數MD1 時採用 Van Albada 限制函數 dm=0.8；
特徵變數MD2 時採用 SUPERBEE 限制函數 dm=2，計算結果如圖 7-202 及 圖7-203 所示，由流線圖可以看到在翼型尾緣處有輕微地不對稱渦流交互地 出現在 y=0 軸的兩側，其翼表面之壓力分佈則仍顯示具有對稱性。翼表面
之壓力分佈與文獻[75]之比較則均相當一致。