孫 巖,江 盟,孟德虹,黃 勇

(中國(guó)空氣動(dòng)力研究與發(fā)展中心 計(jì)算空氣動(dòng)力研究所,綿陽(yáng) 621000)

0 引 言

物面變形或運(yùn)動(dòng)下的非定常流動(dòng)數(shù)值模擬、氣動(dòng)外形優(yōu)化設(shè)計(jì)、靜/動(dòng)氣動(dòng)彈性計(jì)算中外形隨時(shí)間發(fā)生變化,需要對(duì)流場(chǎng)計(jì)算網(wǎng)格進(jìn)行更新以獲取新外形下的流場(chǎng),動(dòng)態(tài)網(wǎng)格生成技術(shù)是解決這一問(wèn)題的有效途徑[1]。動(dòng)態(tài)網(wǎng)格生成主要有兩種實(shí)現(xiàn)方式:網(wǎng)格變形與網(wǎng)格重構(gòu)。相比網(wǎng)格重構(gòu),網(wǎng)格變形具有算法簡(jiǎn)單、計(jì)算量小、自動(dòng)化程度高和不改變網(wǎng)格拓?fù)浼安灰胄碌碾x散誤差等優(yōu)點(diǎn),能夠減小網(wǎng)格變化對(duì)數(shù)值模擬結(jié)果的影響,在氣動(dòng)優(yōu)化和流/固耦合模擬中被廣泛采用[2]。但網(wǎng)格變形技術(shù)也存在一個(gè)明顯的缺點(diǎn),即在變形后網(wǎng)格質(zhì)量會(huì)降低,尤其是大變形情況下,可能會(huì)出現(xiàn)網(wǎng)格單元交叉或負(fù)體積單元的現(xiàn)象,導(dǎo)致流場(chǎng)模擬失敗。另一方面,在非定常計(jì)算中,網(wǎng)格變形方法會(huì)被反復(fù)調(diào)用,網(wǎng)格更新的時(shí)間消耗是影響流動(dòng)模擬效率的關(guān)鍵因素。因此,發(fā)展高效/魯棒/高質(zhì)量的網(wǎng)格變形算法一直是動(dòng)態(tài)網(wǎng)格生成領(lǐng)域的研究熱點(diǎn)[3]。

近十年,一些優(yōu)秀的網(wǎng)格變形算法相繼被提出,如Delaunay背景映射[4]、徑向基插值(Radial Basis Functions,RBF)[5-9]、反向距離加權(quán)(Inverse Distance Weight,IDW)[10]等代數(shù)方法,逐步取代早期采用基于偏微分方程求解的方法,成為目前網(wǎng)格變形采用的主要思路。而一些綜合兩種或多種方法的復(fù)合型網(wǎng)格變形策略也展現(xiàn)出更好的網(wǎng)格變形能力和更新效率,例如RBF-Delaunay[11]、RBFs-MSA[12]、RBF-TFI[13]等,其中RBF-TFI方法因其優(yōu)異網(wǎng)格變形能力在多種問(wèn)題中得到應(yīng)用[14]。

盡管網(wǎng)格變形方法逐步趨向成熟和多樣化,已經(jīng)能夠解決十分復(fù)雜的工程問(wèn)題[15],但目前的算法實(shí)現(xiàn)主要以研究型代碼為主,尚未有高效/魯棒/易用的模塊化網(wǎng)格變形程序供學(xué)術(shù)研究和工程問(wèn)題模擬使用。

國(guó)家數(shù)值風(fēng)洞(National Numerical Windtunnel,NNW)項(xiàng)目是中國(guó)政府在2018年底啟動(dòng)的大型自主軟件研發(fā)項(xiàng)目,致力于發(fā)展以計(jì)算流體動(dòng)力學(xué)(Computational Fluid Dynamics,CFD)為核心的流體數(shù)值模擬軟件群,構(gòu)建CFD 產(chǎn)業(yè)應(yīng)用生態(tài)圈,通過(guò)以點(diǎn)帶面,推動(dòng)中國(guó)計(jì)算機(jī)輔助設(shè)計(jì)行業(yè)的發(fā)展,助力未來(lái)中國(guó)制造業(yè)向著智能化、數(shù)字化方向升級(jí)。作為NNW 項(xiàng)目共性基礎(chǔ)庫(kù)的重要組成部分,網(wǎng)格變形技術(shù)模塊在吸收前期課題組研究代碼的基礎(chǔ)上,通過(guò)邏輯梳理、代碼重構(gòu)、接口規(guī)范、功能擴(kuò)展和算例測(cè)試,首先實(shí)現(xiàn)了多塊結(jié)構(gòu)網(wǎng)格變形程序(Structured Grid Deformation Program,SGDP)V1.0 版 本 的 開發(fā)。目前,SGDP V1.0正在NNW 項(xiàng)目?jī)?nèi)部開展測(cè)試,待成熟完善后,將在全國(guó)范圍內(nèi)推廣使用,作為學(xué)術(shù)研究或工程問(wèn)題模擬的有力工具逐步應(yīng)用于氣動(dòng)外形優(yōu)化、流/固耦合模擬等領(lǐng)域。

