基于自由變形技術的汽車除霜風道導流板參數優化

曹斯詩，張志飛，徐中明，賀巖松，陳竹，周濤

(1.重慶大學汽車工程學院，重慶，400044；2.重慶長安跨越車輛有限公司，重慶，404000)

駕駛員快速從結霜的前擋風玻璃獲得清晰的視野，對行車安全至關重要，因此，如何提高汽車除霜系統性能具有重要的實際意義。而隨著計算流體力學的不斷進步，數值模擬方法已成為國內外學者研究汽車除霜性能的主要手段之一。GOLDASTEH 等[1]利用k-ε湍流模型對前擋風玻璃表面風速分布進行仿真，并與實驗結果進行對比，驗證了仿真結果的準確性。鄧峰等[2]利用半艙模型探討了風速、濕度、溫度等參數對除霜、除霧性能的影響。趙林林等[3]利用實驗和數值模擬的方法，證實整車進行除霜性能穩態仿真可用半艙模型代替。李明等[4]利用k-ε湍流模型研究了出風道結構特征、出風口沖擊角和溫升規律對除霜效果的影響。KUMAR 等[5]優化了前擋風玻璃表面風速分布，證實了風速分布對融霜速度有較大影響。JAHANI 等[6]證實除霜風道結構對前擋風玻璃表面風速分布有較大影響。谷正氣等[7]對除霜風道內的導流板進行參數優化，消除了前擋風玻璃表面的吹風死角，滿足了除霜要求。雖然國內外學者對汽車除霜系統進行了大量研究，證實k-ε湍流模型能準確模擬汽車除霜性能，除霜風道結構會影響前擋風玻璃表面風速分布，從而影響除霜速度，但由于汽車除霜風道形狀往往比較復雜，且不易選取參數進行優化，因此，尋求除霜風道的快速參數優化方法非常重要。自20 世紀80 年代以來，SEDERBERG 等[8]提出自由變形技術(free form deformation，FFD)，由于FFD 技術可對復雜幾何進行簡單、高效參數化，因而在航空、汽車等氣動優化設計領域得到廣泛應用。SAMAREH[9]證實了自由變形技術在氣動優化中能減少參數數量，是一種有效的氣動形狀優化方法。LEE等[10]將FFD技術和B 樣條曲面控制技術應用于3 種優化問題，發現FFD 技術更適合復雜形狀的參數優化。王丹等[11]驗證了FFD 技術對飛機快速參數化的應用效果，使飛機的氣動阻力顯著減小。汪怡平等[12-13]利用自由變形技術對汽車進行快速參數化，并建立近似模型，成功實現了汽車的氣動減阻優化。考慮到自由變形技術能對復雜幾何進行簡單、高效的參數化的優勢和除霜風道結構復雜、不易參數化的缺點，本文作者將自由變形技術應用到汽車除霜系統參數優化中。為了驗證自由變形技術在汽車除霜系統優化中的可行性，以某輕型載貨汽車的除霜系統作為優化對象；綜合考慮該車除霜的穩態和瞬態仿真結果后，選取駕駛員側除霜風道內的3塊導流板放入晶胞型控制框架中進行參數化；利用拉丁超立方試驗設計，搭建可信度高的徑向基模型，最后采用遺傳算法尋得導流板的最優參數，并比較優化前后的除霜性能。

1 汽車除霜性能的數值模擬

1.1 模型的建立與網格劃分

圖1 所示為乘員艙及玻璃分區示意圖。根據GB 11555—2009[14]將圖1 中的前擋風玻璃劃分為A，A′和B區。

圖1 乘員艙及玻璃分區示意圖Fig.1 Schematic diagram of the passenger compartment and glass partition

除霜系統內氣流流動復雜，因此，對除霜系統內網格進行加密處理；對空氣域采用網格質量好、求解穩定性好的六面體網格，玻璃域采用四棱柱網格進行離散；近壁面采用3 層四棱柱網格，總厚度為3 mm，網格厚度增長率為1.2。求解域中截面網格劃分如圖2所示。網格單元數約為629萬個，玻璃厚度為4 mm。

圖2 求解域中截面網格劃分Fig.2 Meshing of the section in the solution domain

1.2 邊界條件與物理條件設置

計算域入口設置為質量流入口，流量為280 m3/h；出口設置為壓力出口，為1×105Pa；空間離散格式為二階迎風格式，計算方法為SIMPLE算法，湍流模型為Realizablek-ε模型。在瞬態計算時，除霜風道入口送風溫度如圖3所示；設置霜層厚度為0.44 mm，環境溫度為255.13 K，設置玻璃、空氣和冰層物性如表1 所示；步進為1 s，最大內部迭代為5步，計算時長為2 400 s。

