基于插值與降維方法的輸電塔線(xiàn)體系隨機(jī)脈動(dòng)風(fēng)場(chǎng)有效模擬

楊雄駿，黃金山，張建國(guó)，雷 鷹

(廈門(mén)大學(xué) 建筑與土木工程學(xué)院，福建 廈門(mén) 361005)

輸電導(dǎo)線(xiàn)在風(fēng)致作用下將產(chǎn)生的低頻率、高振幅的振動(dòng)，因此，風(fēng)荷載是輸電塔線(xiàn)體系結(jié)構(gòu)設(shè)計(jì)必須考慮的重要因素，而隨機(jī)脈動(dòng)風(fēng)場(chǎng)的模擬是輸電塔線(xiàn)體系抗風(fēng)設(shè)計(jì)的首要任務(wù)[1]。關(guān)于風(fēng)場(chǎng)模擬的維度，白海峰等[2]考慮輸電塔線(xiàn)體系結(jié)構(gòu)特征與風(fēng)場(chǎng)特性，將三維風(fēng)場(chǎng)簡(jiǎn)化為多個(gè)一維風(fēng)場(chǎng)。沈國(guó)輝等[3]對(duì)比不同方法的計(jì)算一維風(fēng)場(chǎng)的效率與結(jié)果，對(duì)不同計(jì)算方法的適用性進(jìn)行了討論。樓文娟等[4]采用CFD 數(shù)值模擬與風(fēng)洞試驗(yàn)相結(jié)合研究山地風(fēng)場(chǎng)特性對(duì)輸電線(xiàn)路風(fēng)偏響應(yīng)的影響。韓楓等[5]也對(duì)輸電塔線(xiàn)體系進(jìn)行了三維多變量脈動(dòng)風(fēng)場(chǎng)模擬。

諧波疊加法是模擬脈動(dòng)風(fēng)速時(shí)程的傳統(tǒng)方法，該法雖有理論簡(jiǎn)單清晰等優(yōu)點(diǎn)，但計(jì)算十分耗時(shí)，尤其是當(dāng)需要模擬的點(diǎn)數(shù)較多時(shí)，效率變得更加低下[6-7]。在諧波疊加法中，Cholesky 分解與三角級(jí)數(shù)疊加將耗費(fèi)大量的計(jì)算時(shí)間，因此減小二者的計(jì)算量十分重要。Yang[8]在模擬過(guò)程中引入快速傅里葉變換(fast Fourier transform,FFT)技術(shù)，提高了三角級(jí)數(shù)疊加的計(jì)算效率。Ding等[9]在互譜密度矩陣的 Cholesky 分解中引入三次 Lagrange 插值，減少了 Cholesky 分解的次數(shù)。羅俊杰等[10]采用三次均勻B樣條插值來(lái)擬合分解譜密度曲線(xiàn)。Tao等[11]引入Hermite插值建立的分解譜函數(shù)矩陣在全頻域內(nèi)的近似函數(shù)表達(dá)。

另一方面，諧波疊加法在本質(zhì)上屬于Monte Carlo模擬方法，而在Monte Carlo模擬方法中，處理高維的隨機(jī)變量一直是個(gè)難題。需要進(jìn)行大量的隨機(jī)抽樣才能保證模擬的精度，但是數(shù)量巨大的樣本不僅增加了隨機(jī)模擬的計(jì)算量，更重要的是，結(jié)構(gòu)動(dòng)力響應(yīng)的計(jì)算量也隨之大大增加。Chen等[12]提出隨機(jī)諧和函數(shù)表示方法，并對(duì)頻率選點(diǎn)進(jìn)行了優(yōu)化。劉章軍等[13]通過(guò)引入隨機(jī)函數(shù)的約束條件，對(duì)脈動(dòng)風(fēng)場(chǎng)模擬的譜表示法降維處理。

