曹錦超 項建強(qiáng) 金 磊

（1.三峽新能源山東昌邑發(fā)電有限公司，山東 濰坊 261300；2.中國三峽新能源（集團(tuán)）股份有限公司山東分公司，山東 濟(jì)南 261300）

0 引言

在結(jié)構(gòu)振動分析中，準(zhǔn)確的有限元模型是其關(guān)鍵因素，模型誤差、測量誤差和加工誤差等因素的影響會導(dǎo)致數(shù)值模型與實測結(jié)構(gòu)存在較大差異，無法用其準(zhǔn)確進(jìn)行結(jié)構(gòu)動力響應(yīng)預(yù)測和安全狀態(tài)評估[1]。為了獲取能近似代替實測結(jié)構(gòu)的有限元模型，需要根據(jù)結(jié)構(gòu)的振動測試數(shù)據(jù)對有限元模型進(jìn)行修正。隨著優(yōu)化算法的發(fā)展，基于優(yōu)化算法的模型優(yōu)化算法得到國內(nèi)外學(xué)者的廣泛關(guān)注[2]。與單目標(biāo)優(yōu)化相比，多目標(biāo)問題通常更困難、更復(fù)雜[3]。在單個目標(biāo)中，必須找到最優(yōu)解，它通常是解空間中的1 個點，但是存在多個相等全局最優(yōu)點的情況除外。對多目標(biāo)優(yōu)化問題來說，存在多個最優(yōu)解，而在實際工程應(yīng)用中，不同的目標(biāo)往往相互沖突，有時可能根本不存在真正的最優(yōu)解，因此通常需要運(yùn)用一些近似方法獲取最優(yōu)解的近似值[4]。多目標(biāo)優(yōu)化有許多算法[5-8]，這些算法通常來源于自然界中不同生物的特征。宋鈺等人[9]提出了一種自適應(yīng)布谷鳥搜索算法，研究表明，該方法的全局尋優(yōu)優(yōu)于現(xiàn)有優(yōu)化算法，例如粒子群優(yōu)化算法。在該文中，使用擴(kuò)展布谷鳥算法來解決多目標(biāo)尋優(yōu)問題，并發(fā)展一個多目標(biāo)布谷鳥搜索（MOCS）算法，將應(yīng)用它來解決海上風(fēng)電結(jié)構(gòu)的模型修正問題。

1 多目標(biāo)布谷鳥算法

為了將布谷鳥算法從單目標(biāo)優(yōu)化擴(kuò)展到多目標(biāo)優(yōu)化問題，下面將詳細(xì)論述該文所提出方法的理論和步驟。

1.1 布谷鳥行為

為了繁衍后代，布谷鳥會將自己產(chǎn)的卵放在其他物種的巢穴內(nèi)進(jìn)行孵化，為了模擬布谷鳥算法，通常需要假定3 個理想狀態(tài)：1）布谷鳥每次僅產(chǎn)1 個卵，選擇的巢穴是隨機(jī)的。2）通過選擇，最適合孵化后代的巢穴被保留至下一代。3）可選擇的巢穴數(shù)量是固定的。

1.2 Levy 飛行搜索機(jī)制

Levy 飛行是一種隨機(jī)游動，其步長由Levy 分布得出，通常用簡單的冪指數(shù)公式表示，如公式（1）所示。

式中：L為Levy 分布函數(shù)；s為飛行步長；μ為最小步長，μ<0表示最下步長；γ為尺度參數(shù)。

當(dāng)s→∞時，Levy 飛行為，這通常表示Levy 分布的一個特例。

Levy 分布一般運(yùn)用傅里葉變換來定義，如公式（2）所示。

式中：F（s）為Levy 分布的傅里葉變換；β為levy 分布指數(shù)；α為尺度參數(shù)，其表示的逆的積分一般不容易獲取，除少數(shù)特殊情況外，一般沒有解析表達(dá)式。

當(dāng)β=2 時，可以得到公式（3）。

公式（3）對應(yīng)于高斯分布。

另外一個特例為當(dāng)β=1 時，可以得到公式（4）。

公式（4）對應(yīng)于柯西分布，Levy 分布如公式（5）所示。

式中：p為Levy 分布；x為飛行步長；μ為位置參數(shù)；γ為控制分布尺度的參數(shù)。

對一般情況來說，其逆運(yùn)算的積分如公式（6）所示。

式中：t為時間變量。

當(dāng)s的值足夠大時，可以獲取公式（6）的準(zhǔn)確估計值，如公式（7）所示。

式中：s和k為不同飛行步長。

γ函數(shù)Γ（z）可以為公式（8）。

式中：t為時間變量；z為自變量。

