基于本征正交分解的平流層風場建模與預測

李魁， 鄧小龍， 楊希祥， 侯中喜， 周新

(國防科技大學航天科學與工程學院, 長沙 410073)

臨近空間是介于傳統航空與航天區域之間、尚未為人類所開發利用大氣層，由于其獨特的資源優勢成為人們關注的焦點[1]。臨近空間底部的平流層具有氣流穩定，空氣流動相對緩慢的特點，特別是在一定時間內，存在風速較小的準零風層(Quasi-Zero Wind Layer，QZWL)[2]。準零風層一般是指平流層下層20 km高度附近的大氣層，上下層緯向風風向相反，同時經向風分量亦很小，是由于平流層風場中上下層緯向風逆轉形成的，是存在于平流層風場的一種特殊現象。該區域可用于部署長時間駐空低動態飛行器進行高分辨率對地觀測、通信中繼等任務。平流層風場對低動態臨近空間飛行器總體設計與飛行控制有重要影響[3]。

風場環境的建模與預測，需要對大量歷史數據進行處理和分析，有效提高建模與預測效率和精度值得深入研究[4]。本征正交分解(Proper Orthogonal Decomposition，POD)是風場建模中常用的一種方法，將風場數據通過特征值分解為本征模態和主坐標的組合，其中本征模態主要與空間因素相關，主坐標主要與時間因素相關，物理意義明確，通過模態截斷還能提高計算的效率[5]。

目前，POD方法已被廣泛地應用于解決實際問題數值模擬的降維問題中。Holmes等[6]采用氣動平均的方法對雙坡屋蓋的表面壓力進行測量，指出POD方法是描述鈍體風壓的有效途徑。Fic等[7]將POD方法的離散形式應用于線性和非線性瞬態熱傳導的降階模型中。Tan等[8]應用POD結合三次樣條插值方法對流場壓力系數進行預測。對于風場建模，胡亮[9]利用POD方法對具有橋塔風效應的風場進行了簡化模擬，將POD型譜表示法取代原型譜表示法，物理意義明確且可通過模態截斷節省計算量。陶青秋[10]將POD方法用于分析在大跨屋面和高層建筑物上的風致響應，利用POD與結構模態相結合的雙模態方法來分析結構的動態響應，表明POD方法是一種有效分析多激勵響應的方法。鄒垚等[11]通過對POD原理進行修正后用于風洞實驗風場的重建，主要通過修正第一個本征模態值，對偏差較大測點的風壓場重建進行改進，使重建風場的方差值與實際風場方差值的差距大大縮小。Chen和Kareem[12]利用POD方法對建筑結構的數值風場進行了模擬并提出了狀態空間建模方案，建立了風載荷與建筑結構共振響應在物理意義上的聯系，有助于更好地了解結構隨時空變化的動態負載響應。

本文基于適用于平流層風場環境的POD方法，以長沙地區風場數據為例，研究平流層風場隨高度和時間的變化特點，通過POD方法對風場數據進行降階處理，實現在特定時間范圍內快速準確計算出風場輪廓線。鑒于所選地區平流層20 km附近東西方向風速呈現明顯的空間相關性和周期性，在POD方法的基礎上，提出了一種可以對風場進行預測的Fourier模型。

1 基于POD方法的風場建模

POD方法是指從一組龐大的數據中獲得一組低維最優基，而這組低維最優基是一種能夠表示原復雜系統的降階模型。也就是說POD方法是一種降維的方法，對給定的數據進行最優的低維逼近，用較小的維數將原物理模型的主要特征表現出來[13]。本文通過POD方法獲得風場信息的最優標準正交基，將風場數據中的每個風場信息投影到標準正交基上，獲得每個風場數據在標準正交基上的投影系數[14]。并通過最優標準正交基與相對應的投影系數計算出一組新的數據，將原風場數據的主要特征表現出來。

1.1 瞬像矩陣的生成

瞬像也稱為快照，由Sirovich提出[15]，是指物理場的數值解在不同時刻的空間分布值。建立海拔高度h1至hL的風場模型,取瞬像矩陣為

(1)

式中：V為風速；下標L為瞬像時刻數或快照數。假設風場模型統計天數為M，則A∈RL×M，一般由試驗、經驗、統計等方法獲得，通常M?L。

1.2 POD模態的獲取