本文的目的是通過(guò)對(duì)SGDP V1.0程序的理論方法、組成架構(gòu)、主要功能、運(yùn)行方法、運(yùn)行性能和工程應(yīng)用等進(jìn)行簡(jiǎn)要介紹,讓使用者對(duì)SGDP 程序有更深入的認(rèn)識(shí),促進(jìn)用戶之間進(jìn)行交流和發(fā)現(xiàn)問(wèn)題,推動(dòng)SGDP的持續(xù)良性發(fā)展,為整個(gè)CFD 行業(yè)提供性能更加優(yōu)異的基礎(chǔ)程序庫(kù)。

1 SGDP V1.0簡(jiǎn)介

1.1 理論方法

結(jié)構(gòu)網(wǎng)格變形程序SGDP V1.0采用了徑向基函數(shù)(Radial Basis Functions,RBF)與超限插值(Transfinite Interpolation,TFI)結(jié)合的復(fù)合型動(dòng)態(tài)網(wǎng)格變形方法來(lái)對(duì)物面變形后的網(wǎng)格進(jìn)行更新。RBF方法能夠獲得高質(zhì)量的變形網(wǎng)格,但標(biāo)準(zhǔn)形式的RBF為了保證物面的精確變形,需要將所有物面網(wǎng)格點(diǎn)選做徑向基基點(diǎn),導(dǎo)致插值矩陣龐大、計(jì)算效率低和魯棒性差等問(wèn)題。TFI方法基于棱線的網(wǎng)格點(diǎn)分布特性即可重構(gòu)出新的面和體網(wǎng)格,具有良好的計(jì)算效率和網(wǎng)格質(zhì)量,但TFI無(wú)法對(duì)物面變形后的網(wǎng)格棱線進(jìn)行更新。復(fù)合方法RBF-TFI利用RBF對(duì)物面變形后的空間棱線坐標(biāo)進(jìn)行更新,然后利用TFI方法重構(gòu)得到變形后的空間面網(wǎng)格和體網(wǎng)格。RBF-TFI方法中,變形后的物面網(wǎng)格直接利用物面網(wǎng)格點(diǎn)位移進(jìn)行更新得到,物面的精準(zhǔn)變形得到保證,從而RBF可以選擇少量的徑向基基點(diǎn),能夠極大改善網(wǎng)格變形的效率和魯棒性。RBF-TFI方法具體的網(wǎng)格更新流程如圖1所示,其中為了保證物面和空間棱線交叉位置的一致性,選擇基點(diǎn)時(shí)需要將空間棱線與物面連接的點(diǎn)選做徑向基基點(diǎn)?；c(diǎn)選擇算法首先選擇所有物面網(wǎng)格窗口的角點(diǎn),然后調(diào)用重復(fù)點(diǎn)查找模塊去除重復(fù)的角點(diǎn),并通過(guò)“貪婪”算法選擇其他物面網(wǎng)格點(diǎn)作為剩余的徑向基基點(diǎn)。

圖1 RBF-TFI方法網(wǎng)格更新流程Fig.1 Grid updating procedure of RBF-TFI method

SGDP V1.0中采用了不帶多項(xiàng)式的RBF,可實(shí)現(xiàn)遠(yuǎn)場(chǎng)網(wǎng)格保持靜止,并對(duì)網(wǎng)格塊的棱線進(jìn)行坐標(biāo)更新。

徑向基函數(shù)插值可以表示為:

式中:d 為網(wǎng)格點(diǎn)的位移,x 為網(wǎng)格點(diǎn)的坐標(biāo)矢量,N為徑向基插值基點(diǎn)的數(shù)量,αi(i＝1,2,…,N)為插值系數(shù),φ 為基函數(shù),‖x－xi‖為兩點(diǎn)之間的距離。

兩點(diǎn)之間距離‖x－xi‖可以表示為:

已知N 個(gè)插值基點(diǎn)位置xi(i＝1,2,…,N)的位移di,利用式(1)可以構(gòu)建關(guān)于N 個(gè)關(guān)于插值系數(shù)ai的線性方程組

式中

求解方程組(3)得到插值系數(shù)ai,然后利用式(1)可以插值得到任意位置x 的位移d(x)。

變形網(wǎng)格質(zhì)量與基函數(shù)的類型有關(guān),根據(jù)大量網(wǎng)格算例測(cè)試結(jié)果,本文選擇網(wǎng)格質(zhì)量調(diào)控能力較好的Wendland’s C2函數(shù),具體表示為:

式中,η＝‖x－xi‖/R,其中R 為基函數(shù)的影響半徑,更多有關(guān)RBF的理論可以參考文獻(xiàn)[5]。

SGDP V1.0中采用了Soni等發(fā)展的標(biāo)準(zhǔn)形式TFI[16],能夠更好地處理復(fù)雜外形和極端分布下的網(wǎng)格構(gòu)造,有關(guān)TFI方法的詳細(xì)可以參考文獻(xiàn)[16]。

1.2 組成架構(gòu)

圖2給出了SGDP V1.0程序的文件架構(gòu),整個(gè)程序采用Fortran95語(yǔ)言編寫,其中main.f90定義了程序的主入口,然后通過(guò)launch.f90文件可以選擇調(diào)用不同的網(wǎng)格變形方法,目前僅定義了RBF-TFI方法,通過(guò)launch.f90 可以快速擴(kuò)展其他網(wǎng)格變形方法。文件rbf_tfi.f90定義了網(wǎng)格變形的流程和具體實(shí)現(xiàn),包括物面網(wǎng)格提取、物面網(wǎng)格位移計(jì)算、徑向基基點(diǎn)選擇、插值與系數(shù)矩陣計(jì)算、面/體網(wǎng)格更新等,grid.f90定義了網(wǎng)格數(shù)據(jù)結(jié)構(gòu)、網(wǎng)格數(shù)據(jù)、網(wǎng)格質(zhì)量分析和網(wǎng)格輸入/輸出等模塊,tfi.f90定義了二維/三維TFI的實(shí)現(xiàn)函數(shù),logo.f90定義了一些程序的logo和版權(quán)信息,distance.f90定義了不同距離計(jì)算的函數(shù),base.f90定義了矩陣運(yùn)算、字符串處理、Tecplot文件IO 和重復(fù)點(diǎn)查找等基礎(chǔ)模塊。

SGDP V1.0程序的功能架構(gòu)主要由網(wǎng)格文件輸入/輸出、網(wǎng)格變形、變形結(jié)果輸出三個(gè)大模塊組成。網(wǎng)格文件輸入/輸出模塊用于讀入和寫出結(jié)構(gòu)網(wǎng)格數(shù)據(jù)與邊界信息定義文件,目前能夠處理2種常用結(jié)構(gòu)網(wǎng)格格式:Plot3d 和Gridgen,3 種文件存儲(chǔ)格式:Binary、Unformatted、Ascii,2 種 數(shù) 據(jù) 精 度:Single Precision和Double Precision,組合共計(jì)12種不同的結(jié)構(gòu)網(wǎng)格數(shù)據(jù)文件格式,邊界信息文件格式采用Gridgen中Generic Solver定義的方式,使用Ascii格式進(jìn)行存儲(chǔ)。網(wǎng)格變形模塊根據(jù)不同的變形類型,計(jì)算物面網(wǎng)格點(diǎn)的位移運(yùn)動(dòng),然后插值得到空間棱線網(wǎng)格點(diǎn)的位移,并利用TFI算法對(duì)空間網(wǎng)格數(shù)據(jù)進(jìn)行構(gòu)造更新。變形結(jié)果輸出模塊將變形前后網(wǎng)格質(zhì)量分布、徑向基基點(diǎn)選點(diǎn)結(jié)果、基點(diǎn)插值誤差分布、變形外形等寫出到指定的文件,便于程序運(yùn)行后對(duì)網(wǎng)格變形效果進(jìn)行分析。

圖2 SGDP V1.0程序文件架構(gòu)Fig.2 File architecture of SGDP V1.0

1.3 主要功能