目前對(duì)輸電塔-線(xiàn)系統(tǒng)風(fēng)場(chǎng)模擬的研究中，較少考慮到橫風(fēng)向與豎風(fēng)向的脈動(dòng)風(fēng)的相關(guān)性，或未針對(duì)輸電塔-線(xiàn)系統(tǒng)結(jié)構(gòu)特點(diǎn)對(duì)風(fēng)場(chǎng)進(jìn)行簡(jiǎn)化，導(dǎo)致計(jì)算效率較低。此外，雖然插值方法在改進(jìn)Cholesky 分解方面得到了廣泛的應(yīng)用，但是未有將插值方法與脈動(dòng)風(fēng)場(chǎng)模擬降維方法相結(jié)合的研究。因此本文對(duì)于輸電塔，主要考慮其順風(fēng)向與橫風(fēng)向的脈動(dòng)風(fēng)，而且認(rèn)為其兩個(gè)方向不存在相關(guān)性，將其簡(jiǎn)化為一維多變量進(jìn)行；而輸電導(dǎo)線(xiàn)需考慮橫風(fēng)向與豎風(fēng)向的脈動(dòng)風(fēng)，并認(rèn)為兩者存在相關(guān)性?；谝陨涎芯炕A(chǔ)，本文對(duì)輸電塔進(jìn)行一維多變量脈動(dòng)風(fēng)場(chǎng)模擬，但對(duì)輸電線(xiàn)進(jìn)行考慮二維多變量脈動(dòng)風(fēng)場(chǎng)的模擬，而且將改進(jìn)Cholesky分解的Hermite插值和諧波表示法的降維處理技術(shù)相結(jié)合進(jìn)行輸電塔線(xiàn)體系的脈動(dòng)風(fēng)場(chǎng)模擬，提出了一種較現(xiàn)有的輸電塔線(xiàn)體系脈動(dòng)風(fēng)場(chǎng)數(shù)值模擬方法更有效的模擬方法。

1 輸電塔線(xiàn)風(fēng)場(chǎng)體系

在笛卡兒坐標(biāo)系中空間P點(diǎn)t時(shí)刻的瞬時(shí)風(fēng)速U(P,t)由平均風(fēng)速和脈動(dòng)風(fēng)速組成,即

(1)

式中:u、v和w分別表示縱向、橫向和豎向風(fēng)脈動(dòng)分量；U(z)為縱向平均風(fēng)速，用指數(shù)律分布表示為

(2)

式中:zs為離地參考高度；α為風(fēng)速廓線(xiàn)指數(shù)。

如圖1所示的輸電塔-線(xiàn)系統(tǒng)中，輸電塔是垂直矗立的高聳結(jié)構(gòu)，風(fēng)場(chǎng)的垂向(Z方向)脈動(dòng)能量較低，對(duì)輸電塔可忽略不計(jì)，主要分析其順風(fēng)響應(yīng)(X方向)與橫風(fēng)響應(yīng)(Y方向)；而輸電導(dǎo)線(xiàn)由于受到沿線(xiàn)路方向(Y方向)的約束作用，而只計(jì)縱向(X方向)及垂向(Z方向)脈動(dòng)作用，且同時(shí)考慮這兩個(gè)方向的相關(guān)性，進(jìn)行二維多變量脈動(dòng)風(fēng)場(chǎng)模擬。

圖1 輸電塔線(xiàn)體系示意圖Fig.1 Schematic of transmission tower-line system

由于輸電塔線(xiàn)體系具有高度高的結(jié)構(gòu)特征，采用沿高度變化的風(fēng)場(chǎng)功率譜密度函數(shù)比較符合工程實(shí)際。因此對(duì)縱向風(fēng)速譜可采用的Kaimal譜，而對(duì)橫風(fēng)向和豎風(fēng)向均可采用Panofsky- Kaimal譜，其表達(dá)式分別如下[14]

(3)

(4)

(5)

式中：u*=KU(z)/ln(z/q0)為摩擦風(fēng)速；q0為地面粗糙長(zhǎng)度；K為Karman常數(shù)。

由實(shí)測(cè)數(shù)據(jù)可知，橫向脈動(dòng)風(fēng)速v(t)與縱向脈動(dòng)風(fēng)速u(mài)(t)及豎向脈動(dòng)風(fēng)速w(t)之間不存在相關(guān)性，但縱向脈動(dòng)風(fēng)速u(mài)(t)與豎向脈動(dòng)風(fēng)速w(t)之間具有相關(guān)性，因此引入縱向脈動(dòng)風(fēng)速u(mài)(t)與豎向脈動(dòng)風(fēng)速w(t)之間的交叉譜密度

