中 華 大 學

中 華 大 學 碩 士 論 文

題目：高升力尾翼氣動力初步分析

系 所 別：機械工程研究所 學號姓名：M09608011 許 廷 旭

指導教授：牛 仰 堯 博 士

中華民國 九十八 年 七 月

致 謝

本文之完成首先要感謝指導教授 牛仰堯 博士在這段期間 的督促與指導，再者要感謝台大應力所 張建成 教授、朱錦州 教 授以及謝政達 博士，在Fluent及力源理論的指導，在此致上十二 萬分的敬意。並感謝胡舉軍 博士、張智豪 博士以及 丁大為 博 士 等口試委員對於本拙著之指正。同時感謝國科會計畫[NSC 97-2221-E-216-009-MY2]之經費補助與國家網路與計算中心電 腦設備協助。

此外，也要一併感謝我的家人與親朋好友等在研究所兩年來 所給予的支持與鼓勵，還有實驗室邱勇憲 學長、湯漢威 與 黃 君萍兩位學弟妹，以及我的同窗好友等和學長們的幫忙和協助，

使我的論文可以順利的完成。

摘要

本文的研究目的討論機翼氣動力分析，了解流場變化;以現行 CFD 技術，結合力源理論進行高升力氣動力分析，並分析流場受力 情形，了解高升力機制。分三部份討論，第一部份，以文獻得知單方 程 Splart-Allmaras(S-A)紊流模式的模擬有不錯的成效，與 Fluent 多方 程 Reynolds Stress model(RSM)比較；第二部分，利用力源理論分析 動態失速受力情形；第三部份，分析加入 Gurney flap 裝置改變的升 力與阻力的變化與受力情形，並了解分離出的流場，實際貢獻狀態，

使其有利於應用。

Abstract

In this study, the purpose of this paper to discuss the research wing aerodynamic analysis, to understand changes in the flow field; to the current CFD technology, combined with rising power source of the theory of dynamic analysis efforts and to analyze the flow field where the force to understand the mechanism of rise in force.

Discussion is divided into three sections, the first part to the literature that the single-equation Splart-Allmaras (SA) turbulence model has good simulation results, with more than Fluent equation Reynolds Stress model (RSM) compare with each other; the second part, the use of

Force Element to analysis the dynamic stall force situation; third part of

the analysis by adding devices to change the Gurney flap for lift and drag forces of change and the situation, and learn to isolate the flowfield, the actual contribution to the state, so that it is conducive to the application.

目 錄

誌 謝 ………. I 中文摘要 ………. II 英文摘要 ………. III 目 錄 ………. IV 圖表目錄 ………. V

第一章 緒 論 ……….……….. 1

1.1 前 言 ……….. 1

1.2 文獻回顧 ……… 4

第二章 數值模式 ………. 9

2.1 統御方程式 ……….... 9

2.2 數值方法 ……… 9

2.3 力源理論 ……… 12

第三章 數值結果 ……….… 15

3.1 紊流模式之比較 ………..……… 15

3.2 藉由力源理論在馬赫數(M)=0.1 的動態失速現象 解析 ………..……….. 29 3.3 Gurney flap 流場現象解析 ……….. 40

第四章 總結………. 58

參考文獻 ………. 60

圖與表目錄

圖目錄

圖 1. NACA0012 幾何外型整體網格圖 19 圖 2. NACA0012 局部網格放大圖 19 圖 3. (a) 攻角與升力係數變化圖(實驗值)

(b) 攻角與升力係數變化圖

20 圖 4. (a) 攻角與阻力係數變化圖(實驗值)

(b) 攻角與阻力係數變化圖

21 圖 5. (a) 攻角與力矩係數變化圖(實驗值)

(b) 攻角與力矩係數變化圖

22 圖 6.

S − A

週期壓力係數等值圖與機翼表面壓力係數分

佈圖(a)-(d)

23-24 圖 7.

RSM

週期壓力係數等值圖與機翼表面壓力係數分

佈圖(a)-(d)

25-26 圖 8.

S − A

週期流線圖(a)-(d) 27 圖 9.

RSM

週期流線分佈圖(a)-(d) 28 圖 10. (a) 表示升力元素週期變化圖

(b) 攻角與升力係數的變化圖

33 圖 11. (a)表示阻力元素週期變化圖

(b)攻角與阻力係數的變化圖

34 圖 12. 圖(12) 週期中升力(左)與阻力(右)元素等值線圖

(a)-(h)

35-37 圖 13. 週期中流線圖(a)-(h) 38-39 圖 14. NACA4412 幾何外型整體網格圖 46 圖 15. NACA4412 局部放大網格圖 46 圖 16. NACA4412 加 Gurney flap 局部網格放大圖 47 圖 17. NACA4412 表面壓力分布圖 47 圖 18. 圖(18) NACA4412 與 NACA4412 加 Gurney flap 機翼

表面壓力分佈圖

48 圖 19. 機翼表面壓力分佈等值線圖(a)-(g) 49-50 圖 20. 升力係數與 Gurney flap 尺寸變化圖 51 圖 21. 阻力係數與 Gurney flap 尺寸變化圖 51

圖 22. 升力阻力比值與 Gurney flap 尺寸變化圖 52 圖 23. 升力元素(左)與阻力元素(右)等值線圖(a)-(g) 53-55 圖 24. NACA4412 加 Guenwy flap 流線圖(a)-(g) 56-57

圖目錄

表 1. NACA4412 不同尺寸 Gurney flap 之升力元素貢獻 45 表 2. NACA4412 不同尺寸 Gurney flap 之阻力元素貢獻 45

第一章 緒論

1.1 前言

人類開始有飛行器設計構想的人是李奧納多‧達文西(Leonardo

Da Vinci)。他本人並沒有接受正統的教育培養，憑藉著對大自然的熱 情，透過觀察、自修以及推敲研究一點一滴累積而來的知識，進而呈 現他對科學的認知。達文西一生熱衷思考大自然現象，加上他有一雙 俱繪畫天份的巧手，把分析研究所觀察的東西很有系統的紀錄，他留 下 100 多張關於飛行的畫作可知，最初飛行的夢想，希望可以向鳥類 一般的飛行，但是人類的體力沒辦法勝任，進而言就固定翼的氣體動 力機制，也就是現今的滑翔翼，當時由他弟子進行飛行測試，結果失 敗，他的弟子因此失去一條腿。雖然當時沒辦法使人類有辦法在空中 飛翔，但是達文西留下的寶貴資料，給後繼者無限的啟發。

之後的數百年，飛行技術的構想不斷精進改良。最初的開端由一 位德國人奧特‧李連塔爾(Otto Lilienthal)，他熱衷於飛行不斷試飛改 良，在某次試飛遇到強風落地重傷，次日死亡，當消息傳到萊特兄弟 決定繼承他的遺志。研究李連塔爾留下的飛行和設計數據，不斷精進 改良，最終研發出可自由控制的滑翔機，雖然沒有飛行很遠，但是為 飛行開啟新的里程。

