張冰怡 宋葭禾 何紹文 陳尚賓

（1.華中科技大學(xué) 武漢 430074）

（2.武漢光電國家研究中心 武漢 430074）

1 引言

面對化石燃料緊缺以及污染等諸多問題，發(fā)展可再生能源正成為當(dāng)下社會的重要課題之一。海洋占地球表面積的71%，其中蘊(yùn)藏著大量的可再生能源，波浪能就是其中重要的組成部分，因此波浪能發(fā)電裝置成為了海洋可再生能源研究的一個(gè)重要方向。振蕩浮子式波浪能發(fā)電裝置效率較高，布置靈活，是當(dāng)下一個(gè)研究熱點(diǎn)［1］。為了最大化發(fā)電效率，需盡可能提高能量輸出系統(tǒng)（Power Take-off，PTO）的輸出功率。國內(nèi)外也有許多相關(guān)研究，但大多數(shù)在建模時(shí)僅考慮了浮子在線性波作用下作垂蕩運(yùn)動［2］，本文在此基礎(chǔ)上還考慮了由于海水在水平方向波動所引起的縱搖運(yùn)動，使模型更貼近浮子在海水中的真實(shí)運(yùn)動情況。

波浪能轉(zhuǎn)換裝置形式多樣，但總體可分為三個(gè)環(huán)節(jié)。第一級轉(zhuǎn)換將波浪能轉(zhuǎn)換為實(shí)體所持有的能量，通常表現(xiàn)為隨波浪進(jìn)行起伏的運(yùn)動的機(jī)械能；第二、三級通過中間轉(zhuǎn)換裝置將機(jī)械能傳遞并轉(zhuǎn)換為電能。本文主要研究第一級過程，則可以將振動的浮體簡化為一個(gè)彈簧-質(zhì)量-阻尼系統(tǒng)［2］，如圖1 所示。在波浪的作用下，浮子帶動振子運(yùn)動，再通過兩者的相對運(yùn)動驅(qū)動阻尼器做功并作為能量輸出。本文通過建立浮子動力學(xué)模型，求解使PTO 最大時(shí)的阻尼器阻尼系數(shù)。在基礎(chǔ)模型中僅考慮浮子做垂蕩運(yùn)動；在改進(jìn)模型中考慮垂蕩運(yùn)動和縱搖運(yùn)動。利用三種智能算法求解最優(yōu)阻尼系數(shù)，并將其與現(xiàn)有成果進(jìn)行比較。

圖1 浮子式波浪能發(fā)電裝置示意圖

2 模型建立與求解

為方便分析與建模，本文對模型作如下假設(shè)：1）假定海水無粘、無旋；忽略水面高度變化；2）忽略中軸、底座、隔層與PTO質(zhì)量及各種摩擦；每種情況下，初始時(shí)刻浮子和振子都平衡處于靜水中；3）浮子、振子密度均勻，在繞中軸轉(zhuǎn)動時(shí)視作剛體；4）浮子的錐體部分始終不浮出海面；5）對浮子、振子組成的系統(tǒng)只考慮重力和海水的作用力；6）水對浮子的作用力（矩）都作用于浮子的質(zhì)心。

2.1 基礎(chǔ)模型

基礎(chǔ)模型只考慮系統(tǒng)做垂蕩運(yùn)動，可以將浮體簡化為如圖2 所示的有阻尼強(qiáng)迫振動系統(tǒng)，該系統(tǒng)的強(qiáng)迫力來自浮體所受的垂直波浪力［3～4］。在考慮海水運(yùn)動對浮體的作用時(shí)，通常采用微幅波理論［5～6］，該理論中當(dāng)液體無旋、不可壓縮時(shí)可以用速度勢函數(shù)來描述液體波動情況［7］，本文將其簡化為波浪激勵(lì)力、靜水恢復(fù)力、垂直興波阻力［8］。建立與水面垂直的面上的坐標(biāo)系，規(guī)定豎直向上為正方向，以零時(shí)刻靜止水面為坐標(biāo)原點(diǎn)，對浮子和振子分別做受力分析。

圖2 浮子振子相對位移示意圖及受力分析圖

圖3 退火法求解最優(yōu)解及粒子群法求解最優(yōu)解的算法流程圖

浮子所受的合力寫作：