1.3 初始模型除霜效果分析

圖3 送風口溫升曲線Fig.3 Temperature rise curve of air inlet

表1 材料物理屬性Table 1 Physical properties of materials

對該款輕型載貨汽車模型進行除霜穩態分析，計算所得前擋風玻璃內表面的速度場如圖4 所示。由圖4 可見：在1.5 m/s 和2.0 m/s 以上風速覆蓋區中，A區的覆蓋面積比A′區的覆蓋面積小。圖5所示為不同時刻前擋風玻璃冰層厚度云圖。由圖5可見：融霜速度受風速分布影響較大，即A 區融霜速度比A′區的慢，且在A 區右下側前擋風玻璃融霜速度最慢。

圖4 擋風玻璃內表面速度分布云圖Fig.4 Cloud image of velocity distribution on windshield

根據 GB 11555—2009 規定，A 區、A′區和 B區均滿足除霜要求，但A′區風速分布比A 區的均勻，這使得A 區除霜速度明顯比A′的區慢，且駕駛員正對A 區，使A 區遮擋駕駛員視線的時間較長，不利于行車安全。因此，有必要提高A 區除霜速度，使駕駛員快速獲得清晰的視野，從而提高行車安全性。

2 基于自由變形技術的除霜優化模型

在汽車暖通系統鼓風機功率不變的情況下，A區融霜速度主要受擋風玻璃表面暖風氣流分布情況的影響，而在不改變除霜風道的外形條件下，氣流的分布主要由除霜風道內導流板的幾何形狀決定[4-7]，合理布置駕駛員側除霜風道內的導流板幾何形狀對A 區的除霜性能非常重要，因此，選取駕駛員側除霜風道內主要影響A 區流場的J1，J2和J3導流板進行變形優化，其在風道內的布置如圖6所示。

圖5 不同時刻前擋風玻璃冰層厚度云圖Fig.5 Cloud image of windshield’s ice thickness at different time

圖6 除霜風道內導流板安裝位置示意圖Fig.6 Schematic diagram of installation position of inner deflector in defrosting duct

由于導流板J1，J2和J3幾何形狀較復雜，且難于多參數優化，現引進自由變形技術對導流板J1，J2和J3進行快速參數化，盡量減少參數數量。

2.1 自由變形技術

自由變形技術FFD 采用Bernstein 基函數來建立晶格控制點與晶胞型控制框架內任意點位置之間的函數關系[8]，其表達式為

式中：(s,t,u)為晶胞型控制框架中任一點的局部坐標；Xffd(s,t,u)為對應局部坐標點(s,t,u)變形后的全局坐標；x(s,t,u)和Δx(s,t,u)為對應局部坐標點(s,t,u)的全局坐標和位移量；l，m和n為在局部坐標軸S，T和U上劃分的晶格控制點數量；Pi,j,k和ΔPi,j,k為某晶格控制點變形前的全局坐標和位移。

式中：Bil(s)為第i個l次Bernstein多項式基函數。

以1個半徑為10 mm的球體為例，將該球體嵌入3×3×3 個晶格控制點所組成的晶胞型控制框架中，通過移動頂部中心的晶格控制點5 mm后獲得如圖7所示結果。該自由變形過程主要分嵌入和變形2步。

1) 將球體模型嵌入3×3×3 個晶格控制點所組成的晶胞型控制框架中，設置局部坐標系X0STU。計算球體上各點的局部坐標(s,t,u)和對應點全局坐標x(s,t,u)。

2)拖動晶胞型控制框架頂部中心晶格控制點，移動ΔP1,1,2，并通過式(1)計算球體模型變形后各點全局坐標Xffd(s,t,u)，根據球體模型變形后的全部全局坐標最終獲得圖7所示結果。

圖7 球體自由變形示意圖Fig.7 Schematic diagram of free deformation of sphere

2.2 設計變量與優化目標

在除霜風道中提取出導流板J1，J2和J3的幾何數據，并利用FFD技術將提取的3塊導流板的幾何數據放入如圖8所示的晶胞型控制框架中，該晶胞型控制框架為根據導流板幾何形狀構建的異形六面體。

如圖8(b)所示，選取5個位置，以每個位置的4 個晶格控制點作為1 個變形控制點；將位置1，2，3，4 和 5 的各變形控制點坐標X1，X2，X3，X4和X5和位置點2 的變形控制點坐標Z2作為設計變量，各設計變量初始值為0，取值范圍如表2所示。將剩余控制點固定不動。

圖8 導流板及晶格與參數示意圖Fig.8 Deflector,lattice and parameter diagrams

表2 設計變量取值范圍Table 2 Value range of design variables