一般流體在低速(馬赫數<0.3)時的運動中被視為不可壓縮流 (Incompressible flow)。對於不可壓縮流，在工程與物力的應用上具有 很大的重要性。利用離散數值分法來求解不可壓縮流體方程式已獲得 完整的流場解析。為獲得真實流場的計算結果，必須模擬極為複雜的 幾何外型與物理現象。而其中最主要的關鍵在於網格產的好壞。在計 算流體力學領域中，產生好的網格，才能獲得準確的答案，所以這部 份的處理，占有不可輕忽的重要性。

流體分離(Flow separation)現象，在流體力學上是一個很重要的課 題，且廣泛地出現許多工程上的應用，例如:空中的飛行載具、陸上 的汽車或水中的船舶、潛艇等與流體有相關的運動，皆有流體的分離 現象。一般而言，在接近最佳操作的情況下才會有流體分離產生，例 如:當飛機達到最大升力時，就會發生此現象，在此時阻力會增加，

升力會下降。如果能有效且適當的控制，則會使性能和效率大幅增 加，在現今需要較大的運輸量的情況，可以獲得更大的商業性質，將 有助於減少運輸成本。

一般為了防止與減少流體分離的方法有，例如:使用導葉或不連 續面、除去或增加固體空間…等等，但是這些都是被動方式。須要高 性能的操作要求，就需要主動式的控制方式:開啟隙縫以吹出或者吸 入流體，或利用襟翼改變機翼的幾何外型。

在機翼的流場中，當超過一個臨界攻角的時候，會有一個導致升 力驟降的情形，這個現象稱為失速(stall)。當失速時，一定伴隨流體 分離的現象，於此同時會產生大尺度的渦漩。吾人相信此渦漩結構的 演變，與升力、阻力有很密切的關係。所以，能精確的知道不同升力 與阻力的變化關係，對於流體分離方面的控制方面將更有幫助。實際 工程的應用有:直升機漩翼攻角的改變，渦輪葉片的扭曲震盪，飛機 攻角的急遽變化或遭遇陣風的情況下非穩態的機翼流場，這些細部的 流場變化與升力、阻力關係更有助於應用。

用一個簡單的機械方式來加強升力就是 Gurney flap。大約在 1971 年的時候，這個裝置最初是由一位 Dan Gurney 使用在賽車上。當初 的目的是想要增加賽車的下壓力，使賽車在直線和彎道上速度增加，

為此所以賽車的尾翼採用倒翼。當時的創意，引起他的競爭對手的好 奇，所以他企圖掩飾事實，欺騙他的對手不讓他們他知道真正的目 的。幾年之後他與道格拉斯航空公司的 Bob Liebeck 合作，並討論的 機翼設計的想法與空氣動力學，然後 Liebeck 經過了一連串的測試把 這個襟翼命名為”Gurney flap”。Liebeck 並於 1976 年發表於美國航太 航空研究單位，爾後也激起許多學者的討論，也有許多應用的實例:

飛機上的應用，在起飛與降落的過程中，利用機翼後緣可調式襟翼，

增加空氣動力特性，盡量增大升力，還有在高攻角與有陣風的情況下

不產生失速，增進飛行安全。應用於直升機上面，可以使直升機在橫 向操作上較穩定，並且可以使在起飛時有足夠的升力。另一應用於船 舶，在螺旋槳後緣進行後緣提高，增加船舶的最大推力，同時這樣的 設計也可以減少噪音的產生，達到空蝕現象和鳴音的控制。在風力發 電也是相同的應用，通常安裝在葉片內壓力面的後緣，藉由翼型後緣 的提高，進而增加升力，同時阻力也提高，但是升力係數與阻力係數 的比值是變大的。

1.2 文獻回顧

Metha[1]分析俯仰擺動的 NACA-0012 機翼，經由求解數值上的 二維不可壓縮的 Navier-Stokes 方程式，解析雷諾數(Re)5000 和 10000 的流場變化。Ono[2]等人應用分離渦漩估算法，記算 NACA-0012 機 翼俯仰擺動的流動，獲得相當好的結果。Tassa 和 Sankar[3]使用可壓 縮 Navier-Stokes 方程式，用於求解在高雷諾數下擺動的 NACA-0012 機翼。線性非穩態機翼理論由 Theodorsen、Von Karman 和 Sears[4]

開發，已經廣泛地應用在低於動態的失速角度的小幅度的振動分析，

然而，為了處理動態失速，這樣的理論是很重要的。

Walker、Helin 和 Strickland[5] 以 低 速 風 洞 研 究 NACA0015 研究大角度俯仰運動，在雷諾數(Re)45000 的情形下，減角速度 0.2

速儀量測速度。結果發現，減角速度越大表面流場分離角度就越大，

且機翼前緣產生之渦漩強度就越強。

Francisc 和 Keesee[6]研究 NACA0012 的機翼，它們放置更多的壓 力探針，觀察上下翼面壓力分布狀態，減角速度(reduced angle velocity) 範圍 0.001-0.21，角度變化於 20 度至接近 60 度。結果說明，當減角 速度越大，在運動過程中最大升力係數也就越大。

Hsiau[7]利用拖曳式水槽作研究，觀察 NACA63-018 機翼，減角 速度 0.15 和 0.2，攻角變化 0-60 度，他的結論說明當雷諾數越大，機 翼後緣產生的剪切層渦流強度就越大，越不容易往下游移動，並發 現，減角速度越大，動態失速渦漩就需要越高的攻角才會形成。

Ohmic[8]以流場顯影的技術，觀察 NACA0012 機翼於不同減頻率 及旋轉軸心時，對流場的之影響。使用保角映射的方式，討論流場變 化之改變。Ta Phuoc Loc 和 bouard[9]說明 NACA0012 機翼和橢圓翼 之間的差異，在 NACA0012 的前緣，較容易產生柔順的分離點和分 離方向，由於取率較小，使上下表面渦流的交互作用減小。在機翼尾 緣的部份，因為是尖銳的端面，所以流場分離更容易發展且穩定。

Shih[10]等人利用質點影像位移測速法(PIV)量測，並於數值上使 用離散渦元法研究，NACA0012 機翼以等角速度作俯仰運動之後固定 攻。減角速度 0.131，攻角變化從 0-30 度，之後維持再 30 度，在文

中提到:機翼前緣產生的渦漩發展，會影響到下游的流場行為;並且當 減角速度越大，渦漩分離的狀況會延後。

以上文獻中，不論任何機翼，擺動速度越快，上翼面流場分離與 渦漩的形成、分離，都會發生延後的現象，但是這樣的發展可以使最 大升力係數值提高。

當流場不只向上擺動，且向下擺時會形成一簡諧運動(harmonic

oscillation)，這時流場會有複雜變化。McCroskey[11]等人利用機翼表 面施放煙霧做邊界層分離觀察，以及用熱線測速儀量測表面速度，並 以 壓 力 探 針 量 測 機 翼 NACA0012 正 向 力 和 俯 仰 力 矩 ， 雷 諾 數