SGDP V1.0 程序目前支持4 種類型的網(wǎng)格變形:平動(dòng)、轉(zhuǎn)動(dòng)、彎曲和用戶自定義。平動(dòng)通過(guò)輸入?yún)?shù)文件中X、Y、Z 方向的位移進(jìn)行定義;轉(zhuǎn)動(dòng)為繞固定軸轉(zhuǎn)動(dòng),目前僅支持繞X、Y、Z 三個(gè)軸轉(zhuǎn)動(dòng),通過(guò)輸入?yún)?shù)文件中的轉(zhuǎn)動(dòng)軸、轉(zhuǎn)動(dòng)中心和轉(zhuǎn)動(dòng)角度定義;彎曲主要用于描述沿特定方向的彎曲變形,通過(guò)輸入?yún)?shù)文件中的最大彎曲變形、彎曲變形起始位置、彎曲變形最大位置進(jìn)行定義;用戶自定義通過(guò)讀入文件的方式定義物面網(wǎng)格點(diǎn)的位移,可以描述任意種類的變形,輸入文件的路徑、文件名在參數(shù)文件中指定。

在實(shí)際工程應(yīng)用中,單一的運(yùn)動(dòng)或變形方式非常少見(jiàn),因此平動(dòng)、轉(zhuǎn)動(dòng)、彎曲3種變形方式主要用來(lái)對(duì)程序進(jìn)行測(cè)試和演示算例制作,用戶自定義方式用來(lái)實(shí)現(xiàn)復(fù)雜變形或運(yùn)動(dòng)的輸入。例如,氣動(dòng)彈性計(jì)算或氣動(dòng)外形優(yōu)化設(shè)計(jì)中,流場(chǎng)物面網(wǎng)格點(diǎn)的位移文件可以通過(guò)結(jié)構(gòu)變形插值或者幾何參數(shù)化模塊輸出獲得,然后調(diào)用SGDP V1.0,即可得到物面變形后的流場(chǎng)計(jì)算網(wǎng)格。

此外,SGDP V1.0通過(guò)對(duì)更新前后的網(wǎng)格單元數(shù)據(jù)進(jìn)行分析,可以獲得每個(gè)網(wǎng)格塊的所有網(wǎng)格單元質(zhì)量,并統(tǒng)計(jì)得到每個(gè)網(wǎng)格塊的最差單元質(zhì)量(worst)、最優(yōu)單元質(zhì)量(best)和平均單元質(zhì)量(average),然后輸出變形前后網(wǎng)格單元質(zhì)量在不同塊中的分布特性。

1.4 運(yùn)行流程

SGDP V1.0程序目前可以在Windows和Linux平臺(tái)下運(yùn)行,均采用命令行和配置文件的方式啟動(dòng)。

在Window 下啟動(dòng)CMD 或在Linux 下啟動(dòng)Terminal,定位至工作目錄后,SGDP V1.0可以通過(guò)下面3種方式運(yùn)行:

1)SGDP

2)SGDP project_config

3)SGDP method project_config

其中,第1種運(yùn)行方式,默認(rèn)采用RBF-TFI作為網(wǎng)格變形方法,config.txt文件作為默認(rèn)參數(shù)配置文件;第2中運(yùn)行方式,RBF-TFI作為默認(rèn)網(wǎng)格變形方法,配置參數(shù)文件的名稱project_config可以用戶自己定義,例如SGDP naca0012.txt;第3種運(yùn)行方式,網(wǎng)格變形方法method和配置參數(shù)文件名稱project_config均可以用戶自己定義,例如SGDP RBF_TFI naca0012.txt。

SGDP V1.0程序運(yùn)行前,需要準(zhǔn)備好原始網(wǎng)格和邊界信息文件,并配置好輸入?yún)?shù)文件。SGDP V1.0借鑒了斯坦福非結(jié)構(gòu)流場(chǎng)解算軟件SU2的參數(shù)輸入模式,采用扁平化參數(shù)組織結(jié)構(gòu),參數(shù)讀入沒(méi)有先后順序,從而可以方便地添加或刪除任意的參數(shù),參數(shù)文件中沒(méi)有定義的參數(shù)將使用程序中的默認(rèn)定義值。圖3給出了某個(gè)網(wǎng)格變形算例的輸入?yún)?shù)文件示例,這里采用了用戶自定義的變形輸入方式,只對(duì)相關(guān)的參數(shù)進(jìn)行了定義。

圖3 SGDP V1.0輸入?yún)?shù)文件定義示例Fig.3 Input parameter file example of SGDP V1.0

