Von Karman模型三維大氣紊流仿真理論與方法

高 靜 洪冠新 梁灶清

(北京航空航天大學 航空科學與工程學院，北京 100191)

Von Karman模型三維大氣紊流仿真理論與方法

高 靜 洪冠新 梁灶清

(北京航空航天大學 航空科學與工程學院，北京 100191)

基于三維紊流場相關函數矩陣，用蒙特卡羅法仿真生成Von Karman模型三維空間大氣紊流場，實現了對復雜的Von Karman模型的精確推導和應用，仿真流場的相關特性在理論上可以任意逼近Von Karman理論模型.通過提高白噪聲序列的白化程度及優化蒙特卡羅算法的方法提高了仿真精度和效率.數值仿真結果證明:所生成的空間大氣紊流場符合模型特性，具有較好的統計特性和仿真精度，滿足實時飛行仿真對紊流場數據的要求.

大氣紊流;蒙特卡羅法;Von Karman模型;數值仿真;飛行仿真

大氣紊流對在大氣層中飛行的各種飛行器，包括飛機、導彈、運載火箭、空天飛機以及微小型飛行器、氣球、飛艇等的設計、使用安全等都有較大影響.大氣紊流模型中，與 Dryden模型相比，Von Karman模型的頻譜函數在高頻段的斜率不同，頻譜函數更為復雜，更符合大氣紊流的實際情況［1］.在研究涉及飛機結構振動的飛行品質、飛機的結構疲勞等問題時，由于飛機結構模態頻率通常恰好處在高頻范圍，因而高頻范圍內紊流可能激發飛機的結構振動，只要可行，最好使用Von Karman大氣紊流模型.

對大氣紊流進行建模與仿真研究中，一維大氣紊流序列的仿真一般采用成型濾波器法，此方法要求模型的頻譜函數能夠進行共軛分解，因此無法基于Von Karman模型進行時域仿真，工程應用較少.但在分析空中加油、編隊飛行等多機飛行仿真及飛行品質等研究中，必須使用模型的二維或三維紊流場，且最好采用基于Von Karman模型.文獻［2－3］中生成二維紊流場的方法，在理論上有缺陷，各向同性不滿足.文獻［4］對 Von Karman模型進行有理化逼近，簡化了 Von Karman模型，生成的并不是真正的空間三維大氣紊流場.文獻［5］中提出的空間三維大氣紊流場的理論模型和數值方法，擺脫了以上仿真模型的局限，該方法不用對模型進行簡化，生成的三維大氣紊流場具有良好的均勻性和各向同性.但是文獻［5］只是基于 Dryden模型，沒有實現基于 Von Karman模型的三維紊流場仿真，且由于需要大樣本空間的計算以提高精度，此方法受計算機條件限制，仿真精度和效率都難以令人滿意，無法滿足工程計算需要.

本文應用蒙特卡羅法基于相關函數矩陣，仿真生成了Von Karman模型三維空間大氣紊流場，驗證了該方法的正確性，且在算法上進一步完善，仿真生成了精度較高的可應用飛行動力學仿真的三維大氣紊流場.

1 Von Karman模型

Von Karman模型的能量頻譜函數為

其中，a=1.33;L為紊流尺度;Ω為空間頻率;σ為紊流強度.

由此可推得大氣紊流Von Karman模型的縱向和橫向相關函數如下:

其中，ζ=ξ/(aL);Γ為 Gamma函數;K為 Bessel函數.由于Von Karman模型函數形式復雜，使得常規仿真方法受限.

根據Bachelor公式，得出基于Von Karman頻譜的空間相關函數如下:

基于式(4)可構造三維空間大氣紊流場的空間相關函數矩陣.

2 三維大氣紊流場建模理論

2.1 產生三維高斯白噪聲場并將其離散化

用數學庫函數產生(0，1)上高斯分布白噪聲序列{W}，隨機變量為{W1，W2，…，Wmnl}.

白噪聲的相關函數等于δ函數:

將隨機變量Wi按逐行排列方式排布在劃分為 m ×n×r的空間網格上，m，n，r分別為沿 X，Y，Z軸方向的網格點數.由此生成了如圖1所示的高斯分布三維白噪聲場{Z(x，y，z)}.其中 x=1，2，…，m，y=1，2，…，n，z=1，2，…，r，且 Z(i，j，k)=Wp，其中p=(k－1)×m ×n+(j－1)×m+i.