，減頻率 0.5，攻角變化

15

2.6 10 ×

6 ^D

+ 10 sin

^D

ω t

結果中發現，升阻

力會隨著週期性的擺動，作週期性的變化。McCroskey 和 Pucci [12]

密集的對動態失速方面的研究做實驗，針對 NACA0012 翼型作不同 的減頻率與高雷諾入研究，過程中研究 Navier-stokes 方程式，考慮到 壓縮性透過使用 Beam—Warming-Steger 方程式模擬，在輕微失速的 過程中，有不錯的模擬結果，但在動態失速方面，無法獲得正確的數 值曲線。並於同年 McCroskey [13]針對一些關於此方面的物理現象做 總結說明:減頻率、擺動振幅、以及最大攻角，如果超過某一程度，

對於渦漩的分離，將會是清楚且有規律。Niu [27]針對次音速下動態 失速分析，多元機翼與加入 Gurney flap 的模擬研究。Niu [28]以

AUSMD 方法，對次音速的動態失速分析。

對於商用經濟效益而言，高升力飛機的應用有極大的影響力。當 給定翼型大小時，一套有效率的高升力系統，可以在令其獲得更大更 有效率之負載能力。新一代升力增加裝置亦可以允許較大的攻角，同 時在起飛時能減少機場周圍產生的噪音問題。升阻力比值的增加，使 飛機可以更快更省油的到達巡航高度。最後，簡單的高升力系統將可 以大幅度地減少建造和維護成本，進而提高飛機的營利能力。

其中一個探討的技術就是 Gurney flap。自從 Dan Gurney 將 Gurney flap 應用在賽車上增加下壓力就有多名學者對這個裝置很有 興趣，Neuhart 和 Pendergraft [15]，Storms 和 Jang [16]，Ross 等人[17]

以及 Jeffrey 和 Zhang [18]。這個裝置應用在機翼上，最初是由 Liebeck[19]開始研究，起初 Liebeck 在 Newman 機翼上，加入 1.25%

的弦長研究 Gurney flap，結果顯示它可以提高升力，和減少些微阻 力。然後研究 2%弦長 Gurney flap，發現雖然升力也有增加，但是阻 力也明顯增加，所以 Leibeck 提出假設尾緣的流場和 Gurney flap 的關 係。之後 Lemay 等人[20]也進行根機翼相關的研究，討論升阻力和 力矩方面的影響。Lyon 等人[21]和 Bechert 等人[22]使用 HQ17 和 SD7037 機翼，針對不同雷諾數的大小做實驗分析與模擬結果相互比 較驗證。Jang、Ross 和 Cummings[23]進行 NACA4412 機翼研究，使用

不可壓縮 Navior-Stokes 方程式，計算 Gurney flap 為 0.5-3.0%的弦 長，發現增加了 Gurney flap 之後，儘管在計算上升力的提升不是線 性，但是卻與 Gurney flap 的尺寸大小有關。壓力分布的計算也說明 了增加的負載相較於標準機翼，特別是在接近尾緣的部份會有一個吸 力峰值。計算也發現，在低或中間的升力係數時會有阻力的損失，但 是 L/D 卻是增加的。

Neuhart 和 Pendergraft[24]使 用 NACA0012 的 翼 型 進 行 多 個 Gurney flap 水穴的實驗。透過染色液，顯示 Liebeck 的假設流場會因 為 Gurney flap 產生改變是正確的。於是 Neuhart 等指出 Gurney flap 會影響機翼尾緣的局部流場。Roesch 和 Vuillet[25]用風動實驗研究直 升機的 NACA5414 機翼，尺寸為 1.25%和 5%的弦長，實驗結果顯示 升力係數的大小與 Gurney flap 的大小有關。並提出 5%弦長時，升力 係數增加了 40%，且升力曲線的斜率提高 6%，零升力攻角降為-6 度。

第二章 數值模式

2.1 統御方程式

計算流場變數以求解統御流場特性之 Navier-Stokes 動量方程式為 主，假設把工作流體是為不可壓縮流，則方程式為:

0

∇⋅ =v (1)

( ) 2

v v v p

ρ^∂t +ρ ⋅∇ = −∇ + ∇μ

∂

K vK

(2)

0 T

μ μ= +μ (3) 其中式(3)之μ₀為黏滯性係數，而μ_T則為紊流黏滯性係數，μ_T裡必須 利用紊流模式導入於數值計算過程中。在進行低雷諾數之氣動力特性 分析時，若採用層流模式，此時μ_T可視為 0。

2.2 數值方法

Fluent有兩種數值求解方法，分別為耦合解法(coupled)和分離解 法(segregated)。此兩種解法皆是解積分形式的守恆方程式，包含動 量、質量、能量、和使用者自行定義的標量，且此兩種數值方法皆以 有限體積法(Finite volume method)為基礎，差別在於線性化的方法與 離散方程式的解析方式有所不同。

耦合解法和分離解法的區別在於連續方程式、動量方程式和能量 方程式求解順序的不同，分離解法是按照順序依次解之，而耦合解法

則是合併在同一個線性系統之下同時求解。分離解法主要是使用於不 可壓縮流場和一般可壓縮流場；而耦合解法則主要是用來解高速可壓 縮流場。雖然耦合解能提供較快的收斂速度，但是所需的記憶體容量 為分離解法的1.5至2倍。而在本論文中所探討的不可壓縮流體的物理 現象時，因此採用分離解法來求解物理問題。

壓力與速度耦合關係:

壓力項在 Navier-Stokes 方程式中，是較為不好處理的流場變數，

若壓力為已知，則將壓力梯度計算結果直接倂入源項，在計算過程中 就與一般傳輸方程式是相似的。但是，實際上壓力一般都是未知，若 計算範圍為可壓縮流場我們可以用動量守恒方程式及連續方程式，求 解出速度場(u,v,w)，在透過狀態方程式(State Equation)將壓力場 P 求 解出來，完成一個封閉的方程組，解得完整之流場變數。

目前可以選擇的方法有很多，例如:人工壓縮法、壓力校正法、

SIMPLE法和渦元法...等，然而FLUENT裡使用的方法有:SIMPLE法、

SIMPLEC法和PISO法，以下將簡單敘述使用FLUENT的方法。

SIMPLEC法:

SIMPLEC法是由學者Van Doormal and Raithby依據SIMPLE法再

做修正獲得，不同在於SIMPLE法中的速度修正方程式裡，忽略速度 修正和擴散特性，所以多少會影響到原方程式，因此SIMPLEC將其改 善，其假設速度修正項在網格單元內各位置的值都等於單元中心的 值，再進行簡化降低對原式的影響，所以他們的差別只在於方程式 中，某幾項係數不同。

2.3 力源理論

