基于智能算法优化的浮子式波浪发电装置阻尼设计*

张冰怡 宋葭禾 何绍文 陈尚宾

（1.华中科技大学 武汉 430074）

（2.武汉光电国家研究中心 武汉 430074）

1 引言

面对化石燃料紧缺以及污染等诸多问题，发展可再生能源正成为当下社会的重要课题之一。海洋占地球表面积的71%，其中蕴藏着大量的可再生能源，波浪能就是其中重要的组成部分，因此波浪能发电装置成为了海洋可再生能源研究的一个重要方向。振荡浮子式波浪能发电装置效率较高，布置灵活，是当下一个研究热点［1］。为了最大化发电效率，需尽可能提高能量输出系统（Power Take-off，PTO）的输出功率。国内外也有许多相关研究，但大多数在建模时仅考虑了浮子在线性波作用下作垂荡运动［2］，本文在此基础上还考虑了由于海水在水平方向波动所引起的纵摇运动，使模型更贴近浮子在海水中的真实运动情况。

波浪能转换装置形式多样，但总体可分为三个环节。第一级转换将波浪能转换为实体所持有的能量，通常表现为随波浪进行起伏的运动的机械能；第二、三级通过中间转换装置将机械能传递并转换为电能。本文主要研究第一级过程，则可以将振动的浮体简化为一个弹簧-质量-阻尼系统［2］，如图1 所示。在波浪的作用下，浮子带动振子运动，再通过两者的相对运动驱动阻尼器做功并作为能量输出。本文通过建立浮子动力学模型，求解使PTO 最大时的阻尼器阻尼系数。在基础模型中仅考虑浮子做垂荡运动；在改进模型中考虑垂荡运动和纵摇运动。利用三种智能算法求解最优阻尼系数，并将其与现有成果进行比较。

图1 浮子式波浪能发电装置示意图

2 模型建立与求解

为方便分析与建模，本文对模型作如下假设：1）假定海水无粘、无旋；忽略水面高度变化；2）忽略中轴、底座、隔层与PTO质量及各种摩擦；每种情况下，初始时刻浮子和振子都平衡处于静水中；3）浮子、振子密度均匀，在绕中轴转动时视作刚体；4）浮子的锥体部分始终不浮出海面；5）对浮子、振子组成的系统只考虑重力和海水的作用力；6）水对浮子的作用力（矩）都作用于浮子的质心。

2.1 基础模型

基础模型只考虑系统做垂荡运动，可以将浮体简化为如图2 所示的有阻尼强迫振动系统，该系统的强迫力来自浮体所受的垂直波浪力［3～4］。在考虑海水运动对浮体的作用时，通常采用微幅波理论［5～6］，该理论中当液体无旋、不可压缩时可以用速度势函数来描述液体波动情况［7］，本文将其简化为波浪激励力、静水恢复力、垂直兴波阻力［8］。建立与水面垂直的面上的坐标系，规定竖直向上为正方向，以零时刻静止水面为坐标原点，对浮子和振子分别做受力分析。

图2 浮子振子相对位移示意图及受力分析图

图3 退火法求解最优解及粒子群法求解最优解的算法流程图

浮子所受的合力写作：

考虑F附加惯性力是由于推动浮体摇荡的力还会推动附体周围的流体运动，等效于为浮体加上了一个附加质量，即垂荡附加质量Me［9～10］。可对浮子列出：

再同理对振子做受力分析，如图2 所示。阻尼力做功的平均功率，即PTO的平均输出功率为

求解方程组时运用Matlab（MathoWorks，R2017a）提供的ODE45 求解器求得该微分方程组的数值解；最终要求使平均输出功率最大化，即利用Matlab 求解优化问题。本模型求解属于全局最优问题，常用算法包括全局搜索法（GlobalSearch），粒子群算法（Particle Swarm Optimization，PSO）和模拟退火算法（Simulated Annealing，SA）。我们分别利用三种方法对函数进行求解获取使功率最大化时的阻尼系数，同时比较分析不同算法的结果［11］。利用Matlab可绘出浮子、振子的振幅变化及速度变化（图4）。

图4 浮子、振子作垂荡运动时参数变化图像

规定零时刻静止水面为坐标原点。由图像可知，带入我们假定的数据进行计算，浮子起初受到海水的作用力开始在竖直方向上振动，同时带动振子振动；大约在50s 时二者组成的系统达到平衡，以相同频率振动。在阻尼系数为常数时用三种算法计算最大输出功率和对应的阻尼系数，结果和计算所用时间如表1 所示。所用计算机处理器位Intel Core i7，内存16GB。