SGDP V1.0程序啟動(dòng)后,會(huì)對(duì)讀入的部分特定參數(shù)進(jìn)行解析,并進(jìn)行參數(shù)安全檢查,然后將全部的參數(shù)信息輸出到屏幕,方便用戶了解程序運(yùn)行情況,及時(shí)發(fā)現(xiàn)參數(shù)輸入錯(cuò)誤的情況并進(jìn)行相應(yīng)修正。圖4給出了SGDP V1.0的全部運(yùn)行參數(shù)打印信息。有關(guān)每個(gè)參數(shù)的具體含義,可以參考SGDP V1.0的使用手冊(cè)。

2 性能測(cè)試

計(jì)算時(shí)間消耗是網(wǎng)格變形程序的關(guān)鍵性能指標(biāo),在重復(fù)調(diào)用程序的仿真工程作業(yè)中,對(duì)整個(gè)仿真的時(shí)間和計(jì)算成本有著重要的影響。本節(jié)利用DLR-F6翼身組合體構(gòu)型[17],采用不同規(guī)模的多塊結(jié)構(gòu)網(wǎng)格,對(duì)SGDP V1.0的運(yùn)行時(shí)間性能進(jìn)行測(cè)試。

測(cè)試中采用的不同規(guī)模網(wǎng)格,基于1套基礎(chǔ)網(wǎng)格數(shù)據(jù),利用課題組自主開發(fā)的多塊結(jié)構(gòu)網(wǎng)格自適應(yīng)稀疏/加密程序,成比例的加密得到,不同規(guī)模網(wǎng)格的拓?fù)浣Y(jié)構(gòu)和網(wǎng)格塊數(shù)量是不變的。基礎(chǔ)網(wǎng)格數(shù)據(jù)含有72 個(gè)網(wǎng)格塊,有3 868 996個(gè)網(wǎng)格點(diǎn)、3 566 336個(gè)六面體單元,圖5給出了DLR-F6翼身組合體外形和基礎(chǔ)網(wǎng)格的拓?fù)浣Y(jié)構(gòu)。

圖4 SGDP V1.0全部運(yùn)行參數(shù)輸出信息Fig.4 Printing information of all running parameters of SGDP V1.0

本次測(cè)試,網(wǎng)格加密系數(shù)在0.5～4.0之間變化,網(wǎng)格點(diǎn)數(shù)量從50萬(wàn)增長(zhǎng)到2.3億,選擇DLR-F6 模型繞質(zhì)心旋轉(zhuǎn)20°的算例,徑向基插值最大基點(diǎn)數(shù)目為300。程序在Windows 7系統(tǒng)上運(yùn)行,計(jì)算機(jī)硬件采用Intel i7-8700 CPU 處理器,含有12個(gè)處理器單元,主頻為3.2GHz,內(nèi)存為32Gb。SGDP V1.0目前只有串行版本,采用單CPU 方式運(yùn)行。

圖6給出了SGDP V1.0程序總的CPU 運(yùn)行時(shí)間隨網(wǎng)格點(diǎn)數(shù)目Nnode的變化規(guī)律。隨著Nnode的增加,網(wǎng)格變形消耗的CPU 時(shí)間在逐步增大,但增加的速率是減小的,即每百萬(wàn)網(wǎng)格點(diǎn)消耗的時(shí)間隨著網(wǎng)格點(diǎn)數(shù)的增加逐漸減小。

SGDP V1.0程序總的運(yùn)行時(shí)間主要包括4個(gè)部分:網(wǎng)格文件讀取、徑向基基點(diǎn)選擇、體網(wǎng)格更新和網(wǎng)格文件寫出。圖7給出了SGDP V1.0程序運(yùn)行時(shí)間分布隨網(wǎng)格節(jié)點(diǎn)數(shù)目的變化圖,基于“貪婪”法的徑向基基點(diǎn)選擇占據(jù)程序全部運(yùn)行時(shí)間的絕大部分,比例隨著網(wǎng)格點(diǎn)數(shù)據(jù)的增加逐步減小。這是由于貪婪選點(diǎn)算法是將物面網(wǎng)格點(diǎn)位移誤差最大的點(diǎn)逐步選擇出來(lái),算法的計(jì)算量正比于物面網(wǎng)格點(diǎn)數(shù)量。而隨著網(wǎng)格規(guī)模增加,物面網(wǎng)格點(diǎn)數(shù)量是按照加密比例平方增加,而體網(wǎng)格點(diǎn)數(shù)量是按照加密比例立方增加,所以物面網(wǎng)格點(diǎn)數(shù)量的增加速度要遠(yuǎn)小于體網(wǎng)格點(diǎn)的增加速度。網(wǎng)格更新時(shí)間、網(wǎng)格數(shù)據(jù)讀入/寫出時(shí)間占總運(yùn)行時(shí)間的比例隨著網(wǎng)格點(diǎn)數(shù)目逐步增加的,但網(wǎng)格讀入/寫出時(shí)間占總運(yùn)行時(shí)間的比例非常小,基本在1.0%以下。