一般文獻中，在計算、升阻力與力矩時，最常使用的方法是：先 求得物體表面之壓力分佈與摩擦力，之後再沿著物體表面積分。摩擦 力和表面渦度有關，且物體表面之壓力分佈情形，可經由解動量方程 式而得，優點是：簡單且不需花費太多CPU計算時間，但是，對於整 個流場對物體受力影響無法清楚得知。於是導入力源理論Chang等人 [29][30]用於不可壓縮來幫助流場作分析，此概念源自於D’ Alembert 理論，在均勻入流的條件之下，不可壓縮勢流場中的流體元素對物體 不會有力的貢獻，然而真實流場中，任何具有渦度的流體元素都可視 為 物 體 的 力 源 。 基 於 此 概 念 ， 力 源 理 論 中 ， 引 入 一 輔 助 勢 流 與 Navier-Stokes方程式作內積，再求其體積分，可求得壓力項對物體所 造成的外力，由此引出升、阻力元素的概念來探討流場個別渦旋對 升、阻力的實際貢獻。在非穩態的流場與物體的關係，可以根據這些

做法來更深入了解流場性質。

為了求得作用一個有限物體之各力元元素，首先介紹輔助勢流，

定義

φ φ = ( ) x

在整個流域內的輔助勢流，

φ

可視為對應於物體於某特 定方向運動下所誘發的速度勢能，其解除了滿足

∇ =

²

φ 0

_{外，在遠離}

物體處迅速遞減至無窮遠處趨近於零。

，

φ

在遠域的二維通解可表示為:

( ^A ) ^{log r ... A r ...,}

φ = − ⋅∇ + = − r ^⋅ K +

(7)

而其勢流速度為:

( ) ( )

2

log ... 2 A r r A ...,

A r

φ ^⋅ r ⁻

∇ = − ⋅∇ ∇ + = − +

K K

(8)

其中上述所有方程式中，

r K

為位置向量 為單位向量，向量A與物 體形狀和物體運動型態有關，但與座標位置無關，A可藉由求解勢能 函數 而得。由上述方程式(8)得知，在二維流場中，輔助勢流 和勢流速度則分別以

r

2

2

φ 0

∇ =

1 r

_和

1 r

²的速度遞減。

力源理論升力和阻力可以表示為下面兩式:

v 1 2 2

n v n v n ( j)

2 Re

S S V S

CL dA dA dV dA

φ^∂t φ ω φ ω

= ⋅ + ∇ ⋅ − × ⋅∇ + × ⋅ ∇ +

∫

∂

∫ ∫ ∫

^{φ ( 9 )}

CL= CL_a+ CL_m+ CL_v+ CL_s

其中，CL_a為表面加速度之升力貢獻、CL_m為表面速度之升力貢 獻、CL_v為環境渦度之升力貢獻、CL_s為表面渦度之升力貢獻。

v 1 2 2

n v n v n ( i)

2 Re

S S V S

CD dA dA dV dA

φ^∂t φ ω φ ω

= ⋅ + ∇ ⋅ − × ⋅∇ + × ⋅ ∇ +

∫

∂

∫ ∫ ∫

^φ (10) CD= CD_a+ CD_m+ CD_v+ CD_s

其中，CD_a 為表面加速度產生之阻力、CD_m 為表面速度產生之阻 力、CD_v 為環境渦度產生之阻力、CD_s 為表面渦度產生之阻力。

第三章 數值結果

3.1 紊流模式之比較

在本章第一節的部份，我們採用 Fluent 提供的紊流模式來做比

較，分別是：單方程 Splart-Allmaras(S-A)、多方程 Reynolds Stress model(RSM)。

參考 McCroskey 等人[14]所發表之文章得知 Spalart-Allmaras 在大 部分模擬，成效不錯，McCroskey and Pucci [12]的實驗值與 Shida 等 人[26]所發表的文章，由使用相同的條件來做模擬，以 NACA0012 為 模擬的機翼，馬赫數(M)=0.3，雷諾數( )大約 ，減頻率 (reduced frequency)=0.1，擺動的角度以

1 0

Re 4 10 ×

⁶

1 0 s in ω t +

D D

，藉由模擬 嚴重失速的現象，檢驗紊流模式結果比較。所使用的網格總數為 80793，其幾何外型網格如圖(1)及圖(2)所示。

由於馬赫數是 0.3，再程式中選用可壓縮的條件模擬。圖(3)、

圖(4)、圖(5)為為一個擺動週期中和不同攻角的升力係數、阻力係數 和力矩係數圖。

由升力係數圖中顯示，在這一個週期之變化中，升力在一開始快 速下降至一個最小值，這時候攻角有一個最小值，之後攻角向上在時 間 0.385 左右升力達到最大值，然後又開始下降，之後不斷的循環。

在數值結果方面，一開始只有

RSM

和

S − A

的升力係數值比較接近 與其他比起來有差異，在中段升力提高的部份，紊流方程式所得到的 趨勢接近，在最高點的時候又有差異，最後達到失速下擺至 10 度又 回到開始的階段。在與實驗值比較，兩個紊流模式得到結果，皆與實 驗值不同，發現只有最大與最小值接近，中間變化過程完全不相同。

阻力係數方面，在低攻角的時候阻力會達到一個最低值，隨著攻 角的提高阻力也跟著提高，當達到升力係數最高點的同時，紊流方程 式所得到的值 和 很接近，只有 15 度至 0 度有些微落差。

與實驗值比較的結果也不相同，在最大值的地方就有顯著的差距，中 間變化過程也是差異明顯。

S − A RSM

力矩係數方面， 和 顯示在 0 度至 4 度範圍內，有不 同趨勢發展差距很大，大體上趨勢相近。實驗值比較上，最大最小值 範圍在接近的位置，中間變化過程，如同升、阻力一樣，差異明顯。

S − A RSM

圖(6)與圖(7)分別為

S − A

和 機翼表面壓力分布情形。比 較發現，在三個不同攻角的表面壓力分布狀況，從 Cp 值中，最大和 最小值的壓力係數接近相同，唯一不同在某些表面位置壓力係數不 同，推測原因為:在整體上

RSM

S − A

為單方程方程式， 為多方程 方程式，比較起來多方程的計算考慮，會包含更多物理意義，但在相 對簡單的 模型中與 模擬出來的結果是差不多的，所以

RSM

S − A RSM

RSM

S − A

單方程會比 節省很多時間，但是較複雜問題時，

RS

可 以比 從流場中得到更多的資訊。

M S − A

圖(8)~圖(9)分別為不同紊流模式下，週期循環角度變化之流線 圖。

角度 10 度時，

RSM

和

S − A

的可以看出，再上翼面皆有渦漩 產生，然後 方面有看似在大窩全前面再產生小渦漩的現象。在 20 度時也只有平順的流線流經機翼，觀察之前機翼表面壓力分布圖 得知，在此時壓力變化差距不大所致。當回到 10 度時，又有渦漩產 生。再往下擺動達到 0 度時，紊流模式除了

RSM

S − A

， 只有平順 的流線變化，藉由觀察之前翼面壓力分布狀態，只有在 模式 下，壓力有變化。

RSM S − A

在此模擬中發現，紊流模式所呈現的結果，