由于擋風玻璃的除霜性能可由擋風玻璃表面的流場分布預測，而在工程實際應用中，可通過分析A 和A′區80%的區域是否被2.0 m/s 以上的風速覆蓋以及B 區80%的區域是否被1.5 m/s 以上的風速覆蓋來評估除霜性能是否達到除霜要求[15]，因此，以A和A′區表面2.0 m/s和2.5 m/s以上的風速覆蓋比SA,2.0，SA,2.5，SA′,2.0和SA′,2.5與 B 區表面 1.5 m/s和2.0 m/s 以上的風速覆蓋比SB,1.5和SB,2.0作為目標函數。

3 試驗設計與近似模型

3.1 試驗設計

考慮到設計變量和目標函數較多，在滿足試驗設計要求的前提下，應盡量少采樣。而拉丁超立方方法可以大幅降低多因子、多水平試驗設計的采樣次數[12,16]，因此，選用拉丁超立方方法進行試驗設計。

試驗因子為變形控制點的位置參數，取值范圍見表2。使用拉丁超立方設計進行42組采樣，采樣結果如表3 所示。根據表3，移動變形控制點的相應位置，獲得42組導流板幾何數據。

考慮到乘員艙(全艙模型)體積較大，且全艙模型后部回流的氣流對前擋風玻璃內表面的流場分布影響甚微，且截取B 柱前部的乘員艙模型(半艙模型)和全艙模型的除霜穩態結果幾乎一致[3]，因此，為節省計算時間，將該乘員艙從B 柱中間截斷，保留B柱前部的半艙模型，如圖9所示，并將變形后的導流板J1，J2和J3的幾何數據導入半艙模型中進行數值計算，獲得如圖10所示的42組目標函數值。

表3 試驗設計結果Table 3 Test design results

3.2 近似模型的建立與檢驗

由于設計變量和目標函數較多，采樣數據較少，近似模型關系較復雜。徑向基模型(radial basis function，RBF)是一種神經網絡模型，能很好地逼近復雜非線性函數，在同等條件下，尤其適用于較少樣本點的情況，具有近似模型質量好、擬合速度快和精度高等優點[17]，因此，選用RBF模型作為導流板優化的近似模型。

為檢驗RBF 模型擬合精度，借助決定系數R2進行評估，R2越大，擬合精度越高[13]。經過計算，RBF 模型的R2=1，擬合精度滿足要求(R2＞0.9)。為更直觀地檢驗RBF 模型的擬合精度，另外選取除42組采樣點外的3組新采樣點作為檢驗點，3組檢驗點集如表4所示。根據3組檢驗點對導流板進行自由變形，并將變形后的導流板幾何數據導入半艙模型中進行數值模擬，得到的仿真結果與RBF模型中對應的結果進行比較，對比結果如表5 所示。由表5 可知：3 組檢驗點的仿真結果與對應RBF 模型結果的相對誤差均在±5%以內，擬合精度高，因此，選擇RBF 模型可以作為導流板優化的近似模型。

3.3 近似模型的求解

除霜性能的優化目標是在保證A′區除霜效果的同時，提高A區和B區的除霜效果，屬于典型的多目標優化問題[18]。針對多目標問題，通常是將其轉化為單目標問題進行尋優[19]，因此，設定單目標函數Fmax為

圖9 半艙模型示意圖Fig.9 Diagram of half passenger compartment model

圖10 采樣點目標函數值Fig.10 Objective function value of sample point

表4 檢驗點集Table 4 Checkpoint set

式中：ω1，ω2，ω3，ω4，ω5和ω6為權重系數；0.8＜SA,2.0＜1.0，0.8＜SA′,2.0＜1.0，0.8＜SB,1.5＜1.0。 根 據GB 11555—2009 規定，檢驗 A 區、A'區和 B 區滿足除霜性能要求的時間分別為20，25 和40 min，且各區表面各風速對除霜均有效果，因此，根據1 h 內各區除霜速度，令權重系數分別為0.21，0.21，0.18，0.18，0.11和0.11。

由于遺傳算法為群體搜索，利于搜出全局最優，且效率高[7]，因此，利用遺傳算法在建立的RBF 模型上進行全局尋優。為了得到準確的優化結果，將遺傳算法的初始種群數G按式(4)計算并設為106 個，最大內部迭代為50 步，最小內部迭代為25步，突變率為0.01。

式中：N為設計變量個數。

最終設計變量優化結果如表6所示。按優化后的設計參數對導流板進行自由變形，并與初始導流板J1，J2和J3進行對比，變形前后導流板對比示意圖如圖11 所示。由圖11 可知：優化后的導流板J1和J2較平順，彎曲弧度沒初始的大；導流板J1較初始的長，而導流板J3較初始的短。3塊導流板之間的最小間距增大。