2.2 離散化三維空間隨機場的相關函數矩陣

構造一般相關系數矩陣MU={Cij}={R(ij)}，按照圖1所示，空間劃分為p=m×n×r個網格，則點點相關構造相關系數矩陣的一般表達式為

圖1 隨機變量空間分布示意圖

矩陣中Cij為由式(4)得到的基于Von Karman模型的空間相關函數.

2.3 產生滿足相關函數矩陣的紊流場

由于MU為對稱正定陣，采用Cholesky分解法，可將其分解為一個上三角陣和一個下三角陣的乘積，即MU=AAT.其中，A為下三角陣.設U是需要產生的隨機場，取U=AZ，其中，Z=［Z1，Z2，…，Zn］T為白噪聲隨機變量.

可以證明U是滿足相關系數矩陣要求的隨機場.

以上三維大氣紊流場建模方法，理論完備，可以精確生成任意需要的三維大氣紊流場.在實際數值仿真計算中，仿真精度取決于輸入的白噪聲精度以及紊流的樣本空間大小.可以證明，白噪聲的白化程度越高，紊流場的樣本空間越大，仿真精度就越高，表現在仿真結果上就是仿真曲線和理論曲線貼合得越好.

3 仿真算法

本文采用蒙特卡羅法對三維大氣紊流場進行數值仿真.蒙特卡羅法可以對復雜的隨機性問題進行求解，方法及程序結構簡單，大量簡單重復抽樣，做出的解答能否反映實際，與原來所建立的數學模型是否正確及輸入信息的質量有關，而與本方法無關.由于選用的樣本容量有限，因此總會存在誤差，提高解的精度則需要成百倍的增長樣本的容量.隨著計算機的迅速發展，將蒙特卡羅隨機模擬的方法應用于研究三維大氣紊流場數值仿真中，使得對空間三維大氣紊流速度場的高精度仿真得以實現.

3.1 改進的雙隨機交換法提高輸入白噪聲精度

由于仿真計算中采用數學庫函數生成的偽隨機數存在固有缺陷，導致對紊流場的仿真精度存在很大影響［6］.為了提高偽隨機數的白化程度，本文采用改進的雙隨機交換法.主要步驟如下:

1)隨機抽取序列x(i－1)(n)的兩點數據并互換其位置，得到序列xi(n)，重復執行N次;

2)計算xi(n)的自相關函數 ri和平方和SSi;

3)若SSi＜ε或i達到預定執行運算的最大循環次數Nmax，則停止運行;

4)若 SSi＜ SSi－1，則返回步驟 1)，繼續執行運算，反之，放棄在步驟1)中進行的第i次隨機雙換，并返回步驟1)繼續執行.

5)重復執行步驟1)～步驟4)j次，其中Nj=N/10j，若 Nj＜1，則令 Nj=1.

該方法只是對隨機數序列進行重新組合，因此不改變隨機數樣本的方差、均值以及概率密度特性.由于步驟2)、步驟3)的計算量很大，上述方法用較小的計算量獲得了較大的白化效果.例如處理一個6 000點序列耗時僅幾分鐘，隨機序列白化程度明顯提高.

3.2 相關函數矩陣的計算

Von Karman模型的三維空間大氣紊流場計算，是對計算機速度和內存要求較高的海量數據計算.在具有足夠大樣本容量的計算條件下，計算結果的統計特性將更接近理論值，數據結果的精度也就更高.在計算相關函數矩陣，并Cholesky分解為下三角矩陣的環節，將占用大量計算機內存并耗費大量時間.需采用通過降低程序計算量和優化內存使用的優化算法來提高計算效率，使得在相同的計算條件下可以計算盡可能大的樣本空間.

以10×10×10的空間隨機場為例，構造的基于Von Karman頻譜的縱向相關系數矩陣MU為

其中

C1，1000表示網格空間中第1點與第1000個點間的相關函數，其余類推;Δh是網格間距步長，為簡化表達式，直接按照式(4)計算出了系數;L是大氣紊流縱向尺度.

以上Von Karman模型相關函數表達形式雖然復雜，但是由于本建模方法無需對模型表達式進行簡化，因此可以進行數值計算.