S − A

和 得 到結果與文獻 Shida 等人[25]模擬的結果，失速前有些微的接近，但 是在達到失速的角度之後，仍然有差異，達到失速的角度時，文獻中 皆有明顯的渦漩產生，而吾人模擬之結果，只有在尚未達到失速前接 近。會有如此的結果可能為，Fluent 的數值消散項太強所致，導致渦 漩產生受到消散項的影響，無法得知真正的流場現象。

RSM

然而在與實驗值部份做比較發現，吾人的模擬值與實驗值差異甚 大，然而文獻中模擬的結果也有顯著的差異，証明當達到失速的條件

時，我們所模擬的程式與實際上的現象還是有所差異，在流體力學研 究領域中，討論紊流領域方面的研究，仍然有持續進步的空間。

圖(1) NACA0012 幾何外型整體網格圖

圖(2)NACA0012 局部網格放大圖

(a) 攻角與升力係數變化圖(實驗值)(ref. 11)

(b) 攻角與升力係數變化圖

圖(3) (a)攻角與升力係數變化圖(實驗值)、(b)攻角與升力係數變化圖

(a) 攻角與阻力係數變化圖(實驗值) (ref. 11)

(b) 攻角與阻力係數變化圖

圖(4) (a)攻角與阻力係數變化圖(實驗值)、(b)攻角與阻力係數變化圖

(a) 攻角與力矩係數變化圖(實驗值) (ref. 11)

(b) 攻角與力矩係數變化圖

圖(5) (a)攻角與力矩係數變化圖(實驗值)、(b)攻角與力矩係數變化圖

(a) 10 度

(b) 20 度

圖(6)

S − A

週期壓力係數等值圖與機翼表面壓力係數分佈圖

(c) 10 度

(d) 0 度

圖(6)

S − A

週期壓力係數等值圖與機翼表面壓力係數分佈圖(續)

(a) 10 度

(b) 20 度

圖(7)

RSM

週期壓力係數等值圖與機翼表面壓力係數分佈圖

(c) 10 度

(d) 0 度

圖(7)

RSM

週期壓力係數等值圖與機翼表面壓力係數分佈圖(續)

(a) 10 度 (b) 20 度

(c) 10 度 (d) 0 度 圖(8)

S − A

週期流線圖

(a) 10 度 (b) 20 度

(c) 10 度 (d) 0 度 圖(9)

RSM

週期流線分佈圖

3.2 藉由力源理論在馬赫數(M)=0.1 的動態失速現象解析

由前一節的結果得知，Spalart-Allmaras & Renolds Stress 兩個紊 流模式得到的結果比較接近，由 McCroskey 等人[14]所發表之文章得 知 Spalart-Allmaras 在大部分模擬，成效不錯，且 Spalart-Allmaras 的 計算時間上比較節省，所以在這章節以 Spalart-Allmaras 紊流模式來 模擬。需要設定的條件為: 以 NACA0012 為模擬的翼型，馬赫數 (

M

)=0.1 ， 減 頻 率 (reduced frequency)=0.1 ， 擺 動 的 角 度 以

1 0

^D

+ 1 0 s i n

^D

ω t

，以上述條件模擬嚴重失速的現象，並加入力源 理論，分析流場受力變化情形。

圖(10)表示升力元素對升力的貢獻、攻角與升力係數的變化圖。

CL 為透過力元理論分析所獲得之總升力、CL_m 為機翼表面速度之 升力貢獻、CL_a 為機翼表面加速度之升力貢獻、CL_v 為環境渦度之 升力貢獻、CL_s 為機翼表面渦度之升力貢獻。在數值上有幾個階次 上的差距，但趨勢上是相同的，有一定的參考價值。

觀察力源理論結果，可知升力貢獻情況:整體上 CL_s 在起始時呈 現正值貢獻，但在整體上對升力係數影響微乎其微，對升力變化毫無 影響。換言之對升力有影響的就是 CL_v、CL_a、CL_m。在 CL_m 方面，大部分呈現負值的貢獻，只有一小部分呈現正值，整體上呈負 值貢獻，CL_s 因為因為流體力學的定義在高雷諾數的情況下表面摩

擦力會減小，所以影響不大。然而在 CL_a 貢獻方面，在一開始只有 CL_a 對升力影響是有利的，只有 CL_a 沒有因為攻角的變小而變小，

主要是因為他與運動的加速度有關，無關角度變化且呈現一定的週 期，所以 CL_a 是有對稱性。就 CL_v 而言，整體 CL_v 的主要為正 值貢獻，影響機翼的整體升力，所以可以確切的說環境渦度為主要影 響，在一開始角度下降，升力也跟著下降，表示角度的變化影響著流 場渦度，在 t=0.307 左右有個最低值，之後開始明顯的增加，一直到 最大值約 t=0.32 左右，在此時剛好接近失速角度，升力開始下降至 t=0.353 時，CL_v 開始負值貢獻，直到下一個週期。在此發現 CL_a 與 CL_m 瞬時變化與大小完全取決於機翼的幾何形狀與其所對應的 運動模式(擺動幅度)； CL_v 和 CL_s 之變化由瞬間流場結構決定。

圖(11)表示阻力元素逐時變化圖、攻角與阻力係數的變化圖。由 阻力係數逐時變化圖，可了解阻力來源，同樣的力源理論分析的總阻 力值與其它的阻力元素有階次上差異，但整體趨勢上相同，也具有參 考價值。

同樣一開始也跟升力相同，主要的阻力來源來自於 CD_s，但整 體上影響是次要的。隨著時間的改變，變成由 CD_v 主導整個阻力，

且 CD_a 雖然比例上阻力來不大，但是某個時間點提供的阻力遠比 CD_s 大，短時間成為次要影響。在這邊 CD_m 和 CD_a 相比，CD_m

部分全部為負值，表示他有抵銷阻力的效果。

圖(12)為一週期中升力與阻力元素等值線圖與圖(13)流線圖相互 參照。根據文獻 Chang[26]敘述：任何具有相同渦度方向之流場可能 具有不同方向的升阻力元素(即可提供正或負之貢獻)，而決定某渦度 為正或負貢獻則決定於機翼位置、攻角與瞬時之流場行為。圖中紅色 區域為正值貢獻區、藍色區域為負值貢獻區。

(a)-(c)中，右為阻力元素圖，從 10 度上擺至 20 度，正值區在尾 緣厚度增加迅速，委緣處為主要要主力來源(正值區)，雖然機翼後端 有負值區產生，但影響部分不大，隨著越往下游 影響越小；(a)-(c) 左為升力元素圖，升力元素正值區在後緣增厚，主要升力來源，尾緣 下游處也有負值區產生，與阻力元素同，越遠離影響越小。經流線圖 得知，此處的渦漩對升力與阻力而言都有影響。

(d)下擺至 15 度，從圖(13)流線圖中發現，前緣處又開始產生渦 漩，對升力來說有正值貢獻(升力)，阻力來說負值貢獻(推力)，且從 圖中也發現下游處，升力與阻力元素的正值區，皆會受到負值區影響 推擠，網機翼前緣方向移動，接近尾緣處都有被負值區往上推擠的現 象發生。

