基于AR模型的空間脈動(dòng)風(fēng)速時(shí)程模擬方法研究

趙海霞　（南京航空航天大學(xué)金城學(xué)院，江蘇 南京 211156）

?

基于AR模型的空間脈動(dòng)風(fēng)速時(shí)程模擬方法研究

趙海霞（南京航空航天大學(xué)金城學(xué)院，江蘇 南京 211156）

空間脈動(dòng)風(fēng)速時(shí)程的模擬是橋梁抖振時(shí)域計(jì)算的前提。利用基于自回歸（AR）模型的線性回歸濾波器法，編制了模擬脈動(dòng)風(fēng)速的Matlab計(jì)算程序，模擬了一座大跨連續(xù)剛構(gòu)橋最大懸臂施工階段多點(diǎn)的水平和豎向脈動(dòng)風(fēng)速，并進(jìn)行了功率譜和相關(guān)性檢驗(yàn)，發(fā)現(xiàn)模擬效果較好。

脈動(dòng)風(fēng)速；AR模型；線性回歸濾波器法

1　研究背景

隨著我國橋梁工程的建設(shè)不斷向大跨度方向發(fā)展，橋梁的長細(xì)化使其剛度和阻尼不斷下降，導(dǎo)致結(jié)構(gòu)對風(fēng)的敏感性不斷增加，對橋梁結(jié)構(gòu)進(jìn)行風(fēng)效應(yīng)分析也逐漸成為橋梁結(jié)構(gòu)分析中的重要一環(huán)。由于時(shí)域分析方法可以直接計(jì)算頻域分析中難以處理的非線性響應(yīng)，因而對大跨度橋梁的風(fēng)載作用下的抖振行為進(jìn)行非線性時(shí)域分析日益受到人們的關(guān)注。在抖振響應(yīng)時(shí)域分析中，首先要根據(jù)目標(biāo)功率譜函數(shù)人工模擬空間脈動(dòng)風(fēng)場，從而獲得離散的風(fēng)速時(shí)程和離散的抖振力。

在人工模擬空間脈動(dòng)風(fēng)場方面，基于自回歸（Auto-Regressive，AR）模型的線性回歸濾波器法具有速度快、計(jì)算量小的特點(diǎn)，且精度較高，因此得到了廣泛應(yīng)用[1][2]。該模型將均值為零的白噪聲隨機(jī)系列通過線性濾波器，使其輸出為具有指定譜特征的平穩(wěn)隨機(jī)過程[3]。本文首先詳細(xì)介紹了該方法的計(jì)算原理和過程，利用Matlab軟件編制了相應(yīng)的計(jì)算程序；采用該方法模擬了一座大跨連續(xù)剛構(gòu)橋最大懸臂施工階段多點(diǎn)的水平和豎向脈動(dòng)風(fēng)速，并進(jìn)行了功率譜和相關(guān)性檢驗(yàn)。

2　自回歸（AR）模型

2.1求回歸系數(shù)

采用AR法推廣到模擬多維風(fēng)速過程的技術(shù)，M個(gè)相關(guān)的隨機(jī)風(fēng)過程[u（t）]=[u1（t），…，uM（t）]T可由下式生成：

式中：[u（t-kΔt）]=[u1（t-kΔt），…，uM（t-k Δt）]T；[N（t）]=[N1（t），…，NM（t）]T，Ni（t）為均值為0、具有給定協(xié)方差的正態(tài)分布隨機(jī)過程，i ＝1，…，M；[ψk]為M×M階矩陣，k=1，…，p；p 為AR模型的階數(shù)，一般取4或5。

對任一空間點(diǎn)i（i=1，…，M）具有時(shí)間差的隨機(jī)風(fēng)過程ui（t）與ui（t-kΔt）的協(xié)方差可表示為：

由于ui（t）與ui（t-kΔt）為均值0的平穩(wěn)隨機(jī)風(fēng)過程，其協(xié)方差的值僅為時(shí)間差的函數(shù)，式（2）可改寫為：

在式（1）同時(shí)右乘[u（t-jΔt）]=[u1（t-jΔt），…，uM（t-jΔt）]，j=1，…，p，并兩邊同時(shí)取數(shù)學(xué)期望（均值），考慮到[N（t）]的均值為0，且與隨機(jī)風(fēng)過程ui（t）獨(dú)立，以及協(xié)方差Ru（jΔt）為偶函數(shù)，可得到協(xié)方差Ru（jΔt）與回歸系數(shù)ψk之間的關(guān)系，寫成矩陣形式，有：