考慮F附加慣性力是由于推動浮體搖蕩的力還會推動附體周圍的流體運(yùn)動，等效于為浮體加上了一個(gè)附加質(zhì)量，即垂蕩附加質(zhì)量Me［9～10］?？蓪Ω∽恿谐觯?/p>

再同理對振子做受力分析，如圖2 所示。阻尼力做功的平均功率，即PTO的平均輸出功率為

求解方程組時(shí)運(yùn)用Matlab（MathoWorks，R2017a）提供的ODE45 求解器求得該微分方程組的數(shù)值解；最終要求使平均輸出功率最大化，即利用Matlab 求解優(yōu)化問題。本模型求解屬于全局最優(yōu)問題，常用算法包括全局搜索法（GlobalSearch），粒子群算法（Particle Swarm Optimization，PSO）和模擬退火算法（Simulated Annealing，SA）。我們分別利用三種方法對函數(shù)進(jìn)行求解獲取使功率最大化時(shí)的阻尼系數(shù)，同時(shí)比較分析不同算法的結(jié)果［11］。利用Matlab可繪出浮子、振子的振幅變化及速度變化（圖4）。

圖4 浮子、振子作垂蕩運(yùn)動時(shí)參數(shù)變化圖像

規(guī)定零時(shí)刻靜止水面為坐標(biāo)原點(diǎn)。由圖像可知，帶入我們假定的數(shù)據(jù)進(jìn)行計(jì)算，浮子起初受到海水的作用力開始在豎直方向上振動，同時(shí)帶動振子振動；大約在50s 時(shí)二者組成的系統(tǒng)達(dá)到平衡，以相同頻率振動。在阻尼系數(shù)為常數(shù)時(shí)用三種算法計(jì)算最大輸出功率和對應(yīng)的阻尼系數(shù)，結(jié)果和計(jì)算所用時(shí)間如表1 所示。所用計(jì)算機(jī)處理器位Intel Core i7，內(nèi)存16GB。

表1 阻尼系數(shù)為常數(shù)時(shí)三種算法的計(jì)算結(jié)果

由表中數(shù)據(jù)可看出，三種算法計(jì)算得到的最大功率為226W 左右，對應(yīng)阻尼系數(shù)為37200N·s/m。而三種算法中粒子群算法所需時(shí)間最少。

2.2 優(yōu)化模型

本問題中浮子不僅做垂蕩運(yùn)動還做縱搖運(yùn)動，因此在轉(zhuǎn)軸上還安裝了旋轉(zhuǎn)阻尼器和扭轉(zhuǎn)彈簧。此時(shí)浮子與振子在二維平面上運(yùn)動。

浮子在做垂蕩運(yùn)動的基礎(chǔ)上還做縱搖運(yùn)動，所以在之前的受力基礎(chǔ)上還需考慮扭轉(zhuǎn)彈簧和旋轉(zhuǎn)阻尼器的作用以及海水的各作用力對應(yīng)的旋轉(zhuǎn)力矩。考慮浮子在豎直平面上的剖面，其二維質(zhì)心與O點(diǎn)的距離為

其中，Xi為各部分質(zhì)心與錐頂?shù)木嚯x；Mi為各部分的質(zhì)量。浮子受力在基礎(chǔ)模型類似，但需考慮到部分力與豎直方向存在夾角。由牛頓第二定律，對浮子可列出：

為分析縱搖運(yùn)動時(shí)的角位移與角速度，分析浮子力矩，由剛體轉(zhuǎn)動定律可列出：

其中I為浮子的總轉(zhuǎn)動慣量；correction為對浮子轉(zhuǎn)動慣量的扭矩修正系數(shù)。這是因?yàn)榕まD(zhuǎn)彈簧與旋轉(zhuǎn)阻尼器對浮子的作用力作用在轉(zhuǎn)軸上，而我們在計(jì)算力矩時(shí)需將其轉(zhuǎn)移到浮子的質(zhì)心。扭轉(zhuǎn)修正系數(shù)的表達(dá)式為

再同理分析振子。分別求解上述微分方程，獲得浮子、振子的運(yùn)動規(guī)律。

此時(shí)的輸出功率由直線阻尼器和旋轉(zhuǎn)阻尼器兩部分做功組成。對于直線阻尼器和旋轉(zhuǎn)阻尼器，分別有