(e)右為阻力元素等值線圖，機翼表面分佈情形來看，前緣的負值 貢獻區為主要推力來源(負阻力)，紅色正值區有部分與機翼尾緣負值

區抵銷； (e)左為升力元素等值線圖，參照圖(14)流線圖(e)，得知機 翼產生之前元渦漩對升力是有幫助的，因為角度 10 度上下翼面流速 度不同，正值區貢獻在上翼面明顯，負值區只在尾緣或更下游區增 加，當更遠離機翼時，負值影響越小，顯示尾緣區只有部份影響。

(f)下擺至 5 度，上翼面前緣渦漩，參照圖(14)流線圖(f)發現有稍 微往後移動的趨勢，觀察圖右阻力元素圖，前緣渦漩還是負值貢獻，

但是在近翼面處有較小正值區產生，蔣負值貢獻區往上推擠，渦漩之 後上翼面仍然是正值區，主要阻力來源。

(g)右為阻力元素等值線圖，在下擺至 0 度時，前面負值區往後 移動，正值區佔領整個機翼後緣部份，此時的阻力較大。圖左顯示比 對流線圖(14)-(g)，發現上翼面渦漩向後移動，此渦漩對升力還是正 值貢獻居多，至尾緣處上翼面正值升力貢獻居多。

(h)上擺至 5 度，渦漩往後移動，圖右阻力元素圖中，在上百過 程中，渦漩被視為祖力來源，且往後延伸，使尾緣部份正值區(阻力) 加大；圖左升力元素中，只要渦漩產生，仍舊對升力是有利的，且正 值區網後增厚，推擠負值區(負升力)至下游處。

由此可知，力源理論可以幫助了解流場中的受力情形。將有助於 分析哪些部份影響升阻力情形最大。

(a) 表示升力元素週期變化圖

(b) 攻角與升力係數的變化圖 圖(10)

(a)表示阻力元素週期變化圖

(b)攻角與阻力係數的變化圖 圖(11)

(a) 10 度(上擺行程)

(b) 15 度(上擺行程)

(c) 20 度

圖(12) 週期中升力(左)與阻力(右)元素等值線圖

(d) 15 度(下擺行程)

(e) 10 度(下擺行程)

(f) 5 度(下擺行程)

圖(12) 週期中升力(左)與阻力(右)元素等值線圖(續)

(g) 0 度

(h) 5 度(上擺行程)

圖(12) 週期中升力(左)與阻力(右)元素等值線圖(續)

(a) 10 度(上擺行程) (b) 15 度(上擺行程)

(c) 20 度 (d) 15 度(下擺行程) 圖(13) 週期中流線圖

(e) 10 度(下擺行程) (f) 5 度(下擺行程)

(g) 0 度 (h) 5 度(上擺行程) 圖(13) 週期中流線圖(續)

3.3 Gurney flap 流場現象解析

接下來將討論在機翼尾緣處加 Gurney flap，並討論其現象。參考 文獻 Jang、Ross 和 Cummings[23]發表之文章，並使用相同條件模擬。

使用 NACA4412 機翼模擬，馬赫數(M)=0.085，雷諾數(Re)大約

，使用 Gurney flap 的尺寸 0.5%-15%的弦長。網格參考文 獻建構方式，圖(14)~圖(16)所示，網格數約 40000。

1.64 10 ×

6

圖(17)為本文與文獻的比對結果，結果非常接近，只有在最大值 的部份與接近 0.1 附近的轉折處略有不同。在機翼表面壓力分布情 形，如圖(18)所示，為四種不同尺寸 Gurney Flap 與 NACA4412 機翼 作比較，以 NACA4412 機翼情形顯示，壓力分布在尾緣的部份會有 閉合情形，表示在尾緣最尖端的部份，壓力分佈相同，有 Gurney flap 之後，在尾緣部份就還是能保有壓力差；在不同尺寸的比較，當尺寸 越大，上下翼面壓力差值就會增大，這樣升力就會增加。圖(19)表示，

不同長度 Gurney flap 與 NACA4412 的表面分壓力分布情形，(a)為在 沒有 Gurney flap 的情況時，因為翼型的關係，上下表面流速不同造 成的壓力差；(b)-(g)為不同尺寸弦長 Gurney flap 的情形，同樣會因為 機翼形狀的關係造成壓力差，但是加入 Gurney flap 之後，在尾緣的 下游處也有低壓產生，同時在 Gurney flap 內側會有一個較大的壓力 區產生，尾緣下游處的低壓和內側的高壓，都會因為 Gurney flap 的

長度而改變，隨著長度的增加，內側高壓區與外側低壓區分佈會越來 越大。

圖(20)表示，不同 Gurney flap 尺寸下升力係數的改變，在文獻 中只有從 0.5-3%部分的模擬，與目前模擬的值比較有接近。文獻中 最大部份只有作到 3%弦長 Gurney flap，我們考慮，如果飛機在巡航 狀態下，在 Gurney flap 長度增加多長的情況下會有最大升力產生，

於是繼續增加長度，發現長度不斷加長升力也不斷提高，在大約 10%

的長度時，對於 NACA0012 的機翼，升力有一最大值，之後就呈現 下降趨勢。

圖(21)為不同 Gurney flap 尺寸下，阻力係數的改變，阻力在 0.5%-1.25%之間的弦長，有上升但是幅度不大，從 1.25%開始阻力突 然增大，到 2%又平緩，之後長度的阻力就一直不斷增加。圖(22)為 升阻力比與不同 Gurney flap 尺寸關係。圖中發現在 2%時，升阻力比 會有一個峰值，在長度比較小的條件下 2%的弦長會是個最佳尺寸，

有足夠大的升力，同時阻力也小；直到 6%-7%之間升阻力比有個較 小的值，之後到 10%又有一極值，但是升阻力比與在 2%的情形下小，

再觀看升力圖(20)，在 10%弦長時，升阻力比值較小，但是他有最大 升力，此時需要最大升力值，在不考慮阻力增加的情形下，10%的襟 翼長度會有最大升力。

由於計算結果發現在 10%長度時會有最大升力的情況，所以針 對七個狀態 1%、2%、4%、7%、10%、15%和 NACA4412 機翼繼續 討論。接下來使用力源理論分析，圖(23)表示在 1%、2%、4%、7%、

10%、15%弦長的 Gurney flap 與 NACA4412 機翼升力與阻力元素等 值線圖；並配合圖(24)流線圖。

(a)為 NACA4412 機翼，觀察圖右阻力元素部分:藍色部分為負值 區，在機翼前緣處上下部分與在後緣下游區產生，紅色正值的部份，

包覆機翼上翼面中間至後緣部份，在從下翼面後緣部份到接近前緣；