(6)

綜上所述，對(duì)于輸電塔，可將風(fēng)場(chǎng)簡(jiǎn)化為兩個(gè)相互獨(dú)立的一維風(fēng)速場(chǎng)(X方向與Y方向)；對(duì)于輸電線(xiàn)導(dǎo)線(xiàn)，則需要模擬一個(gè)X、Z方向相關(guān)的二維風(fēng)速場(chǎng)。

2 輸電塔-線(xiàn)風(fēng)場(chǎng)模擬方法

2.1 基于傳統(tǒng)諧波疊加法及Cholesky分解

2.1.1 輸電導(dǎo)線(xiàn)的風(fēng)場(chǎng)模擬

一組有n個(gè)點(diǎn)的二維多變量的平穩(wěn)高斯隨機(jī)風(fēng)場(chǎng)f(t)=[f1(t)f2(t) …fn(t)]T，節(jié)點(diǎn)的風(fēng)速向量空間fj(t)=[uj(t)wj(t)]T。n個(gè)點(diǎn)的互功率譜密度矩陣S(ω)可以表示為

(7)

式中，考慮X、Z方向脈動(dòng)風(fēng)速相關(guān)，每個(gè)功率譜矩陣Sij(ω)為 2×2的子矩陣，Sij(ω)可以表示為

(8)

相干函數(shù)用來(lái)描述脈動(dòng)風(fēng)在空間上的相關(guān)性。采用Davenport相干函數(shù)表達(dá)式[15]

(9)

式中:r=x,y,z；s=u,v,w；zi、zj分別為離地高度；Crs為指數(shù)衰減系數(shù)。

S(ω)為半正定矩陣，進(jìn)行Cholesky分解可得

S(ω)=H(ω)·H*T(ω)

(10)

式中，H*T(ω)為H(ω)的轉(zhuǎn)置共軛矩陣。

S(ω)為實(shí)對(duì)稱(chēng)矩陣，故H(ω)為實(shí)數(shù)矩陣。對(duì)某一頻率點(diǎn)ω的譜密度矩陣進(jìn)行分解，可得

(11)

隨機(jī)過(guò)程fj(t)的樣本可由下式來(lái)模擬[16]

φ2m,l))

(12)

(13)

式中:N為頻率等分?jǐn)?shù)，取為足夠大的正整數(shù)，一般取2的整數(shù)次冪;Δω=ωu/N為頻率增量；ωu為上界截止頻率，即當(dāng)ω>ωu時(shí)S(ω)=0；φ為均勻分布于[0,2π]區(qū)間的隨機(jī)相位角。ωml為雙索引頻率，其表達(dá)式為

m=1,2，…，n；l=1,2，…，N

(14)

2.1.2 輸電塔的風(fēng)場(chǎng)模擬

輸電塔的風(fēng)場(chǎng)模擬的方法與上節(jié)輸電導(dǎo)線(xiàn)模擬風(fēng)場(chǎng)方法類(lèi)似，但縱向風(fēng)場(chǎng)與橫向風(fēng)場(chǎng)是分別獨(dú)立模擬。分解的譜密度矩陣S(ω)降為n階(輸電塔線(xiàn)風(fēng)速譜密度矩陣為2n階)。以豎向風(fēng)vj(t)(j=1,2,…,n)為例，對(duì)于一維n變量零均值的平穩(wěn)隨機(jī)過(guò)程vj(t)，同樣將其進(jìn)行Cholesky分解

S(ω)=H(ω)·H*T(ω)

(15)

而隨機(jī)過(guò)程vj(t)的樣本可由下式來(lái)模擬

(16)