表1 阻尼系数为常数时三种算法的计算结果

由表中数据可看出，三种算法计算得到的最大功率为226W 左右，对应阻尼系数为37200N·s/m。而三种算法中粒子群算法所需时间最少。

2.2 优化模型

本问题中浮子不仅做垂荡运动还做纵摇运动，因此在转轴上还安装了旋转阻尼器和扭转弹簧。此时浮子与振子在二维平面上运动。

浮子在做垂荡运动的基础上还做纵摇运动，所以在之前的受力基础上还需考虑扭转弹簧和旋转阻尼器的作用以及海水的各作用力对应的旋转力矩。考虑浮子在竖直平面上的剖面，其二维质心与O点的距离为

其中，Xi为各部分质心与锥顶的距离；Mi为各部分的质量。浮子受力在基础模型类似，但需考虑到部分力与竖直方向存在夹角。由牛顿第二定律，对浮子可列出：

为分析纵摇运动时的角位移与角速度，分析浮子力矩，由刚体转动定律可列出：

其中I为浮子的总转动惯量；correction为对浮子转动惯量的扭矩修正系数。这是因为扭转弹簧与旋转阻尼器对浮子的作用力作用在转轴上，而我们在计算力矩时需将其转移到浮子的质心。扭转修正系数的表达式为

再同理分析振子。分别求解上述微分方程，获得浮子、振子的运动规律。

此时的输出功率由直线阻尼器和旋转阻尼器两部分做功组成。对于直线阻尼器和旋转阻尼器，分别有

最终的系统输出功率为两种阻尼器的功率之和。利用问题二中的求全局最优解的算法求解P最大时ζ1与ζ3的值。

三种算法求解最优阻尼系数的计算结果和用时见表2。全局搜索法和粒子群法计算的结果相近，求得的最大输出功率为325 W 左右，对应直线阻尼系数和旋转阻尼系数分别为64200N·s/m 和32200N·s/m；三种算法种粒子群法用时最短。由图5和图6可以比较得出，在此问题中，粒子群算法相较于退火法和全局搜索法具有显著优势，原因可能来自于粒子群算法的迭代次数比退火法少了一个数量级，因此会更快收敛。

表2 做垂荡、纵摇运动时三种算法求解最优阻尼系数的计算结果

图5 粒子群算法迭代次数

图6 模拟退火算法

图5 为粒子群算法运行过程中，最优结果随迭代次数的变化趋势；图6 为模拟退火算法运行时的参数情况，图6（a）为最优结果随迭代次数的变化趋势，图6（b）为实时结果随迭代次的变化趋势，图6（c）为最优结果对应的自变量取值，图6（d）为当前结果对应的自变量取值。

3 模型分析

根据文献［13］，估算波浪在单位时间内传递到浮子所在地的能量为

我们假定波幅与浮子垂荡位移约为同一数量级，根据第三问计算的浮子位移，取平衡时的波幅H 为0.5m，得到单位时间内波浪传递到浮子所在地的能量E 为1255.63 J，由于海浪上浮子的运动不止垂荡和纵摇，所以在我们模型中最为贴近实际的是运动更为复杂的改进模型，故取P 为325.12 W。计算得到效率为25.9%。

在浮子的振动已经稳定时，此时浮子做稳定的有阻尼受迫振动，任意时刻的能量不等于耗散的能量，系统的能量处于随时变化的状态［14］，波浪传递来的能量分别被兴波阻力、附加转动惯量、附加质量等系统外阻力耗散掉一部分，因此此效率比较合理。

4 结语

本文对传统的浮子波浪能模型进行修正改进，考虑了浮子作纵摇运动并建立动力学模型。利用全局搜索算法、粒子群算法与退火算法计算最优阻尼系数，通过对三种算法的运行时间、迭代次数、收敛结果的对比可以直观地观察到三种算法的表现，同时还可以互相印证，得出更为可信的结果。其中粒子群算法依靠在定义域内随机产生大量粒子完成搜索，搜索速度最快，且结果相对准确，是研究本问题的最优算法。但是粒子群算法理论上容易陷入局部最优解，所以还需在后续研究中不断优化、验证。

粒子群算法计算改进后模型得到发电效率为25.9%。由文献［2］的结果可知，10%～30%是目前浮子式波浪能装置发电的合理范围，而我们的模型效率落在此范围之间，从侧面印证我们建立了一个较合理的模型。