最終的系統(tǒng)輸出功率為兩種阻尼器的功率之和。利用問題二中的求全局最優(yōu)解的算法求解P最大時(shí)ζ1與ζ3的值。

三種算法求解最優(yōu)阻尼系數(shù)的計(jì)算結(jié)果和用時(shí)見表2。全局搜索法和粒子群法計(jì)算的結(jié)果相近，求得的最大輸出功率為325 W 左右，對應(yīng)直線阻尼系數(shù)和旋轉(zhuǎn)阻尼系數(shù)分別為64200N·s/m 和32200N·s/m；三種算法種粒子群法用時(shí)最短。由圖5和圖6可以比較得出，在此問題中，粒子群算法相較于退火法和全局搜索法具有顯著優(yōu)勢，原因可能來自于粒子群算法的迭代次數(shù)比退火法少了一個(gè)數(shù)量級，因此會更快收斂。

表2 做垂蕩、縱搖運(yùn)動時(shí)三種算法求解最優(yōu)阻尼系數(shù)的計(jì)算結(jié)果

圖5 粒子群算法迭代次數(shù)

圖6 模擬退火算法

圖5 為粒子群算法運(yùn)行過程中，最優(yōu)結(jié)果隨迭代次數(shù)的變化趨勢；圖6 為模擬退火算法運(yùn)行時(shí)的參數(shù)情況，圖6（a）為最優(yōu)結(jié)果隨迭代次數(shù)的變化趨勢，圖6（b）為實(shí)時(shí)結(jié)果隨迭代次的變化趨勢，圖6（c）為最優(yōu)結(jié)果對應(yīng)的自變量取值，圖6（d）為當(dāng)前結(jié)果對應(yīng)的自變量取值。

3 模型分析

根據(jù)文獻(xiàn)［13］，估算波浪在單位時(shí)間內(nèi)傳遞到浮子所在地的能量為

我們假定波幅與浮子垂蕩位移約為同一數(shù)量級，根據(jù)第三問計(jì)算的浮子位移，取平衡時(shí)的波幅H 為0.5m，得到單位時(shí)間內(nèi)波浪傳遞到浮子所在地的能量E 為1255.63 J，由于海浪上浮子的運(yùn)動不止垂蕩和縱搖，所以在我們模型中最為貼近實(shí)際的是運(yùn)動更為復(fù)雜的改進(jìn)模型，故取P 為325.12 W。計(jì)算得到效率為25.9%。

在浮子的振動已經(jīng)穩(wěn)定時(shí)，此時(shí)浮子做穩(wěn)定的有阻尼受迫振動，任意時(shí)刻的能量不等于耗散的能量，系統(tǒng)的能量處于隨時(shí)變化的狀態(tài)［14］，波浪傳遞來的能量分別被興波阻力、附加轉(zhuǎn)動慣量、附加質(zhì)量等系統(tǒng)外阻力耗散掉一部分，因此此效率比較合理。

4 結(jié)語

本文對傳統(tǒng)的浮子波浪能模型進(jìn)行修正改進(jìn)，考慮了浮子作縱搖運(yùn)動并建立動力學(xué)模型。利用全局搜索算法、粒子群算法與退火算法計(jì)算最優(yōu)阻尼系數(shù)，通過對三種算法的運(yùn)行時(shí)間、迭代次數(shù)、收斂結(jié)果的對比可以直觀地觀察到三種算法的表現(xiàn)，同時(shí)還可以互相印證，得出更為可信的結(jié)果。其中粒子群算法依靠在定義域內(nèi)隨機(jī)產(chǎn)生大量粒子完成搜索，搜索速度最快，且結(jié)果相對準(zhǔn)確，是研究本問題的最優(yōu)算法。但是粒子群算法理論上容易陷入局部最優(yōu)解，所以還需在后續(xù)研究中不斷優(yōu)化、驗(yàn)證。

粒子群算法計(jì)算改進(jìn)后模型得到發(fā)電效率為25.9%。由文獻(xiàn)［2］的結(jié)果可知，10%～30%是目前浮子式波浪能裝置發(fā)電的合理范圍，而我們的模型效率落在此范圍之間，從側(cè)面印證我們建立了一個(gè)較合理的模型。