整體上，正值區包覆區域較大，負值區包覆於前緣部份，對機翼來說 影響阻力的因素在於機翼後緣，所以對正值部分加以控制能有效降低 阻力；觀察圖左升力元素部分:因為機翼形狀的關係，上翼面流速較 快，正值區域提供的升力較下翼面負值區大；且下游處有對稱正負值 區相互抵銷，因此對機翼而言，升力元素等值圖呈現此分佈型式。

(b)-(g)圖中，發現 NACA4412 尾緣部份加入 Gurney flap 之後，

觀察圖右阻力元素等值線圖中，發現受力情形更複雜。(b)-(c)圖右，

為 1%弦長的 Gurney flap 阻力元素分布情形，在 Gurney flap 後方有 產生兩個小渦璇(流線圖(24)-(b)與(c)所示)，觀看力源元素分佈情形，

在 Gurney flap 內側有一較小的負值阻力產生，與 NACA4412 機翼相 比，因為受到 Gurney flap 的影響，在靠近 Gurney flap 的下翼面阻力

有分離的情形，內側額外產生負值貢獻區，對機翼而言為推力，下游 區兩對稱正負元素會互相抵消，使的 Gurney flap 下游區產生的阻力 小。(b)-(c)圖左，在機翼下游區變化不大，主要差異在於下翼面，加 入 Gurney flap 與未加入的 NACA4412 相比，再下翼面負值貢獻區有 減少的趨勢，表示增加之後升力增大的趨勢。

(d)-(g)圖右阻力元素等值圖部分，會在下游區產生額外一組較小 的對稱正負值的區域，同時比對流線圖，產生的渦漩比 1%與 2%弦 長時還大，且跟 1%與 2%弦長比起來，使的原來較大的兩個正值與 負值區，有分離的現象發生，隨著 Gurney flap 的長度增加尾緣部份 較大兩個正負值區域，分離的越明顯，所以受力有分離的現象，分離 之後受力不平衡使的阻力增加；(d)-(g)圖左升力元素等值圖部分，當 Gurney flap 尺寸增加，同樣和阻力一樣在下游區額外產生一組對稱正 負值區，使得下游區原來正負值區同樣有分離的現象，觀察圖(24)流 線圖中，推測在 15%時尾緣渦流作用劇烈，使的升力下降。

在 2%弦長時，有較大的升阻力比，是由於渦漩較小，產生的正 負值區較集中，有互相抵銷的結果，於此產生的阻力就減少；在 15%

弦長時，受力的情形，因為 Gurney flap 後緣渦漩較大，會使尾緣下 游處正負值區有分離的現象，反而不利減小阻力，所以 15%弦長時阻 力有增加的趨勢。

表(1)和表(2)顯示 NACA4412 不同尺寸 Gurney flap 之升力和阻 力元素貢獻。CL 為透過力元理論分析所獲得之總升力、CL_v 為環境 渦度之升力貢獻、CL_s 為機翼表面渦度之升力貢獻；CD 為透過力元 理論分析所獲得之總阻力、CD_v 為環境渦度產生之升力、 CD_s 為 機翼表面渦度產生之阻力。

由表(1) 中知道，主要的升力來源，都來自於 CL_v；表中也發 現在沒有加入 Gurney flap 的情況下，總升力與 2%弦長相比小了很 多，在加入之後，發現 CL_s 和 CL_v 增加升力幅度很高，到 10%弦 長時，都還是成長的，15%時升力沒有與長度增加成正比，反而減小。

表(2)中得知，沒有加入 Gurney flap 與 2%弦長相比，阻力比 2%弦長 時還大，且發現 2%弦長時，CD_s 還負值貢獻；但是沒有像升力一樣 的表現方式，長度增加，阻力也明顯增加，2%-7%與 10%-15%增加 的幅度不小。從表中也驗證了 Jang 等人[22]發現，增加 Gurney flap 之長度，升力成長不是線性，且到達一定長度之後會有下降趨勢；同 時在 10%弦長時，有最大升力值表現，但是由圖(23)了解，L/D 並不 是最大值。

Gurney flap size CL CLs CLv

clean 0.270607074 0.050826037 0.219781037 1.00% 0.68236871 0.165385491 0.516983219 2.00% 0.84558283 0.222246815 0.623336015 4.00% 1.0108742 0.269628958 0.741245242 7% 1.1279781 0.28597005 0.84200805 10% 1.8262916 0.598466974 1.227824627 15% 1.341582 0.3027745 1.0388075 表(1) NACA4412 不同尺寸 Gurney flap 之升力元素貢獻

Gurney flap size CD CDs CDv

clean 0.065288564 0.021775429 0.04351314 1.00% 0.039091713 -0.010941536 0.050033239 2.00% 0.045996289 -0.011701722 0.057698011 4.00% 0.076345883 -0.005856243 0.082202126 7% 0.10685768 0.02470202 0.082155661 10% 0.11791696 -0.012889392 0.130806352 15% 0.17178639 -0.014140805 0.185923195 表(2) NACA4412 不同尺寸 Gurney flap 之阻力元素貢獻

圖(14) NACA4412 幾何外型整體網格圖

圖(15) NACA4412 局部放大網格圖

圖(16) NACA4412 加 Gurney flap 局部網格放大圖

圖(17) NACA4412 表面壓力分布圖

圖(18) NACA4412 與 NACA4412 加 Gurney flap 機翼表面壓力分佈圖

(a)NACA4412 機翼 (b)1%弦長 Gurney flap

(c)2%弦長 Gurney flap (d)4%弦長 Gurney flap

(e)7%弦長 Gurney flap

(f)10%弦長 Gurney flap (g)15%弦長 Gurney flap 圖(19)機翼表面壓力分佈等值線圖比較

圖(20) 升力係數與 Gurney flap 尺寸變化圖

圖(21) 阻力係數與 Gurney flap 尺寸變化圖

圖(22) 升力阻力比值與 Gurney flap 尺寸變化圖

(a)NACA4412 機翼

(b)NACA4412 加入 1%弦長 Gurney flap

(c)NACA4412 加入 2%弦長 Gurney flap 圖(23) 升力元素(左)與阻力元素(右)等值線圖

(d)NACA4412 加入 4%弦長 Gurney flap

(e) NACA4412 加入 7%弦長 Gurney flap

圖(23) 升力元素(左)與阻力元素(右)等值線圖(續)

(f) NACA4412 加入 10%弦長 Gurney flap

(g) NACA4412 加入 15%弦長 Gurney flap 圖(23) 升力元素(左)與阻力元素(右)等值線圖(續)

(a) NACA4412 機翼 (b) 1%弦長 Gurney flap

(c) 2%弦長 Gurney flap (d) 4%弦長 Gurney flap

圖(24) NACA4412 加 Guenwy flap 流線圖

(e) 7%弦長 Gurney flap

(f) 10%弦長 Gurney flap (g) 15%弦長 Gurney flap

圖(24) NACA4412 加 Guenwy flap 流線圖(續)

第四章 總結

經由以上的模擬，對機翼的週期性俯仰擺動與 Gurney flap，升阻 力變化情形有了初步的認識，在此歸納出以下幾點結論:

1. 對於 Fluent 的四種紊流模式比較，準確度上 Splart-Allmaras 和 Reynolds Stress model 可以比較準確的描述出渦漩的現象，且只有 在失速的前模擬，較接近文獻，失速之後的模擬 就失去準確性。

2. 對於在動態失速結果，升力與阻力，在數值上有幾個階次上的差 距，但趨勢上是相同的，有一定的參考價值。在力源理論上，升 力方面為環境渦度主導整體升力貢獻，加速度之升力貢獻為次要 影響；阻力方面，同樣為環境渦度主導整體阻力的影響，在速度 之阻力影響與加速度之阻力為次要。

3. 在 NACA4412 加 Gurney flap 力源理論解釋方面，主要的升力與阻 力來源，都來自環境渦度的貢獻與影響，表面渦度影響雖小，但 至於到能忽略的地步。尺寸大小方面，隨著增加長度，升力皆上 升趨勢，至 10%之後，反而升力下降；阻力方面，加入襟翼與未 加入相比，在小尺寸時有降低趨勢，隨著尺寸增加，阻力也不斷

上升。

4. 在力源理論上，他對流場可以表現出受力情形狀況，參考文獻中 發現，在低雷諾數的情形力源理論可以準確的預測，在高雷諾數 的情形下，則有差異。需要更深入了解，解決其中差異性。

5. 力源理論目前使用的範圍只介於不可壓縮的範圍，假使可以討論 至可壓縮流場範圍，將可以更深入的突破，對於流場受力情形分 析，將有助於航空工程方面的突破，於此將更具貢獻，是未來值 得探討的地方。

參考文獻

[1]Metha, U.B., "Dynamical Stall of an Oscillating Airfoil." AGARD paper 23., 1977

[2]Ono, k., Kuwahara,k., and Oshima, K., "Numerical Analysis of Dynamic Stall Phenomena of an Oscillating Airfoil by the Discrete Vortex Approximation."

Proceeding of 7th ICNMFD, New York., 1981

[3]Tassa, Y. and Sankar, N.L., "Dynamic stall of NACA0012 Airfoil." AIAA Paper 81, 1289., 1981

[4]Theodorsen, T., "General Theory of Aerodynamic Instability and the Mechanism of Flutter." NACA Rept.496, 1935

[5]Walker, J.M.,Helin H. E.,and Strickland, J.H."An Experimental Investigation of an Airfoil Undergoing Large-Amplitudel Pitching Motions", AIAA Journal, Vol 23, 1141-1142., 1985

[6]Francisc, M.S. and Keesee, J.E., "Airfoil Dynamic Stall Performance with Large-Amplitude motions", AIAA Journal, 23, No.11., 1985

[7]Hsiau,Y.C., and Chu,C.C., "Starting Vortex and Lift on an Airfoil", Physic of Fluids A, 5, No. 3, 662-676, 1993

[8]Ohmic,K., Coutanceau, M., Daube, O. and Ta Phuoc Loc., "Further Experiments on Vortex Formation Around an Oscillating and Translating Airfoil at Large Incidences", J. Fluid Mech., 225, 607-630., 1991

[9]Ta Phuoc Loc and Bouard, R., "Numercial Solution of the Early Stage of the Unsteady Viscous Flow around a Circular Cylinder: A Comparison With Experimental Visualization and Measurements", J. Fluid Mech.,160, 93-117., 1985

[10]Shih, C., Lurenco, L., Van Dommelen, L., and Krothapalli, A."Unsteadt Flow Past an Airfoil Pitching at a Constant Rate", AIAA Journal, Vol.30, No. 5, 1153-1161., 1992

[11]McCroskey, W.J., Carr, L.W. and McAlister, K.W.,"Dynamical Stall Experiments on Oscillating Airfoil", AIAA Journal, 14, 57., 1976

[12]McCroskey, W.J. and Pucci, S.L., "Viscous-Inviscid Interaction on Oscillating Airfoil." AIAA Paper85-0125. 1985.

[13]McCroskey, W.J., "The Phenomenon of Dynamic Stall." NACA Teachnical Memorandum 81264., 1981

[14]Srinvasan, G.R., Ekaterinaris, J.A. and McCroskey, W.J., "Dynamic Stall of An Oscillating Wing Part:1 Evalution of Turbulience Models." AIAA-93-3403., 1993

[15]Neuhart H, Pendergraft O.C. A water tunnel study of gurney flaps., NASA-TM-4071, 1988

[16]Storms BL, Jang CS, Lift enhancement of an airfoil using a Gurney Flap and Vortex Generators. AIAA 0647., 1993

[17]Jang CS, Ross JR, Cummings RM. Numerical investigation of an airfoil with a Gurney flap. Aircraft Design,75-88 ,1998.

[18]Jeffrey D,Zhang S. Aerodynamics of Gurney flaps on a single-element high-lift wing. Journal of Aircraft, 2000.

[19]Leibeck, R. H., "Design of Subsonic Airfoils for High Lift." J. Aircrast. Vol. 15.

No. 9, Sept. 1978, 547-561

[20]Lemay J, Giguere P, Dumas G, Gurney Flap effects and scaling for low-speed airfoils, AIAA 1995-1881

[21]Lyon CA, Broeren AP, Giguere P, Gopalarathnam A, Selig MS. Summary of low-speed airfoil data, Soartech publication, Virginia, 1997

[22] Bechert DW, Meyer R, Hage W. Drag reduction of airfoils with miniflaps. AIAA paper 2000-2315., 2000

[23]Jang, C.S., Ross. J.C.,and Cummings, R.M.,"Computational Evaluation of an Airfoil with a Gurney Flap." AIAA Paper 92-2708-CP, June 1992

[24]Neuhart, D.H., and Pendergraft, Jr., O.C., "A Water Tunnel Study of Gurney Flaps," NASA TM 4071, 1988

[25]Roesch, P., and Vuillet, A., "New Designs for Improved Aerodynamic Stability on Reccent Aerospatiale Helicopters." Vertica, Vol. 6., pp.145-164. 1982

[26]Shida, Y., Kunio, K.,Ono K., and Takami, H.,"Computation of Dynamic Stall of a NACA-0012 Airfoil.", AIAA Journal VOL. 25. No.3., 1987

[27]Niu, Y.Y., "MUSCL Type Limiters for flux Spling Methods and Applications."

Dissertation of Philosophy, 1995

[28]Niu, Y.Y., and Liou, M.S., "Numerical Simulation of Dynamic Stall Using an Improved Advection Upwind Splitting Method." AIAA Journal, Vol.37, No. 11.1999

[29]Chang, C.C., "Potential Flow and Forces for Incompressible Viscous Flow."

JSTOR Vol 437: 517-525.1992

[30] Hsieh,C.T., Hsieh., Chang,C.C., and Chu, C.C., "Revisiting the aerodynamics of hovering flight using simple models." J. Fluid Mech. Vol 623, 121-148. 2009