摘要
以自製差分法求解器求解二維不可壓縮 Navier–Stokes 方程式，重現圓柱尾流的卡門渦街。在筆記型電腦（Intel Core i7-8550U / RAM 8 GB）上，完成雷諾數 3 水準 × 空間網格 3 水準 × 時間步長 2 水準的層流案例，合計約 28.5 分鐘。
本文為經同儕審查之研究筆記的部落格版（僅含不需重新計算的輕微文件修正）。
數值方法
| 項目 | 內容 |
|---|---|
| 控制方程式 | 2D 不可壓縮 Navier–Stokes（無因次） |
| 對流項 | 二階精度上風差分 |
| 擴散項 | 二階精度中心差分 |
| 時間推進 | Chorin 分數步法（投影法）+ 一階顯式 Euler |
| 壓力 | SOR 法（ω=1.7） |
| 圓柱邊界 | 直接強制型嵌入式邊界法（體積率平滑約 1 格寬） |
| 流體力 | 完整動量虧損形式（非穩態、對流、黏性、壓力） |
計算域依設計為 13D × 6D（上游 3D / 下游 10D / 上下各 3D），圓柱直徑 D=1，均勻來流 U∞=1。
計算條件
空間網格
| 水準 | 類型 | 分割 | 最小 Δx/D | 單元數 |
|---|---|---|---|---|
| M1 | 均勻 | 130×60 | 0.100 | 7,800 |
| M2 | 均勻 | 208×96 | 0.0625 | 19,968 |
| M3 | 伸縮（圓柱近旁加密） | 200×100 | 0.042 | 20,000 |
時間步長
- DT1: Δt = 0.008（CFL ≈ 0.13）
- DT2: Δt = 0.004
各案例於卡門渦約釋放 5 個時自動停止。
主要結果
1. 雷諾數與阻力／升力
圖：平均阻力係數（左）與升力 RMS（右）對雷諾數
建議設計 M1+DT1:
| Re | 平均 Cd | Cl RMS | Strouhal |
|---|---|---|---|
| 100 | 1.617 | 0.437 | 0.185 |
| 400 | 1.800 | 0.844 | 0.185 |
| 800 | 1.844 | 0.912 | 0.247 |
Cl RMS 隨 Re 增大，與渦激振動風險上升一致。Cd 絕對值略高於文獻（Re=100 約 1.3–1.4），可歸因於阻塞率 16.7%（域高 6D）。
2. 阻力／升力時間歷程
圖：建議穩定設計 M1+DT1 的 Cd（左）與 Cl（右）（上至下 Re=100, 400, 800）
初始過渡後，Cl 趨近零均值正弦波，Cd 以約兩倍頻率微振，與卡門渦街物理一致。
3. 流速場五階段時間演進
圖：Re=100（M1+DT1）流速 |u|/U 時間演進
圖：Re=400（M1+DT1）
圖：Re=800（M1+DT1）
可見：雙生渦 → 對稱性破壞 → 規則卡門渦街的典型過程。
4. 速度向量（圓柱近旁）
圖：Re=100 最終時刻速度向量。未見網格尺度不自然振盪。
穩定設計結論
「穩定」定義為同時滿足：
- 計算不發散
- 速度向量無網格引起之不自然擾動
- 計算時間與記憶體在筆電上可接受
穩定空間網格
- 建議：M3（近旁最小 Δx/D ≈ 0.042，約 24 分割/D）
- 與 M2（全域均勻、約 16 分割/D）相當的 Cd 精度，成本更合理
- M1（Δx/D = 0.10）供初步／敏感度分析
穩定時間步長
- 建議：DT1（Δt = 0.008）
- 與 DT2 相比平均 Cd 差數 %，St 與 Cl RMS 幾乎一致
- 以成本效益選擇 DT1
圖：各案例壁鐘時間
計算環境
| 項目 | 值 |
|---|---|
| CPU | Intel Core i7-8550U |
| RAM | 8 GB |
| 實作 | Python + Numba JIT |
| 層流案例總時間 | 約 28.5 分鐘 |
| 發散案例 | 0 |
結語
差分法、分數步法與嵌入式邊界法的組合，可在有限資源下穩定重現卡門渦街。建議設計為 M3（伸縮網格）+ DT1（Δt=0.008）。
後續工作：時間二階精度、壓力多重網格、三維擴展、以及 k-ω SST 等紊流模型。
本文於同儕審查通過後自動發表。審查判定：通過（僅需不重新計算之輕微文件修正）。
Top comments (0)