3.3 仿真大氣紊流場的特性檢驗

三維大氣紊流場仿真計算結果是否合理，精度是否滿足要求，需要對仿真結果進行相關特性檢驗.檢驗的準則就是看仿真紊流序列的相關特性是否符合該紊流模型的相關函數理論表達式［7］.Von Karman模型的縱向和橫向相關函數的理論表達式見式(2)和式(3).仿真序列的相關函數如下式:

將仿真的相關函數曲線和理論相關函數曲線畫在同一幅圖上，來衡量它們的貼合程度.

3.4 抽樣-檢驗方法

蒙特卡羅方法提高解的精度需要成百倍地增加樣本的容量，而增加樣本容量受內存容量和CPU計算速度的限制，這是蒙特卡羅方法通過增加采樣點數提高精度的主要限制因素.另一方面，結果的誤差難以估計，更是蒙特卡羅法固有的特點.注意到蒙特卡羅法是一個純粹的抽樣方法，本文提出了“抽樣-檢驗方法”:抽樣方法采用蒙特卡羅法，檢驗則是對抽樣結果的誤差進行統計計算，進而篩選誤差達到要求的仿真紊流場.本文對檢驗長度分段逐段檢驗來實現，對每段分別設定誤差，以實現精細的誤差控制，從而得到滿足相應要求的仿真結果.按照抽樣-檢驗方法的要求改進后的算法，可以在短時間內完成同一尺度的成千上萬個紊流場計算，按照數理統計理論的大數定律，可以從中挑選到相關特性與Von Karman模型吻合度很高的紊流場，這樣就不需要成百倍地增加樣本的容量而得到更高精度的仿真結果.圖2為仿真算法流程示意圖.

圖2 仿真算法流程示意

4 仿真算例和分析

選取兩個算例，用本文方法計算基于Von Karman模型的三維紊流場.這兩個算例的紊流場尺度，在常見的紊流場的大尺度范圍內，具有典型的工程應用意義.選取紊流強度σ=1m/s，縱向紊流尺度為L=530m.根據反復計算調試，設定比較合理的誤差率控制條件，以達到精度和計算效率的雙贏.

算例1 網格點數:6000(1500×2×2).檢驗長度:100點，分3段(40點，30點，30點)逐段檢驗.U序列的控制誤差率為:第1段，0.05;第2段，0.05;第3 段，0.1.V，W 序列的控制誤差率為:第1段，0.1;第 2段，2;第 3段，1.計算時間:48min.相關性檢驗結果見圖3.

圖3 三維大氣紊流場(1500×2×2)的相關性檢驗

算例2 網格點數:8000(2×2000×2).檢驗長度:100點，分3段(40點，30點，30點)逐段檢驗.U序列的控制誤差率為:第1段，0.05;第2段，0.05;第3 段，0.1.V，W 序列的控制誤差率為:第1 段，0.1;第 2 段，1.5;第 3 段，0.75.計算時間:0.25min.相關性檢驗結果見圖4.

圖4 三維大氣紊流場(2×2000×2)的相關性檢驗

圖5是算例1的Von Karman模型紊流速度序列仿真結果，圖6是仿真紊流場(30×30)剖面示意圖.

圖5 Von Karman模型紊流速度序列仿真結果

本方法未對模型進行任何簡化，因此仿真生成的三維空間大氣紊流場完全符合所選用的Von Karman大氣紊流模型.由圖3～圖6的相關性檢驗圖、仿真紊流速度圖及紊流場剖面圖等表明，仿真相關性曲線與理論曲線吻合程度高，與以往其他研究結果相比，仿真精度有了大幅度提高，非常令人滿意.另外本文通過具體算法的改進，在計算條件、耗時和結果逼真度之間找到了一個比較理想的平衡方法，計算時間大大縮短，使仿真效率得到提高.

圖6 仿真三維空間大氣紊流場剖面示意圖

5 結論

通過大量計算分析表明，應用蒙特卡羅法基于相關函數矩陣，仿真生成Von Karman模型三維空間大氣紊流場的模型方法正確，仿真結果合理，計算精度和效率高.