當(dāng)z=n（n為正數(shù)）時，Γ（n）=（n-1）!。

在解決多維尋優(yōu)問題時，Levy 飛行明顯優(yōu)于布朗隨機(jī)漫步，效率高的原因如公式（9）所示。

式中：δ為方差；t為時間變量。

1.3 多目標(biāo)布谷鳥尋優(yōu)算法

對具有K個不同目標(biāo)的多目標(biāo)優(yōu)化問題來說，修改布谷鳥第一個和第三個假設(shè)條件，以滿足多目標(biāo)尋優(yōu)需求：1）每只布谷鳥1 次產(chǎn)K個蛋，并將其傾倒在隨機(jī)選擇的巢中。第K個蛋對應(yīng)第K個目標(biāo)的解決方案。2）每個巢穴都有可能被棄置，并根據(jù)蛋的相似性或差異性產(chǎn)生新的巢穴，同時運(yùn)用新的隨機(jī)組合產(chǎn)生分集。當(dāng)產(chǎn)生新的解時，對第i 個未知數(shù)來說，位置更新如公式（10）所示。

式中：x為修正系數(shù)；t為算法的第t次循環(huán)；α為求解問題的步長，取值大于0。

一般情況下取a=0（1）。為了使α能夠解決各種問題，須滿足公式（11）。

式中：α0為常數(shù)；xi和xj為前后2 次循環(huán)產(chǎn)生的解。

公式（11）表示當(dāng)前后2 次產(chǎn)生的解比較接近時，獲取全局最優(yōu)解的難度較大，因此需要運(yùn)用前后2 次迭代運(yùn)算產(chǎn)生的解之間的差異，以產(chǎn)生新的解。⊕表示點乘，Levy 飛行提供了一個隨機(jī)行走，而隨機(jī)步長是從Levy 分布中提取的較大步長。隨機(jī)分布如公式（12）所示。

式中：Levy～u為隨機(jī)分布；t為時間變量；β為指數(shù)參數(shù)。

當(dāng)公式（12）所示的方程的方差無窮大時，其平均值也無窮大。布谷鳥的連續(xù)跳躍步長基本上形成了一個隨機(jī)行走過程，該過程遵循具有重尾的冪指數(shù)規(guī)律步長分布。對Lévy 飛行來說，產(chǎn)生步長為s的樣本難度較大，如公式（13）～公式（15）所示。

式中：N為正態(tài)分布；Γ為γ分布。

1.4 非劣最優(yōu)目標(biāo)域

對優(yōu)化問題的解向量來說，u=（u1，…，un）T∈R，當(dāng)每次循環(huán)且ui≤vi時，解向量可以控制其他向量v=（v1，…，vn），如公式（16）所示。

對目標(biāo)函數(shù)最大值問題來說，可以通過定義?和?代替<和>。因此π，點x*∈R 稱為非支配解，該值不受任何條件支配。多目標(biāo)優(yōu)化問題的非劣最優(yōu)目標(biāo)域PF 可以定義為非支配解，如公式（17）所示。

式中：PF為非支配解；R 為實數(shù)域；s'為飛行步長。

根據(jù)搜索空間非劣最優(yōu)目標(biāo)域如公式（18）所示。

式中：f（x）為自變量x的函數(shù)；f（x'）為自變量x'的函數(shù)。

其中，f=（f1，…，fk）T。

2 基于改進(jìn)布谷鳥算法的模型修正方法

2.1 適應(yīng)度函數(shù)構(gòu)建

對基于優(yōu)化算法的模型修正方法來說，在計算結(jié)構(gòu)修正系數(shù)時，通常運(yùn)用實測結(jié)構(gòu)和有限元模型的振動特征參數(shù)構(gòu)建目標(biāo)函數(shù)，并以此作為優(yōu)化算法的適應(yīng)度函數(shù)。為了準(zhǔn)確獲取結(jié)構(gòu)的修正系數(shù)，該方法運(yùn)用實測結(jié)構(gòu)和數(shù)值模型的固有頻率和振型構(gòu)建目標(biāo)函數(shù)，如公式（19）所示。

式中：ωu，i和ωe.i分別為有限元模型和實測結(jié)構(gòu)第i階固有頻率；?iju和?ije分別為有限元模型和實測結(jié)構(gòu)第i階振型第j個自由度振型幅值。

2.2 修正系數(shù)計算

對基于優(yōu)化算法的模型修正方法來說，可以通過計算適應(yīng)度函數(shù)值最小值獲得結(jié)構(gòu)修正系數(shù)，即在每次循環(huán)時需要最小化適應(yīng)度函數(shù)，如公式（20）所示。

2.3 算法計算流程