表5 RBF模型精度檢驗Table 5 Accuracy test of RBF model

表6 設計變量優化結果Table 6 Optimization results of the designed variables

圖11 變形前后導流板示意圖Fig.11 Schematic diagrams of the deflector before and after deformation

將優化后的導流板J1，J2和J3的幾何數據導入半艙模型中進行數值模擬。將得到的半艙模型仿真結果與RBF 模型中對應的尋優結果進行比較，如表7所示。由表7可知：A區表面2.5 m/s以上風速覆蓋比的相對誤差最大，為-5.98%，其余區域風速覆蓋比的相對誤差均在±3%以內，說明所建RBF近似模型可信度高。

表7 風速覆蓋比仿真與尋優結果對比分析Table 7 Comparison and analysis of simulation and optimization results of wind speed coverage ratios

4 除霜效果分析

將優化后的導流板J1，J2和J3的幾何數據導入全艙模型中進行數值模擬，并與半艙模型的數值模擬結果進行對比，如表8所示。由表8可知：半艙模型與全艙模型數值模擬結果相似，誤差均在±5%以內，進一步說明對體積較大的乘員艙進行優化研究時，可用半艙模型進行研究。

表8 半艙和乘員艙模型各區風速覆蓋比Table 8 Wind speed coverage ratios in each section of the model of half cabin and crew cabin

4.1 穩態結果分析

為了更直觀地說明優化后模型的風速分布情況，將優化后的全艙模型與初始的全艙模型進行對比。圖12 所示為除霜風道內流線圖，表9 所示為優化前后各區風速覆蓋比。由圖12和表9可知：

1)由于優化后的導流板J1被拉長，阻礙了除霜風道內駕駛員側氣流流向副駕駛員側，使A′區2.0 m/s 和2.5 m/s 以上風速覆蓋比從初始的99.63%和90.41%分別下降到99.47%和86.96%；

2)與初始導流板J2相比，優化后的導流板J2彎曲弧度減小，更平滑。流經導流板J2附近區域的氣流更順暢、發散，使A 區2.0 m/s 和2.5 m/s 以上風速覆蓋比從初始的82.79%和64.09%提升到99.34%和86.91%，提升比例均在19%以上。

圖12 除霜風道內流線圖Fig.12 Internal flow diagram of defrosting duct

表9 優化前后各區風速覆蓋比Table 9 Wind speed coverage ratios of each district before and after optimization

3) 優化后的導流板J3更短，且在導流板J3上部，靠近前擋風玻璃側，無渦流作用，使氣流流動更順暢；而且導流板J3與上部導流板J2和J4協同作用，加速氣流流動。這使B 區的1.5 m/s 和2.0 m/s 風速覆蓋比均提升4.5%以上，從初始的93.98%和79.11%提升到98.48%和84.81%。

4.2 瞬態結果分析

圖13 所示為前擋風玻璃不同時刻各區域除霜覆蓋比，由圖13 可見：由于優化后的A′區表面各風速覆蓋比相比初始模型的略有下降，但A區和B區表面風速分布得到明顯改善，增強了A區和B區的對流換熱能力，使A 區除霜面積達到80%和100%的時間分別比初始模型提前17 s和40 s；B區除霜面積達到95%的時間比初始模型提前6 s；而A′區除霜面積達到80%和100%的時間僅比原模型分別延后了5 s和7 s。對比結果表明，優化后前擋風玻璃表面的除霜總時間縮短，前擋風玻璃A 區融霜速度比A′區的快。因此，對除霜風道內的導流板進行優化，能有效提高汽車的除霜性能。

圖13 前擋風玻璃不同時刻各區域除霜覆蓋比Fig.13 Defrosting coverage ratio of each area in front windshield at different time

5 結論

1)對某輕型載貨汽車除霜系統進行數值模擬，發現A 區表面流場分布沒有A′區的均勻。為此，本文作者以除霜風道內的3 塊導流板為優化對象，利用自由變形技術對其進行快速參數化，減少了優化的參數數量。

2)通過抽樣的方法搭建了高精度的近似模型，并利用遺傳算法尋優得到了最優結果。優化后A區2.0 m/s 和2.5 m/s 以上風速覆蓋比均提高19%以上，A區除霜面積達到80%和100%的時間分別提前 17 s 和 40 s；使 B 區 1.5 m/s 和 2.0 m/s 風速覆蓋比均提升4.5%以上，B 區除霜面積達到95%的時間提前6 s，有效改善了前擋風玻璃表面風速分布和對流換熱能力。

3)自由變形技術能簡單、高效地對除霜風道內復雜形狀的導流板進行快速參數化，有效改善汽車除霜系統性能。