雖然運(yùn)用FFT技術(shù)可以加快上述諧波疊加法中三角級(jí)數(shù)的疊加，然而在隨機(jī)過(guò)程的模擬中，仍需進(jìn)行大量的Cholesky分解，尤其是對(duì)輸電導(dǎo)線(xiàn)的風(fēng)場(chǎng)模擬，需考慮縱向脈動(dòng)風(fēng)與豎向脈動(dòng)風(fēng)的相關(guān)性，所以譜密度矩陣階數(shù)更大，Cholesky分解更加困難，使得風(fēng)場(chǎng)模擬效率大大降低。由于Hermite插值函數(shù)具有較高的準(zhǔn)確性，因此本文采用Hermite插值來(lái)減少諧波疊加法中Cholesky分解次數(shù),有效進(jìn)行對(duì)輸電塔-導(dǎo)線(xiàn)系統(tǒng)風(fēng)場(chǎng)的數(shù)值模擬。

2.2 基于Hermite插值改進(jìn)Cholesky分解

對(duì)輸電導(dǎo)線(xiàn)或輸電塔的數(shù)值模擬，僅需在k個(gè)插值節(jié)點(diǎn)頻率及其鄰點(diǎn)處進(jìn)行對(duì)互譜矩陣的 Cholesky 分解，而對(duì)其他頻率點(diǎn)的分解譜函數(shù)則采用Hermite插值擬合。這樣將傳統(tǒng)諧波合成法中的Cholesky分解次數(shù)由n×N次縮減為n×2k次(2k<

采用三階Hermite函數(shù)近似表達(dá)H矩陣，需要計(jì)算插值點(diǎn)處的函數(shù)值與其一階導(dǎo)數(shù)值。插值點(diǎn)處的一階導(dǎo)數(shù)值H′(ωi)，可利用插值點(diǎn)處的函數(shù)值與鄰近點(diǎn)ωi+Δω的差商近似。建立區(qū)間[ωi,ωi+1]內(nèi)的Hermite插值函數(shù)

a3(ω-ωi)2(ω-ωi+1)

(17)

a3=

(18)

將上述Hermite插值函數(shù)引入諧波疊加法，將顯著減少互譜密度矩陣S的Cholesky分解次數(shù)，插值間距可采用前密后疏指數(shù)分布的形式。

2.3 基于降維的風(fēng)場(chǎng)譜表示法

要得到可以完備描述n×N維隨機(jī)變量空間概率信息的樣本集，需要從n×N維隨機(jī)變量空間中采樣足夠多的樣本點(diǎn)。顯然，隨機(jī)變量空間的維度決定了采樣點(diǎn)的數(shù)量。但是這里存在兩個(gè)問(wèn)題，一是從高維空間中采樣出具有完備性的樣本集相當(dāng)困難；二是采樣效率十分低下。本文采用一種降維方法，可以先在低維隨機(jī)空間中采樣一個(gè)完備樣本集，再用隨機(jī)函數(shù)映射到高維空間中去，從而達(dá)到在高維隨機(jī)空間中獲得完備樣本集的目的。

2.3.1 輸電塔線(xiàn)風(fēng)場(chǎng)模擬的諧波疊加法的降維表示

將式(12)、(13)及(16)中的三角函數(shù)展開(kāi)得

cos(ωmlt+φml)=

cos(ωmlt)cos(φml)-sin(ωmlt)sin(φml)

(19)

式(12)、(13)及(16)中的隨機(jī)變量φml并沒(méi)有獨(dú)立性的要求，它僅存在以下的約束關(guān)系[18-19]

E[cos(φml)]=0,E[sin(φml)]=0

(20a)

E[cos(φml)sin(φm′l′)]=0

(20b)

(20c)

l=1,2,…,nsel

(20d)

為此，可構(gòu)造一個(gè)僅由兩個(gè)獨(dú)立的基本隨機(jī)變量表示所有φml并滿(mǎn)足式(20)的隨機(jī)函數(shù)映射

φml(Θ)=m×Θ1+l×Θ2

m=1,2,…,n;l=1,2,…,N

(21)

式中，Θ={Θ1,Θ2}為(0,2π)上均勻分布的獨(dú)立實(shí)值隨機(jī)變量。

由于隨機(jī)函數(shù)映射是確定性映射，為此通過(guò)隨機(jī)函數(shù)得到n×N個(gè)變量后，還需要隨機(jī)打亂其排列。這樣，每計(jì)算一個(gè)風(fēng)場(chǎng)樣本均只需要首先從2維的隨機(jī)變量空間中采樣一個(gè)點(diǎn){Θ1,Θ2}，然后通過(guò)式(22)的變換得到n×N個(gè)隨機(jī)變量，然后打亂排列，再代入式(12)～式(13)或式(16)即可得到一個(gè)風(fēng)場(chǎng)樣本。可見(jiàn)，降維表示后，為得到完備的樣本集所需要的采樣空間維度將大大降低。

2.3.2 基本隨機(jī)變量的完備抽樣

為了獲得一個(gè)完備的樣本集，數(shù)論選點(diǎn)法[20]是一個(gè)很好的選擇。主要包括：

(1) 選一個(gè)合適的素?cái)?shù)nsel作為樣本個(gè)數(shù)