圖5 DLR-F6翼身組合體外形及基礎(chǔ)網(wǎng)格拓?fù)浣Y(jié)構(gòu)Fig.5 Shape and basis grid topology of DLR-F6 wing-body configuration

圖6 SGDP V1.0總運(yùn)行時(shí)間隨網(wǎng)格節(jié)點(diǎn)數(shù)目的變化Fig.6 Total CPU time cost of SGDP V1.0 with different grid nodes

圖7 不同網(wǎng)格節(jié)點(diǎn)數(shù)目下SGDP V1.0運(yùn)行時(shí)間分布Fig.7 CPU time cost percent of SGDP V1.0 with different grid nodes

對(duì)于網(wǎng)格變形程序多次調(diào)用的算例,由于物面精確變形不是通過(guò)徑向基函數(shù)插值RBF來(lái)保證,所以基點(diǎn)貪婪選擇算法可以利用最大變形只執(zhí)行一次,然后將基點(diǎn)在物面網(wǎng)格點(diǎn)上的索引編號(hào)通過(guò)數(shù)組或文件記錄下來(lái),下一次程序調(diào)用直接使用此編號(hào)去獲取徑向基基點(diǎn)的坐標(biāo),然后更新系數(shù)和插值矩陣,將極大地改善程序的運(yùn)行效率。對(duì)于這類情況,例如氣動(dòng)彈性多次迭代或優(yōu)化中不同參數(shù)下氣動(dòng)外形變化,SGDP V1.0程序的運(yùn)行時(shí)間主要依賴于網(wǎng)格讀寫時(shí)間和體網(wǎng)格更新時(shí)間。而對(duì)于利用內(nèi)存進(jìn)行網(wǎng)格數(shù)據(jù)傳遞的例子,網(wǎng)格讀取時(shí)間也將是一次性的。

圖8給出了網(wǎng)格讀寫時(shí)間、網(wǎng)格變形更新時(shí)間隨總的網(wǎng)格點(diǎn)數(shù)目的變化規(guī)律?？梢钥闯?網(wǎng)格讀寫時(shí)間、網(wǎng)格更新時(shí)間隨網(wǎng)格點(diǎn)數(shù)目基本是線性變化的趨勢(shì),每百萬(wàn)網(wǎng)格點(diǎn)的讀取時(shí)間約為0.01 s,每百萬(wàn)網(wǎng)格點(diǎn)的寫出時(shí)間約為0.008 s,每百萬(wàn)網(wǎng)格點(diǎn)的更新時(shí)間約為0.1 s。因此,對(duì)于千萬(wàn)量級(jí)的網(wǎng)格,SGDP V1.0程序可以在數(shù)秒內(nèi)完成一次變形網(wǎng)格生成,綜合來(lái)說(shuō)是十分高效的。

圖8 網(wǎng)格讀寫時(shí)間、更新時(shí)間隨網(wǎng)格節(jié)點(diǎn)數(shù)目的變化Fig.8 CPU time cost of grid reading,writing and updating with different grid nodes

3 變形算例

本節(jié)提供一些典型情況下的網(wǎng)格變形算例,用來(lái)驗(yàn)證SGDP V1.0網(wǎng)格變形能力和網(wǎng)格質(zhì)量維持特性。

圖9給出了NACA0012翼型繞軸旋轉(zhuǎn)90°前后的變形網(wǎng)格,其中旋轉(zhuǎn)中心位于弦線1/4位置。旋轉(zhuǎn)情況下,變形向遠(yuǎn)場(chǎng)的衰減過(guò)程會(huì)使網(wǎng)格單元產(chǎn)生明顯的剪切變形,對(duì)網(wǎng)格單元質(zhì)量有著顯著的影響,如圖9(b)。