POD模態又稱為POD基。POD方法的目的是尋找一組最優的POD基φ1,φ2,…,φr(φi為列向量)，使得數據模型中任意一天特定時刻的風場可以表示為POD基的線性組合[16]，即

V(t)=α1(t)φ1+α2(t)φ2+…+αr(t)φr

(2)

式中：αi(t)為POD基的系數，i=1,2,…,r。為了尋求滿足式(2)的最優POD基，定義相關矩陣：

R=ATA

(3)

則R∈RL×L，求取相關矩陣R的非零特征值以及特征向量φi：

Rφi=λiφi

i=1,2,…,n;λ1≥λ2≥…≥λn>0;n≤L

(4)

將特征向量化為最優POD基：

(5)

且φi(i=1,2,…,n；n≤L)是一組標準正交基，即

(6)

只選取前r階模態，POD模態矩陣[17]記為

(7)

投影系數矩陣記為

(8)

則降階模型矩陣B=ΦΩ∈RL×M。

1.3 降階模型誤差評估

POD模態表示捕獲物理場的主要特征，截斷后的前r階POD模態較全階模態所捕獲的能量比為

(9)

故定義相對截斷誤差[18]為

(10)

2 基于Fourier級數的風場預測模型

通過第1節的POD方法夠獲得風場數據在最優POD基上的投影系數，若具有明顯的空間相關性和周期性，可以使用Fourier級數直接進行表示。通過最優POD基與Fourier級數所表示的投影系數計算出特定高度的風速隨時間變化情況[19]。其方程表示如下：

(11)

式中：t為時間，d；m為Fourier級數展開的項數；aki和bki為Fourier系數；T為投影系數的振蕩周期。

3 實例分析與結果

以長沙地區為例，采用2005年1月1日—2009年12月31日每天00：00所測量的風場數據，風場數據由國防科技大學氣象海洋學院提供，數據記錄了長沙地區00：00和12：00(世界時)的風速、濕度、溫度、氣壓等氣象要素，數據覆蓋海拔高度0～60 km之間每隔0.5 km所對應的天氣數據。考慮到臨近空間低動態飛行器駐空期間的工作高度范圍，本文針對海拔高度10～30 km的范圍進行研究，并通過3次樣條對風場數據進行了插值處理。

3.1 風場建模與預測

相對模態能量分布主要用于識別數據降階處理所需的模態數量。圖1(a)給出了東西方向風場各階POD模態的相對模態能量，前5種模態能量所占的比重相對于其他模態能量具有明顯的區別，可以捕獲到風場總特征性的98.9%。同時，由圖1(b)可知，南北方向風場的相對模態能量分布除前2種模態之外，其余模態能量呈逐漸下降的趨勢，若只選取前2種模態進行降階處理，不足以表現出風場的主要特征。因此，南北方向風場不能高效地進行降階建模。

圖1 相對模態能量分布Fig.1 Relative mode energy distribution

通過選用多種POD模型分別對2007年11月6日0：00時刻東西方向和南北方向的風場進行建模。由圖2(a)可知，10階和15階POD模型較好地表示出東西方向風場的實際數據，5階POD模型能體現出東西方向風場的主要特征。圖2(b)中實際風場模型清楚地表明南北方向風速隨空間的變化是極其不規則的，導致南北方向風場的低階重建性較差。因此，低階的POD模型(比如5階)在捕獲東西方向風場數據的總體趨勢上是有效的。

經過POD方法降階處理后可以得到相對應的投影系數。以第一階POD模態的投影系數Ck1為例，通過Fourier級數擬合，可以得出其Fourier預測模型，分別從東西方向(見圖3)和南北方向(見圖4)進行分析。從圖3(a)和圖4(a)可看出，兩者的投影系數Ck1都具有周期性和規律性，都能夠通過Fourier級數將其擬合(分別見圖3(b)和圖4(b))。但相對于南北方向風場，東西方向風場的投影系數更緊密，周期性更強，Fourier級數的擬合程度更高，說明東西方向風場更適用于建立Fourier預測模型。

圖2 采用POD方法重建風場Fig.2 Reconstruction of wind field using POD method