(2) 計(jì)算nsel的一個(gè)原根a

(3) 計(jì)算整型向量h

h=([1,a])modnsel

(22)

(4) 計(jì)算最優(yōu)晶格點(diǎn)集

l=1,2,…,nsel

(23)

根據(jù)數(shù)論選點(diǎn)法得到樣本集具有完備性，且每個(gè)樣本的賦得概率均為1/nsel。由于φml與基本隨機(jī)變量之間存在確定性映射，所以模擬的風(fēng)場(chǎng)樣本的完備性可以得到繼承，即風(fēng)場(chǎng)樣本也是完備的，且各樣本的賦得概率也為1/nsel。完備的風(fēng)場(chǎng)樣本集可用于風(fēng)場(chǎng)隨機(jī)特性和結(jié)構(gòu)可靠性方面的研究。

3 數(shù)值算例與分析

3.1 數(shù)值算例模型

為了驗(yàn)證所提方法的準(zhǔn)確性，進(jìn)行圖2所示模型的風(fēng)場(chǎng)模擬。輸電塔上模擬點(diǎn)位置取在豎向軸線(xiàn)上，每個(gè)塔模擬4個(gè)點(diǎn)，分別模擬X,Y兩個(gè)方向的風(fēng)速，不考慮輸電塔之間相干性的影響。輸電線(xiàn)模擬點(diǎn)等間距布置，每層每跨均布4點(diǎn)，模擬X,Z方向的風(fēng)速，考慮各點(diǎn)相干性影響。設(shè)計(jì)風(fēng)速U(z10)=30 m/s，風(fēng)速廓線(xiàn)指數(shù)α=0.16，Karman常數(shù)K=0.4，地面粗糙長(zhǎng)度q0=0.01,上限截止頻率ωup=2π，頻率等分?jǐn)?shù)N=1 024，M=2N，插值間距Δ=64Δω，時(shí)間間隔Δt=0.5 s。

3.2 模擬結(jié)果分析

圖2 輸電塔線(xiàn)體系示意圖 (m)Fig.2 Schematic of transmission tower-line system (m)

圖3 H矩陣元素插值精度對(duì)比Fig.3 Comparison of interpolation accuracy of H matrix elements

選取輸電塔第4點(diǎn)與輸電線(xiàn)第16點(diǎn)為例，進(jìn)行統(tǒng)計(jì)分析。脈動(dòng)風(fēng)速時(shí)程曲線(xiàn)如圖4所示，模擬風(fēng)場(chǎng)功率譜與目標(biāo)功率譜如圖5所示，可以看出脈動(dòng)風(fēng)能量主要集中于低頻，且縱向風(fēng)譜能量大于橫向風(fēng)譜能量，橫向風(fēng)譜能量大于豎向風(fēng)譜能量，模擬出的縱向脈動(dòng)風(fēng)速也最大，橫向脈動(dòng)風(fēng)速次之，豎向脈動(dòng)風(fēng)速最小，模擬的脈動(dòng)風(fēng)的功率譜密度函數(shù)與目標(biāo)值擬合較好。模擬風(fēng)場(chǎng)互相關(guān)函數(shù)與目標(biāo)相關(guān)函數(shù)如圖6所示。在時(shí)滯等于0時(shí)，相關(guān)性達(dá)到最大，且縱向風(fēng)速的相關(guān)程度要大于橫向風(fēng)速的相關(guān)程度。模擬的脈動(dòng)風(fēng)的互相關(guān)函數(shù)與目標(biāo)值變化趨勢(shì)基本一致，都隨著時(shí)滯增大相關(guān)性迅速減小。對(duì)功率譜密度與相關(guān)函數(shù)的誤差進(jìn)行定量分析，以功率譜密度為例，定義誤差如下