變形前,NACA0012網(wǎng)格最差單元質(zhì)量為0.800,變形后,NACA0012 網(wǎng)格最差單元質(zhì)量為0.597。90°旋轉(zhuǎn)后,盡管網(wǎng)格單元質(zhì)量有了顯著的降低,但依然維持較高的水平(＞0.30),說(shuō)明SGDP V1.0在顯著旋轉(zhuǎn)變形下仍然具有很好的網(wǎng)格質(zhì)量維持能力。從圖9(b)還可以看出,物面附近網(wǎng)格隨物面一起運(yùn)動(dòng),能夠保證物面附近黏性流場(chǎng)計(jì)算區(qū)域網(wǎng)格單元保持良好的連續(xù)性和正交性。

圖9 NACA0012旋轉(zhuǎn)90°的變形網(wǎng)格Fig.9 Deformed grid of NACA0012 airfoil with 90°rotation

圖10 為某大展弦比柔性機(jī)翼模型在氣動(dòng)載荷作用下變形前后的流場(chǎng)計(jì)算網(wǎng)格。該大展弦比機(jī)翼半展長(zhǎng)為1542,翼梢最大彎曲變形達(dá)到402,約為半展長(zhǎng)的26%。在如此大幅度變形下,機(jī)翼附近網(wǎng)格單元仍然能夠隨著機(jī)翼一起運(yùn)動(dòng),保證了變形后的網(wǎng)格質(zhì)量,說(shuō)明SGDP V1.0具有處理大變形情況的能力。

圖10 大展弦比柔性機(jī)翼模型變形前后網(wǎng)格剖面Fig.10 Grid section of high aspect-ratio flexible wing model before and after deformation

圖11 為NASA 共同研究模型翼/身/掛架/短艙組合體構(gòu)型CRM-WBPN 彎曲變形前后的流場(chǎng)計(jì)算網(wǎng)格。CRM-WBPN 模型計(jì)算網(wǎng)格有153 個(gè)網(wǎng)格塊和約1200萬(wàn)網(wǎng)格點(diǎn)(半模),網(wǎng)格拓?fù)湓跈C(jī)翼/掛架/短艙連接位置十分復(fù)雜。從圖11 可以看出,SGDP V1.0能夠很好地處理大彎曲變形下CRM 的網(wǎng)格變形,表明SGDP V1.0具有非常好的復(fù)雜外形網(wǎng)格變形處理能力。

圖11 NASA CRM-WBPN 模型變形前后網(wǎng)格剖面Fig.11 Grid section of NASA CRM-WBPN model before and after deformation

圖12 為NASA CRM-WBPN 大彎曲變形前后單個(gè)網(wǎng)格塊最差和最優(yōu)單元質(zhì)量隨網(wǎng)格塊的分布。變形前后,最差、最優(yōu)單元質(zhì)量基本沒(méi)有變化,依然具有較高的網(wǎng)格質(zhì)量,甚至某些塊變形后最差單元質(zhì)量還有所改善(第88塊),說(shuō)明SGDP V1.0在大變形情況具有十分優(yōu)異的網(wǎng)格質(zhì)量維持能力。

4 工程應(yīng)用

目前,結(jié)構(gòu)網(wǎng)格變形程序SGDP V1.0已經(jīng)被完整地集成到了亞跨超流場(chǎng)模擬平臺(tái)(Trisonic Platform,TRIP)中,用于開展不同類型的靜/動(dòng)氣動(dòng)彈性耦合數(shù)值模擬,如機(jī)翼載荷分配、副翼效率、風(fēng)洞模型結(jié)構(gòu)靜變形影響修正[18-19]、戰(zhàn)斗機(jī)垂尾抖振、機(jī)翼顫振[15]、空調(diào)管道流/固耦合振動(dòng)等。

圖12 NASA CRM-WBPN 模型變形前后最差、最優(yōu)網(wǎng)格單元質(zhì)量分布Fig.12 Worst and best grid cell quality distribution of NASA CRM-WBPN model before and after deformation

圖13 HIRENASD機(jī)翼變形計(jì)算結(jié)果Fig.13 Deformation calculation results of HIRENASD wing