利用東西方向風場的Fourier預測模型，可以獲得某一天特定時刻下的投影系數，能對該時刻下平流層20 km附近的東西方向風速隨高度的變化情況進行預測。如圖5所示，以2008年11月7日為例，5階POD模型能夠捕獲到實際風場的主要特征，Fourier預測模型能夠體現出5階POD模型的整體變化趨勢狀況。

圖3 系數擬合(東西方向)Fig.3 Coefficient fitting (east-west direction)

圖4 系數擬合(南北方向)Fig.4 Coefficient fitting (north-south direction)

圖5 各種方法下風矢量圖的比較Fig.5 Comparison of wind vector plots for various methods

將Fourier預測模型向外延伸和擴展，可以對未來特定時刻下的風場進行預測。通過2005—2009年風場數據的統計和分析，對2010年10月30日0:00平流層20 km附近的東西方向風速進行預測，結果如圖6所示，Fourier預測模型能預測到風速隨高度的變化情況，并且能夠體現出平流層低速風帶所在的位置高度。

圖6 實際風矢量圖與預測風矢量圖的對比Fig.6 Comparison of actual wind vector with predicted wind vector plots

圖8 不同高度的風速變化情況Fig.8 Change of wind speed at different altitudes

圖9 不同高度的殘差分析Fig.9 Residual analysis at different altitudes

3.2 Fourier預測精度分析

由圖5、圖6可知，針對不同的高度，Fourier預測模型與實際風場模型之間存在不同的偏差，說明Fourier預測精度與空間分布有關，下面通過模型之間的比較，對Fourier預測精度進行分析。

將圖5中3種風場模型分別進行比較，如圖7(a)所示，其中POD降階殘差表示5階POD風場模型與實際風場之差，Fourier擬合殘差表示Fourier級數進行擬合時與5階POD模型之差，Fourier預測殘差表示Fourier預測模型與實際風場之差。由圖7(a)可知，3種殘差的絕對值都小于5，其中POD降階殘差曲線的差波動性最小，說明5階POD模型與實際風場模型數據很接近。由于Fourier預測模型是基于POD方法降階處理后建立的，所以Fourier預測殘差等于POD降階殘差與Fourier擬合殘差兩者之和，說明Fourier預測精度不僅與POD降階模型所捕獲到的能量有關，而且和Fourier級數的擬合程度有關。將圖6中Fourier預測模型與實際風場模型進行比較，如圖7(b)所示，Fourier預測殘差曲線的波動范圍較小，在海拔高度19～22 km之間的殘差接近于0，說明在此范圍內，Fourier預測精度較高，而在15 km附近波動最大，說明預測精度較低。

由圖7可知，在20 km處Fourier預測精度相對較高，在15和24 km處的預測精度較低。分別對15、20和24 km，3個高度的東西方向風速隨時間的變化情況進行分析，并分別計算出實際風場模型、5階POD模型與Fourier預測模型兩兩之間的殘差，結果如圖8、圖9所示。由圖9可知，在這3個高度位置上，POD降階殘差都接近于零，說明對于5階POD模型，在任意高度上都能夠捕獲到實際風場的主要特征；對于Fourier預測模型，如圖8所示，在20 km處東西方向的風場數據相對于其他2個高度的數據更加緊湊、平穩，Fourier級數的擬合程度也更高，說明Fourier預測模型與實際風場隨時間變化的緊湊性與規律性有著密切聯系，風場數據越緊湊，周期性越明顯，Fourier預測精度就越高。

圖7 殘差分析Fig.7 Residual analysis

4 結 論

本文通過POD方法將平流層風場數據進行高精度高效率降階處理后，為Fourier預測模型對風場進行預測節省了大量的計算時間，大大提高了計算效率，得到：

1) 提出了一種對平流層風場數據進行降階處理的POD方法，在POD方法的基礎上，提出了一種可以對平流層風場進行預測的Fourier預測模型。

2) 以長沙地區為例，選取海拔高度10～30 km的5年風場數據，采用提出的POD方法與Fourier預測模型對風場進行了建模與預測，并對Fourier預測精度進行了分析。

3) 研究結果表明：采用POD方法可以對東西方向風場進行高效高精度的降階建模，由于南北方向風場變化極其不規則，不可以采用低階的POD模型進行建模；通過Fourier預測模型能對東西方向風場進行準確預測，預測精度與實際風場隨時間變化的規律性有關，風場數據越緊湊，周期性越明顯，Fourier預測精度越高。