(a)

(a)

(a)

(24)

表1 樣本的功率譜相對(duì)誤差Tab.1 Comparisons of PSD relative errors

在降維處理方面，本文模擬樣本數(shù)選取為113，233，383。以樣本數(shù)113為例，其原根可取為33，由式(22)可知整型向量h=[1，33]，隨后基于式(23)可獲得113×2維的最優(yōu)晶格點(diǎn)集Θ(133)={Θ1,Θ2}(133)。將Θ(133)={Θ1,Θ2}(133)通過(guò)隨機(jī)函數(shù)式(22)變換得到4×1 024個(gè)變量后，并隨機(jī)打亂其排列得到所需完備的隨機(jī)變量φml(Θ)，最后將φml(Θ)代入譜表示式(16)即可得到模擬樣本。以第4點(diǎn)為例，表2為生成相同樣本數(shù)時(shí)傳統(tǒng)方法和采用插值與降維譜表示法的均值誤差及耗時(shí)的比較。計(jì)算機(jī)為：CPU是Intel i5-4 200 H@2.80 GHz、4 G內(nèi)存，采用MATLAB軟件編程。均值誤差定義如下

表2 不同數(shù)目樣本的均值相對(duì)誤差以及耗時(shí)比較Tab.2 Comparison of relative errors of mean and computation time of different numbers of samples

(25)

根據(jù)模擬結(jié)果的統(tǒng)計(jì)分析，其統(tǒng)計(jì)特征檢驗(yàn)符合假定條件。并且和傳統(tǒng)譜表示法相比較，在生成相同數(shù)目樣本時(shí)，插值與降維譜表示法不僅可以達(dá)到更好的精度，同時(shí)耗時(shí)也明顯減少，表明所采用方法的合理性與高效性。

4 結(jié) 論

本文基于輸電塔線(xiàn)體系的結(jié)構(gòu)特征和風(fēng)場(chǎng)特性，對(duì)輸電塔與輸電線(xiàn)分別進(jìn)行風(fēng)場(chǎng)模擬。將輸電塔風(fēng)場(chǎng)簡(jiǎn)化為兩個(gè)獨(dú)立的一維風(fēng)場(chǎng)，而對(duì)輸電線(xiàn)，引入縱向與豎向交叉譜的影響，模擬一個(gè)二維風(fēng)場(chǎng)。由于傳統(tǒng)諧波疊加法中Cholesky需要較大計(jì)算量，根據(jù)函數(shù)曲線(xiàn)隨頻率的變化規(guī)律，采用Hermite插值來(lái)擬合Cholesky分解后譜密度曲線(xiàn)，從而提高計(jì)算效率。

進(jìn)一步應(yīng)用隨機(jī)函數(shù)的約束條件，將傳統(tǒng)的譜表示法中的隨機(jī)變量個(gè)數(shù)大大減少，實(shí)現(xiàn)隨機(jī)脈動(dòng)風(fēng)場(chǎng)的譜表示降維模擬，與傳統(tǒng)譜表示法相比，達(dá)到相同精度所需生成的代表性時(shí)程數(shù)量減少。

這樣，與目前相關(guān)的輸電塔線(xiàn)體系脈動(dòng)風(fēng)場(chǎng)數(shù)值模擬方法相比，本文提出了更高效的模擬方法。通過(guò)數(shù)值算例分析，表明了本文采用方法能快速準(zhǔn)確地模擬出符合要求的多變量隨機(jī)脈動(dòng)風(fēng)場(chǎng)，可為輸電塔線(xiàn)體系動(dòng)力響應(yīng)分析提供基礎(chǔ)。