本文的仿真理論，實現了對復雜的Von Karman模型的精確推導和應用，使得仿真結果在理論上可以無限逼近理論值.在仿真方法上，改進的雙隨機交換法提高了白噪聲的精度，優化了蒙特卡羅方法，采用抽樣-檢驗方法，實現了精細的誤差控制，從而解決了蒙特卡羅方法提高解的精度需要海量樣本容量和誤差控制的困難.本文提出的仿真方法，大大提高了仿真精度和效率.

本方法可直接應用于飛行仿真，使飛行仿真中的風場模擬更加真實，對模擬器的研制具有實際意義，完全符合工程應用需要.在實時飛行仿真中，可對白噪聲和確定紊流空間的相關函數矩陣進行預計算，以滿足實時仿真需要.如計算條件有限，又需要較大紊流場，可采用基于對稱原理的紊流場擴展方法，將生成的局部空間紊流場通過拼接的辦法有效擴展為大范圍連續紊流場.

References)

［1］ Etkin B.Turbulent wind and its effect on flight［J］.Journal of Aircraft，1981，18(5):327 －345

［2］肖業倫.用于飛行仿真的二維紊流場的數字生成法［J］.航空學報，1990，11(4):124 －129 Xiao Yelun.Digital generation of two-dimensional field of turbulence for flight simulation［J］.Acta Aeronautica et Astronautica Sinica，1990，11(4):124 －129(in Chinese)

［3］陸宇平，胡亞海.基于空間相關函數的二維紊流場數值生成法［J］.南京航空航天大學學報，1999，31(2):139－145 Lu Yuping，Hu Yahai.Digital generation of two-dimensional field of turbulence based on spatial correlation function［J］.Journal of Nanjing University of Aeronautics＆ Astronautics，1999，31(2):139－145(in Chinese)

［4］張峰，汪沛，王沖，等.基于Von Karman模型的三維大氣紊流仿真［J］.計算機仿真，2007，24(1):35 －38 Zhang Feng，Wang Pei，Wang Chong，et al.Simulation of threedimensional atmospheric turbulence based on Von Karmanmodel［J］.Computer Simulation，2007，24(1):35 － 38(in Chinese)

［5］洪冠新，肖業倫.蒙特卡羅法仿真生成三維空間大氣紊流場［J］.航空學報，2001，22(6):542 －545 Hong Guanxin，Xiao Yelun.Monte Carlo simulation for3-D-field of atmospheric turbulence［J］.Acta Aeronautica et Astronautica Sinica，2001，22(6):542 －545(in Chinese)

［6］馬樹峰，岳曉奎.大氣紊流的數字仿真［J］.計算機輔助工程，2008，17(3):54 －56 Ma Shufeng，Yue Xiaokui.Digital simulation on atmospheric turbulence［J］.Computer Aided Engineering，2008，17(3):54 －56(in Chinese)

［7］肖業倫.大氣擾動中的飛行原理［M］.北京:國防工業出版社，1993:166－173 Xiao Yelun.Flight theory in atmosphere turbulence［M］.Beijing:National Defense Industry Press，1993:166 － 173(in Chinese)

(編 輯:李 晶)

Theory and method of numerical simulation for 3D atmospheric turbulence field based on Von Karman model

Gao Jing Hong Guanxin Liang Zaoqing
(School of Aeronautic Science and Engineering，Beijing University of Aeronautics and Astronautics，Beijing 100191，China)

atmospheric turbulence;Monte Carlo method;Von Karman model;numerical simulation;flight simulation

V 321.2

A

1001-5965(2012)06-0736-05

2011-03-28;網絡出版時間:2012-04-01 11:59

http://www.cnki.net/kcms/detail/11.2625.V.20120401.1159.003.html

凡舟青年科研基金資助項目(20080502)

高 靜(1981－)，女，陜西子洲人，助理研究員，gaojing@buaa.edu.cn.

Three-dimension-field atmospheric turbulence based on Von Karman model was generated with Monte Carlo simulation，in view of correlation function matrix.According to the derivation of the Von Karman model for mulae，the correlation curves of the simulated field can highly approach the theory results.The precision and efficiency of atmospheric turbulence simulation were improved by using double stochastic interchange minimization algorithm to improve Gauss white noise sequence，and also by the optimization of Monte Carlo method.The simulated atmospheric turbulence field，is identical with the Von Karman model by numerical test，and has perfect statistical characteristic and simulation precision.The simulation method could be used in real-time flight simulation.