其中

根據(jù)隨機(jī)振動(dòng)理論[3]，功率譜密度與相關(guān)函數(shù)（協(xié)方差）之間符合維納-辛欽（Wiener-Khintchine）公式，即：

通過rjk（n）考慮風(fēng)速時(shí)程的空間相關(guān)特性，其三維表達(dá)式為：

式中：Cx、Cy和Cz分別表示空間任意兩點(diǎn)左右、上下和前后的衰減系數(shù)，一般分別取為16、8和10[4]；（xi，yi，zi）、（xk，yk，zk）分別為空間i、k點(diǎn)的三維坐標(biāo)，i、k=1，…，M；x、y、z分別為垂直于來流的水平方向、來流方向和豎向；分別表示第i點(diǎn)和第k點(diǎn)的平均風(fēng)速。

求解式（4）給出的線性方程組，可以得到回歸系數(shù)矩陣[ψ]。

2.2求給定方差的隨機(jī)過程[N（t）]

對式（1）同時(shí)右乘[u（t）]=[u1（t），…，uM（t）]，并取期望，有

求出[RN]后，對其作喬利斯基（Cholesky）分解[RN]=[L][L]T，則

式中：[L]為下三角矩陣；[n（t）]=[n1（t），…，nM（t）]T，ni（t）是均值為0、方差為1且彼此相互獨(dú)立的正態(tài)隨機(jī)過程，i=1，…，M。

2.3求最終的M個(gè)隨機(jī)過程

求出回歸系數(shù)矩陣[ψ]及[RN]后，可按式（1）求解出M個(gè)空間相關(guān)的隨機(jī)風(fēng)過程。

將式（1）按時(shí)間間隔Δt離散化，分別考慮三種情況，即ui（t）為①偶函數(shù)；②奇函數(shù)；③當(dāng)t<0時(shí)，ui（t）＝0。由此可得出不同的矩陣方程形式。其中，以第三種假設(shè)計(jì)算起來最為方便，其遞推的矩陣表達(dá)式為：

從而得到M個(gè)具有時(shí)間、空間相關(guān)、時(shí)間間隔Δt的離散脈動(dòng)風(fēng)速時(shí)程向量。

為了避免模擬結(jié)果失真，Δt必須滿足以下條件：

式中：ωup是截?cái)囝l率（rad/s）。

3　風(fēng)場模擬的Matlab程序

根據(jù)線性回歸濾波器法的自回歸（AR）模型的思路，所編制的Matlab程序流程如圖1所示。

圖1　自回歸模型生成風(fēng)速時(shí)程流程圖

4　算例

以某連續(xù)剛構(gòu)橋最大懸臂施工階段為對象，選定的模擬點(diǎn)即風(fēng)荷載加載節(jié)點(diǎn)共23點(diǎn)，其中主梁上17點(diǎn)，墩上6點(diǎn)。圖2所示為部分模擬點(diǎn)位置示意圖。主梁上的模擬點(diǎn)模擬水平脈動(dòng)風(fēng)速和豎向脈動(dòng)風(fēng)速，墩上的模擬點(diǎn)僅模擬水平脈動(dòng)風(fēng)速。

圖2　風(fēng)速模擬點(diǎn)示意圖（單位：m）

順風(fēng)向風(fēng)速模擬目標(biāo)譜采用Kaimal水平脈動(dòng)風(fēng)速譜[5]：

豎直向風(fēng)速模擬目標(biāo)譜采用Panofsky譜：

式中：Su（n）、Sw（n）分別為水平向和豎向脈動(dòng)風(fēng)速功率譜；n為脈動(dòng)風(fēng)頻率（Hz）；Z為有效高度；U（Z）為高度Z處的平均風(fēng)速；u*為氣流摩阻速度（亦稱剪切速度）；K是無量綱常數(shù)，K≈0.4；Zd為零平面高度；為周圍建筑物的平均高度；kd是地面阻力系數(shù)；Z0是地面粗糙長度。

考慮上述風(fēng)速模擬的目標(biāo)譜為非均方規(guī)一的單邊功率譜，通過與實(shí)際脈動(dòng)風(fēng)速的均方差比較，可以換算得到圓頻率表達(dá)的風(fēng)譜形式：