圖13 給出了HIRENASD 機(jī)翼靜氣動(dòng)彈性計(jì)算的變形結(jié)果。HIRENASD 機(jī)翼是由德國(guó)亞琛工業(yè)大學(xué)設(shè)計(jì),用于研究跨聲速飛行雷諾數(shù)條件下氣動(dòng)彈性現(xiàn)象與機(jī)理的風(fēng)洞試驗(yàn)?zāi)Ｐ?曾被選為第一屆氣動(dòng)彈性預(yù)測(cè)研討會(huì)(AePW1)的標(biāo)準(zhǔn)算例,詳細(xì)信息見(jiàn)文獻(xiàn)[20]。本文算例的計(jì)算條件取自文獻(xiàn)[21],用于驗(yàn)證氣動(dòng)彈性計(jì)算程序的精準(zhǔn)度,具體計(jì)算參數(shù)為馬赫數(shù)Ma＝0.80、雷諾數(shù)Rec＝1.4×107、載荷因子q/E＝4.7×10－7、迎角α＝3°,其中載荷因子q/E 為風(fēng)洞來(lái)流速壓q 與風(fēng)洞模型材料彈性模量E 的比值。圖13中LE、TE分別表示機(jī)翼前后緣,TRIP表示采用SGDP V1.0的氣動(dòng)彈性計(jì)算模塊,TAU 為德國(guó)宇航院的CFD 計(jì)算軟件,EXP表示風(fēng)洞試驗(yàn)得到的結(jié)果。可以看出,采用SGDP V1.0的TRIP 獲得和TAU、EXP比較一致的彎曲和扭轉(zhuǎn)變形結(jié)果,在靠近翼梢位置TRIP 與TAU 計(jì)算的扭轉(zhuǎn)變形存在略微的差異,可能是TAU 僅采用了前4階模態(tài)導(dǎo)致扭轉(zhuǎn)變形計(jì)算截?cái)嗾`差增大引起的。

圖14 HIRENASD機(jī)翼壓力系數(shù)計(jì)算結(jié)果Fig.14 Pressure coefficients calculation results of HIRENASD wing

圖14 給出了HIRENASD 機(jī)翼3個(gè)展向剖面位置的壓力系數(shù)計(jì)算結(jié)果,圖中Rigid表示機(jī)翼剛性條件下獲得的壓力分布,Elastic表示采用SGDP V1.0計(jì)算的彈性狀態(tài)下的壓力分布?？梢钥闯?彈性狀態(tài)下的壓力分布與試驗(yàn)值十分吻合,在靠近翼梢位置與剛性結(jié)果差異明顯,側(cè)面驗(yàn)證了SGDP V1.0程序的可靠性,能夠提供高質(zhì)量的變形網(wǎng)格。

5 結(jié) 論

對(duì)國(guó)家數(shù)值風(fēng)洞工程N(yùn)NW 標(biāo)準(zhǔn)共性基礎(chǔ)庫(kù)的第1款結(jié)構(gòu)網(wǎng)格變形程序SGDP V1.0的特性進(jìn)行了簡(jiǎn)要的概述,并通過(guò)DLR-F6不同規(guī)模網(wǎng)格對(duì)SGDP V1.0的性能進(jìn)行了測(cè)試,以及開展了SGDP V1.0網(wǎng)格變形能力和網(wǎng)格質(zhì)量維持特性驗(yàn)證、靜氣動(dòng)彈性數(shù)值模擬等工作。從測(cè)試與計(jì)算結(jié)果可以得到:

1)SGDP V1.0擁有非常高的網(wǎng)格變形效率,可以在數(shù)秒內(nèi)完成千萬(wàn)量級(jí)結(jié)構(gòu)網(wǎng)格更新,能夠滿足目前絕大多數(shù)工程實(shí)際中的網(wǎng)格變形需求。

2)針對(duì)不同類型結(jié)構(gòu)網(wǎng)格變形問(wèn)題,SGDP V1.0均能獲得較高質(zhì)量的變形網(wǎng)格,具有良好的網(wǎng)格變形處理和質(zhì)量維持能力,程序展示了優(yōu)良的魯棒性。

3)SGDP V1.0采用了扁平化和模塊化的參數(shù)輸入方式,可以通過(guò)少量幾個(gè)輸入?yún)?shù)的修改應(yīng)用于不同類型的網(wǎng)格變形算例,展示了良好的通用性。

下一步SGDP V1.0將結(jié)合使用經(jīng)驗(yàn)和反饋意見(jiàn),對(duì)程序不斷地進(jìn)行迭代完善,并結(jié)合未來(lái)超大規(guī)模網(wǎng)格使用需求,開展并行版本的開發(fā)。SGDP V1.0將根據(jù)單位下一步的知識(shí)產(chǎn)權(quán)政策,在逐步完善后確定是否開源,如有使用需求,可以通過(guò)郵箱y.sun＠cardc.cn與作者聯(lián)系,研究團(tuán)隊(duì)將提供執(zhí)行程序版本。

致謝:感謝中國(guó)空氣動(dòng)力研究與發(fā)展中心洪俊武研究員在多塊網(wǎng)格自適應(yīng)稀疏/加密方面提供的建議和幫助。