為了便于理解該文所提出的模型修正方法，該方法的計算步驟如下所示：1）初始化改進(jìn)布谷鳥算法的參數(shù)，即種群大小Popsize 和循環(huán)次數(shù)Maxiter，根據(jù)特征值分析方運(yùn)用布谷鳥算法產(chǎn)生的系數(shù)計算有限元模型的固有頻率ωu和振型?u。2）運(yùn)用模態(tài)參數(shù)識別方法計算實測結(jié)構(gòu)固有頻率ωe和振型?e，并根據(jù)公式（20）構(gòu)建目標(biāo)函數(shù)。3）根據(jù)公式（10）～公式（18），不斷更新算法中布谷鳥的位置和速度，并在每次循環(huán)時通過最小化公式（19）所示的適應(yīng)度函數(shù)計算結(jié)構(gòu)的修正系數(shù)，算法達(dá)到終止條件后，通過輸出算法全局最優(yōu)值就可以獲取結(jié)構(gòu)的修正系數(shù)。

3 單樁式海上風(fēng)電結(jié)構(gòu)數(shù)值算例

為了驗證該文所提方法的有效性，選用單樁式風(fēng)電結(jié)構(gòu)數(shù)值模型開展相關(guān)研究，為了模擬海上風(fēng)電結(jié)構(gòu)的真實服役環(huán)境，以山東濰坊某風(fēng)電場單樁式風(fēng)電結(jié)構(gòu)為例設(shè)計數(shù)值模型。

3.1 風(fēng)電結(jié)構(gòu)數(shù)值模型

風(fēng)電結(jié)構(gòu)數(shù)值模型縮尺比為1 ∶70，模型由基礎(chǔ)、塔筒和機(jī)艙3 個部分組成，模型結(jié)構(gòu)尺寸和示意圖如圖1 所示?；A(chǔ)高度為1 200 mm，由等直徑鋼管組成，截面面積為?10 mm×2 mm，塔筒高度為800 mm，由變直徑鋼管組成，直徑變化范圍為500 mm～900 mm，壁厚為1 mm，機(jī)艙運(yùn)用集中質(zhì)量模擬。材料選用Q235 鋼材，彈性模量為210 GPa，密度為7 850 kg/m3，泊松比為0.3。為了便于后續(xù)的模型修正，將數(shù)值模型劃分為10 個單元，每個單元具有2 個節(jié)點和4 個自由度。

圖1 風(fēng)電結(jié)構(gòu)數(shù)值模型

3.2 荷載與工況

為了模擬實際海上風(fēng)電結(jié)構(gòu)振動特性的改變，假設(shè)僅單元1、單元6 和單元7 的剛度和質(zhì)量發(fā)生變化，單元剛度變化系數(shù)分別為α1=-0.25、α6=-0.3 和α7=-0.35，單元質(zhì)量變化系數(shù)分別為β1=0.3、β6=0.3 和β7=0.25。根據(jù)是否考慮噪聲和空間不完備條件的影響，共設(shè)計3 個工況，具體參數(shù)見表1。為了獲取結(jié)構(gòu)的振動響應(yīng)信號，在數(shù)值模型頂部節(jié)點施加周期荷載，荷載如圖2 所示，運(yùn)用Newmark-β 方法計算結(jié)構(gòu)的振動響應(yīng)。

圖2 周期荷載時程曲線

表1 海上風(fēng)電結(jié)構(gòu)模擬工況

3.3 模型修正結(jié)果

先運(yùn)用結(jié)構(gòu)固有頻率和振型構(gòu)建目標(biāo)函數(shù)，再運(yùn)用改進(jìn)的布谷鳥算法計算不同工況下結(jié)構(gòu)的修正系數(shù)，結(jié)果如圖3 所示。為分析空間不完備和測量噪聲對算法準(zhǔn)確性的影響，運(yùn)用特征值分析方法計算不同工況下修正后模型的固有頻率，結(jié)果見表2。由表2 可知，在5%噪聲和空間不完備條件的影響下，修正后模型前四階固有頻率誤差最大值為-2.244%，說明該文提出的模型修正方法的準(zhǔn)確性和魯棒性較高，也驗證了算法的有效性。

表2 不同工況下修正后模型固有頻率

圖3 風(fēng)電結(jié)構(gòu)修正系數(shù)計算結(jié)果

4 結(jié)語

該文基于改進(jìn)的布谷鳥算法提出了一種新的模型修正方法，并運(yùn)用風(fēng)電結(jié)構(gòu)數(shù)值模型驗證了方法的有效性。結(jié)果表明，在5%噪聲和空間不完備條件的影響下，該文提出的模型修正方法能夠準(zhǔn)確計算結(jié)構(gòu)的修正系數(shù)。