模擬的其他參數(shù)有：場地類別為Ⅱ類；場地基本風(fēng)速為24.1m/s[6]，其他高度處的風(fēng)速按照抗風(fēng)規(guī)范的指數(shù)分布，風(fēng)速剖面無量綱冪指數(shù)α=0.16；地面粗糙長度Z0取為0.05；零平面高度Zd取為0；截至頻率取為2Hz；頻率等分?jǐn)?shù)取為12000；采樣時(shí)距取為0.1s，共取12000步，1200s。分別選取左懸臂端（點(diǎn)1）、距左懸臂端8m處（點(diǎn)2）、主梁0號塊中點(diǎn)（點(diǎn)9）、右懸臂端處（點(diǎn)17）和墩左肢中部（點(diǎn)20）共五點(diǎn)作為展示對象。圖3是用自回歸（AR）模型模擬的這五點(diǎn)的水平風(fēng)速時(shí)程，圖4是其豎向風(fēng)速時(shí)程（不含點(diǎn)20）。

圖3　各點(diǎn)水平脈動(dòng)風(fēng)速

圖4　各點(diǎn)豎向脈動(dòng)風(fēng)速

圖5～圖6是分別用改進(jìn)的平均周期圖法和最大熵法[5]對所模擬的點(diǎn)1風(fēng)速時(shí)程做的譜估計(jì)；圖7所示為點(diǎn)1的自相關(guān)函數(shù)檢驗(yàn)，圖8所示點(diǎn)1和點(diǎn)2的互相關(guān)函數(shù)的檢驗(yàn)。從圖中我們可以看出，模擬值與目標(biāo)值之間相差很小，可見模擬的結(jié)果良好。

圖5　點(diǎn)1功率譜密度函數(shù)（最大熵法）

圖6　點(diǎn)1功率譜密度函數(shù)（多周期圖法）

圖7　點(diǎn)1自相關(guān)函數(shù)

圖8　點(diǎn)1、2互相關(guān)函數(shù)

4　結(jié)語

本文詳細(xì)介紹了基于自回歸（Auto-Regressive，AR）模型的線性回歸濾波器法的計(jì)算原理和過程，利用Matlab軟件編制了相應(yīng)的計(jì)算程序。采用該程序模擬了一座大跨連續(xù)剛構(gòu)橋最大懸臂施工階段多點(diǎn)的水平和豎向脈動(dòng)風(fēng)速，并進(jìn)行了功率譜和相關(guān)性檢驗(yàn)，發(fā)現(xiàn)模擬效果較好。因此本文的計(jì)算方法和程序具有一定的推廣價(jià)值。

[1]舒新玲，周岱.風(fēng)速時(shí)程AR模型及其快速實(shí)現(xiàn)[J].空間結(jié)構(gòu)，2003，9 （4）：27-32.

[2]李元齊，董石麟.大跨度空間結(jié)構(gòu)風(fēng)荷載模擬技術(shù)研究及程序編制[J].空間結(jié)構(gòu)，2001，7（3）：3-11.

[3]俞載道，曹國敖.隨機(jī)振動(dòng)理論及其應(yīng)用[M].上海:同濟(jì)大學(xué)出版社，1988.

[4]邊建烽，魏德敏.大跨空間結(jié)構(gòu)風(fēng)速時(shí)程的數(shù)值模擬理論[J].暨南大學(xué)學(xué)報(bào)（自然科學(xué)版），2005，26（1）：87-90.

[5]JTG T D60-01-2004，公路橋梁抗風(fēng)設(shè)計(jì)規(guī)范[S].

[6]張文明.大跨連續(xù)剛構(gòu)橋最大懸臂施工階段風(fēng)致響應(yīng)分析[D].武漢：華中科技大學(xué)，2007.

[7]黃文梅.信號分析與處理：MATLAB語言與應(yīng)用[M].國防電子科技大學(xué)出版社，2000.

TU442.5+9

A

1007-7359(2016)03-0058-04

10.16330/j.cnki.1007-7359.2016.03.020

趙海霞（1981-），女，畢業(yè)于西安建筑科技大學(xué)，碩士；講師，南京航空航天大學(xué)金城學(xué)院土木工程教研室。