我国光伏发电产业发展迅速，前景十分广阔[1]，2018年底全国光伏发电规模已经达到1.74亿kW，居世界第一位。然而，土地资源的紧缺限制了其进一步发展，传统陆上光伏电站要占用巨大土地面积，而水上光伏以其发电效率高，不占用土地资源的优势受到重视[2]。水面漂浮光伏电站可以安装在水塘、水库、湖泊等水面上，缓解土地资源日益紧张的问题[3]。我国拥有广阔的近海区域及丰富的湖泊、水库等资源，具有发展水上光伏的天然优势。水上光伏以其固定形式分为固定式及漂浮式两种。固定式光伏电站采用桩基固定，出于成本考虑只适用于浅水水域。相对而言，漂浮式光伏电站具有更加广阔的发展空间。
目前水面漂浮光伏电站已经在法国、美国、英国和日本等国家相继建成[4－5]，我国于2016年在安徽省淮南市建成全球第一个单体容量达到40 MW的漂浮光伏电站——潘阳光伏40 MW水面光伏发电项目。漂浮光伏电站由浮体提供浮力漂浮在水面上，浮体采用高密度聚乙烯制作而成，目前主要有立方体式、浮体加钢支架式和浮体加小角度倾斜式三种形式，小角度倾斜式水面漂浮光伏电站如图1所示。图中左下水域为方阵西侧。为防止光伏电站在湖泊风浪流等环境载荷作用下发生漂移继而产生碰撞，需要系泊系统限制漂浮方阵的移动范围，因此，漂浮方阵所受的风、浪、流等环境载荷对其系泊设计以及方阵整体及局部的结构强度设计十分重要。
图1 漂浮式光伏电站
Fig.1 Floating photovoltaic power plant
一般而言，风载荷是光伏电站的主要环境载荷，国内外已有学者针对陆上光伏电站风载荷展开研究。王建勃等[6]对不同数量光伏组件阵列风荷载所受风载荷进行数值模拟计算，根据所得风载荷分布规律提出了设计建议。黄张裕和阎虹旭[7]通过数值模拟对光伏板风载荷系数群体遮蔽效应进行了研究，并给出了设计建议。阮辉等[8]对一低矮细长大跨度光伏阵列进行了数值模拟，分析了不同风向角下电池板表面不同分布分块区域的净载荷体形系数。Shademan和Hangan[9]采用数值模拟的方法对单体光伏板的阻力系数及其分布进行了研究，给出了不同风向角下的光伏板不同区域的阻力系数，并以三个单体模型为对象，研究了不同风向角下的遮蔽效应。Bitsuamlak等[10]采用大涡模拟方法研究了固定式光伏板的风荷载，并与试验结果进行了对比，同时分析了光伏板阵列的遮蔽效应。现有研究成果表明，采用CFD数值模拟是研究光伏板风载荷的有力手段，风载荷的分布及遮蔽效应是该问题研究的重点方向。
与陆上的固定式光伏电站不同，漂浮式光伏电站需要系泊系统来确保整体方阵的稳定安全，因此方阵整体的风浪流等环境载荷是系泊设计必不可少的载荷参数。但是，迄今对方阵整体风载荷尚未有专门的研究，对大规模漂浮方阵的流载荷及波浪载荷的预报研究也处于空白状态。
本文以前述淮南潘阳40 MW水面光伏发电项目170 m×170 m的漂浮方阵为研究对象进行环境载荷的数值预报。漂浮方阵的浮体单元数以万计，这使得通过数值模拟预报载荷变得十分困难。本文通过数值计算仿真、风洞试验测量分析等技术手段研究漂浮式光伏电站漂浮方阵在风载荷、流载荷、波浪载荷作用下漂浮方阵的受力情况，导出了载荷分布的规律和数值预报方法，为水上漂浮式电站的设计提供基础性数据。
1 漂浮方阵风载荷计算
本文采用CFD方法[11]计算漂浮方阵风载荷。在进行整体方阵计算之前，首先对单体光伏板+主浮体模块进行相关研究，包括两方面内容，其一是进行精细模型与简化模型的计算对比，其二是完成简化模型风洞试验并与计算结果对比。
对整体方阵则采用大规模计算，得到风载荷分布规律，针对性地提出2.5 D的计算策略提高计算精度。根据风载荷分布规律设计并完成方阵模型的风洞试验，对比试验和计算结果，得到漂浮方阵风载荷的计算方法。
1.1 单体模型计算及验证
主浮体的外形较复杂，采用精细模型会显著增加实尺度方阵建模难度并且需要相当多的网格单元来描述模型，因此，建模时将主浮体简化为梯形六面体，首先进行简化模型与精细模型的CFD计算对比，考察采用简化模型的可行性，然后通过风洞模型风载荷试验来验证CFD计算的精度。两种模型的形状如图2所示。
两种模型的计算方法完全相同。本文CFD计算采用通用软件FLUENT，如图3所示，计算域横截面6 m×3.5 m与风洞试验条件一致。采用旋转模型的方式完成不同风向的建模和风载荷计算，风向角定义见图4。入口边界条件为速度入口，出口为压力出口，四周边界为无滑移边界条件。计算中湍流模型采用realizable k-ε模型[12]，壁面处采用非平衡壁面函数进行模拟，速度压力耦合采用SIMPLE算法，空间离散采用标准型，其他离散项均采用一阶迎风格式。
图2 主浮体精细模型与简化模型对比
Fig.2 Comparison between the detailed model and the simplified model of the main floating body
图3 单体计算域示意图
Fig.3 Computational domain of single model
图4 单体简化模型
Fig.4 Single simplified model
计算得到各风向角下的水平风载荷，将其换算成按式(1)定义的无量纲的风载荷系数Cf形式进行对比：
式中各量均采用国际单位制，其中：FIN是水平方向受力；ρ/(kg/m3)为空气密度；A/m2为特征面积(由于从模型到实体换算得到载荷时在面积方面仅和缩尺比有关，故将特征面积取为1，换算成实体时乘以缩尺比的平方即可)；v/(m/s)为风速。本文中风、浪、流三种环境载荷系数均采用此计算方式。
图5对比了简化模型与精细模型的计算结果，二者趋势完全一致，各风向角下，简化模型的计算结果比精细模型的计算结果略大。由此表明，采用简化模型计算在不同风向角下的受力规律与精细模型基本一致，且采用简化模型进行计算偏安全，可以采用简化模型进行该问题的研究。
图5 单体风载荷系数
Fig.5 Wind load coefficients of the single model
为验证CFD计算的准确性，进行了单体简化光伏板+主浮体模块的模型试验，并与单体简化模型计算结果进行对比。单体模型按照1∶2的缩尺比进行制作，基本尺寸见表1，模型示意图及试验风向如图6所示。
表1 模型基本尺寸
Table 1 Basic dimensions of the model
单体模型 长/mm 宽/mm 短边高/mm 长边高/mm主浮体 485 420 92 189光伏板 829 496 12 20
图6 单体简化模型
Fig.6 Single simplified model
该单体模型在上海交通大学多功能风洞[13]中进行试验，试验风速9 m/s。考虑到对称性，在0°～180°区间，每隔15°进行测量，180°～360°区间设置较少测量点，试验测得各风向角下的模型所受风载荷。
试验结果见图5。对比试验结果与简化模型计算结果，二者吻合良好，趋势一致，多数风向角、特别是大载荷风向时数值接近。计算与试验的结果对比充分验证了采用CFD模拟该类型流场的准确性，为后续的CFD工作提供了基础支持。
1.2 实体方阵风载荷计算
本文研究对象为170 m×170 m的水上漂浮光伏电站方阵。该方阵包含众多单体模型，其布置方式为南北向120行，东西向96列。由于南北方向相邻两副浮体的布置方式如图7所示存在一定差异，左侧(北侧)副浮体相比右侧缺少了与主浮体相邻的部件，使相邻两主浮体的迎风面积有所不同，从而导致南北相邻两单元受力的不同。本研究中南北方向共120个单元之多，其整体载荷分布规律是本文的研究重点，而相邻两单元的载荷变化并不影响整体规律，故无单独处理的必要。因此，将南北相邻两单元合并成一个单元进行分析，后文图表中南北向为60行，按从北向南顺序排序1～60。方阵中最西端与最东端各有一列为浮体走道，没有光伏板的主浮体，在本文分析中将其各成一列合并到整体方阵中，东西向合计98列，按照从西向东的顺序排序1～98。
图7 单元示意图
Fig.7 Schematic of a unit
要得到整体方阵的风载荷，直接对其进行CFD计算是最理想的方法，但是对170 m×170 m方阵进行高精度的CFD计算所需网格单元数至少为10亿级，计算量巨大。为解决此问题，研究中采用了以下方法：首先进行粗网格简略计算，再进行局部高精度计算，最后对简略计算结果进行修正得到较为准确的数值解。
1.2.1 实体方阵简略计算
计算域如图8所示。流域四周边界距方阵外围550 m、流域顶部距水面200 m，计算网格单元总数约为4800万。主浮体与副浮体凹凸不平的表面在数值计算中通过设置表面粗糙度进行模拟，粗糙度高度设置为5 mm。
图8 实体方阵计算域
Fig.8 Computational domain of array
流域顶部及底部水面设为滑移壁面，四周边界根据风向角的变化，迎风面设为速度入口，风速为15 m/s，背风面设为压力出口。考虑到方阵东西方向的对称性，本研究进行了0°、30°、45°、60°、90°、120°、135°、150°、180°共9个风向角下的风载荷计算。其中0°为北风，90°为西风，180°为南风，以此类推。数值模拟方法与前述单体模型计算一致。
计算完成后提取方阵所有光伏板、主浮体、副浮体的总体风载荷进行对比分析，计算结果如图9所示。
图9 风载荷曲线
Fig.9 Wind load curves
整体风载荷曲线基本呈现出沿90°对称的形状，90°风向角下风载荷最小。0°与180°风载荷远大于90°，0°略大于180°，且为各风向下最大值，因此，北风为最大载荷风向。
结合光伏板、主浮体、副浮体的受力曲线，光伏板所受风载荷按风向角的变化规律与整体风载荷规律一致；主浮体与副浮体风载荷随风向角变化规律不大。90°风向角下光伏板风向投影面积极小，风载荷最小。除90°附近的其他风向角下，光伏板贡献了大部分风载荷，为主要受风构件。
1.2.2 风载荷分布规律及2.5D计算方法的提出
尽管实体方阵3D简略计算可以得到不同风向角下的风载荷变化规律及分布情况，但由于计算网格不够精细，数值解的准确度难以得到保障。根据简略计算得到的载荷分布规律，从整体方阵中选取一部分进行局部计算，再推广至整体方阵是解决此问题的核心思想。
图10为北风时方阵中风载荷分布云图，从图10可以直观感受到风载荷由迎风位置到末端的衰减趋势，并且，在中部位置风载荷沿行数的变化规律一致。为更直观分析此规律，将方阵中每一列单元的载荷进行累加，命名为列风载荷，绘制列风载荷分布曲线如图11所示。列风载荷分布曲线基本呈现对称形状，其中中间大部分区域数值接近，曲线平缓，两边缘处受力逐渐减小。
图10 0°风向角风载荷分布
Fig.10 Wind load distribution at 0° incident angle
在北风作用下，列风载荷分布在中部位置基本是一致的，只要知道了中部位置的列载荷，整个方阵上的风载荷就容易计算。基于这一规律，本研究提出2.5 D计算的方法，即针对南北向一列基本单元进行计算，其东、西两侧边界采用周期边界条件，此方法可认为是无限多列单元的计算。该方法在简化计算量的同时，保证了光伏板与主副浮体的形状完整，与实尺度方阵情况更为接近。本项目中方阵较庞大，共98列，可近似为无限列单元的计算，图11也反映出中间位置的列载荷均匀一致，理论上2.5 D的计算结果应与3 D计算中中间位置一列的结果较为接近。
图11 0°风向角列风载荷分布
Fig.11 Column wind load distribution at 0° incident angle
1.2.3 2.5D计算方法的数值验证
为验证上述设想，进行了2.5 D与3 D的计算对比，两个计算的网格密度与计算设置完全一致，从而排除了由于网格差异与计算方法不同而产生的误差。2.5 D计算域可参考图8。图12显示了3 D计算第49列的结果与2.5 D的结果，两条曲线的形状与数值大小均十分吻合，从而验证了计算北风风载荷时用2.5 D计算推广到3 D计算的可行性。2.5 D计算可以采用较精细的网格，提高风载荷的数值预报精度。
图12 2.5D与3D计算结果对比
Fig.12 Comparison between 2.5D and 3D calculations
1.3 模型方阵风载荷计算与验证
本研究对漂浮式光伏方阵的风载荷进行了风洞试验测量。根据CFD计算揭示的规律，结合风洞试验条件，确定了18行×11列共计198块单元模块(太阳能板及主浮体)的1∶6缩尺模型方阵的风洞试验方案，测量了不同风向下的特定位置单元的风载荷。与相同条件的CFD计算进行对比，可以验证针对该类型阵列风载荷计算数值方法的精度，同时可以为高精度的2.5 D计算提供计算设置及网格划分的参考。
1.3.1 模型方阵试验
试验方阵布置方式为南北向18行(北端为第1行，南端为第18行)，东西向11列(从西向东列编号为1～11)，如图13所示。试验中用测力天平对第1行、第18行和第6列各单元模块进行风载荷测量。试验风速为15 m/s，测量了该方阵模型在0°～180°、间隔22.5°共9个方向角下的各选定单元上的风载荷。
图13 试验模型方阵
Fig.13 Test model array
1.3.2 模型方阵计算
选取北风试验条件下的工况进行CFD计算。网格的划分方法与实体方阵一致，网格单元总数约1700万，如图14所示，计算域主体部分与风洞试验段尺寸一致，另外附加了一定长度的来流段与去流段。
图14 试验方阵计算域
Fig.14 Computational domain of test array
试验方阵数值模拟中，湍流模型采用可实现的k-ε模型，壁面处采用非平衡壁面函数进行模拟，速度压力耦合采用SIMPLE算法，压力项离散采用PRESTO，其他离散项均采用一阶迎风格式。
1.3.3 计算结果对比
计算完成后提取试验测试位置的浮体与光伏板模型的受力，与试验值进行对比，绘制曲线如图15。
图15 模型方阵计算与试验结果对比
Fig.15 Comparison between computed results and measured data of model array
图15(a)为中间一列的风载荷系数对比曲线，两条曲线的整体趋势及形状基本一致，两条曲线的后半部分在形状上基本吻合，且数值更为接近。相比于后半部分，前半部分两条曲线则有一定程度的偏差，尤其是迎风位置的首个单元，计算值比试验值偏小。分析其原因，本研究的太阳能模型方阵由钝体组成，空气流过时会产生复杂的旋涡结构，进而导致流场中钝体所受摩擦阻力及压差阻力情况均十分复杂，采用基于雷诺平均的RANS方程求解存在一定的局限性，采用更加精细的网格使用LES或DES方法或许可以得到较好的局部结果。
图15(b)为南北两行的风载荷系数对比结果，北第一行两条曲线的偏差已在图15(a)中有所反映，该曲线显示出这种偏差在迎风的第一行是共性的；南第一行的两条曲线则更为接近。同时，图中曲线反映的方阵载荷沿东西向的分布与初步计算的结果是一致的，即载荷分布东西向对称，中间大部分区域载荷基本一致，两边缘略有不同。
尽管在某些局部位置计算值与试验值有一定偏差，但作为关注重点的列整体风载荷系数，试验值为0.0516，计算值为0.0539，二者仅有约4.4%的相对偏差，且计算预报结果稍大于试验值，从偏于保守安全的角度考虑，这样的偏差可以接受。因此，本计算方法可以应用于整体载荷的数值预报。
1.3.4 应用2.5D计算结果修正漂浮方阵风载荷数值预报结果
采用与上述试验工况计算一致的数值求解方法及网格划分形式对2.5 D工况进行高精度计算。计算后的2.5 D整体风载荷系数为5.103，与实体方阵简略计算得到的中间第50列风载荷系数9.386做比，得到0.547的比例系数。假定高精度计算和实体方阵简略计算得到的载荷分布规律类似，用2.5 D计算得到的比例系数对简略计算的结果进行修正，即可得到较为准确的整体方阵风载荷数值结果。
2 漂浮方阵流载荷计算
流载荷的计算方法与风载荷完全一致。建模时将浮体的水下部分简化为六面体，不同部分吃水分为17 mm与50 mm两种。一个基础单元如图16所示。
图16 模型单元示意
Fig.16 Unit model grid
计算区域中浮体方阵四周计算域各向外扩张600 m，网格单元总数约为2233.6万。
图17和图18分别给出0°来流时列流载荷和行流载荷的分布，其中，列载荷分布绘出东、西两端补充浮体流载荷，行载荷分布中首行数值包含北端补充浮体流载荷。方阵96列列流载荷分布规律表现为中间部分较高且在一段区域内平稳，边缘较中间部分略微降低。行流载荷分布总体上可分为两段，前20行组成流载荷下降段，下降段之后流载荷趋于平稳，为平稳段。
图17 0°来流列流载荷分布
Fig.17 Column load distribution at 0° incident current
图18 0°来流行流载荷分布
Fig.18 Row load distribution at 0° incident current
图19为以各单元编号坐标及对应流载荷大小绘出的流载荷分布图，中间曲线为2.5 D计算结果。可见，2.5 D计算结果可以较准确地预报0°流向角下方阵中间部分的行载荷分布规律。
图19 0°来流时流载荷三维分布及与2.5 D计算结果对比
Fig.19 Distribution of current loads of 3D computation and comparison with 2.5 D computed results
对0°来流0.183 m/s、0.5 m/s、1 m/s、1.5 m/s流速下2.5 D南北一列单元流载荷系数进行计算，结果见表2。
表2 不同流速2.5D流载荷系数
Table 2 2.5D current load coefficient at various current speeds
流速/(m/s) 2.5D流载荷系数0.183 1.901 0.5 1.713 1.0 1.687 1.5 1.673
不同流速下的2.5 D流载荷系数变化基本平稳，仅在0.183 m/s的低流速下受到摩擦阻力特性的影响，流载荷系数稍有增大。不过，此时流载荷数值上是小量，对工程问题的研判影响不大。因此，流速大于0.5 m/s即可认为流载荷与流速平方成正比，流速较低时，适当提高流载荷系数也可快速预报流载荷。
3 漂浮方阵波浪载荷计算
目前，国内外针对水面漂浮电站的波浪载荷计算未有先例。本文基于势流理论方法[14]，采取船舶与海洋工程通用软件SESAM对漂浮方阵波浪载荷进行计算，在鄱阳湖环境条件下，针对170 m×170 m大型漂浮方阵开展波浪载荷数值分析研究，研究单位波幅规则波作用下漂浮方阵的波浪载荷随行列的变化规律，根据此规律获得了整体漂浮方阵的波浪载荷，并据此计算了50年一遇的极端条件下漂浮方阵所受的波浪载荷。
以漂浮方阵几何中心为原点，x轴由北指向南，y轴由西指向东，z轴铅垂向上。系泊系统在垂直方向可以靠重力与浮力的平衡提供恢复力，在水平方向内方向必须通过系泊系统的锚链提供的锚力来限制方阵的位移。漂浮方阵水平方向内所受的波浪力会使其发生漂移，而绕z轴的力矩会使漂浮方阵发生旋转。为了给系泊设计提供基础性数据，对漂浮方阵在波浪作用下所受的最大波浪力Fx、Fy以及波浪力矩Mz进行计算及分析研究，其中波浪力矩Mz以漂浮方阵几何中心为力矩中心。
参考鄱阳湖历年风浪统计[15]，选取的波浪周期T范围为1.5 s～5.2 s(圆频率ω范围为1.2 rad/s～4.2 rad/s)，平均水深为6m，计算的浪向为0°、45°和90°。
3.1 网格依赖性研究
同计算风载荷和流载荷一样，将两个主浮体及相连的副浮体作为一个基本单元，通过简单复制就可得到整个漂浮方阵。为了能够保证计算精度的同时减少计算量，首先研究一个基本单元的网格数对波浪力Fx的影响，计算结果见图20。
图20 网格依赖性研究
Fig.20 Study on grid resolution
由图20可知，布置38个面网格的计算结果偏大，布置67和97个面网格的结果很接近，表明布置67个面网格计算精度已足够，于是采用一个基本单元上布置67个面网格的网格方案进行计算。
3.2 漂浮方阵随行列的变化规律
在波浪的周期范围内计算在单位波幅条件下的最大波浪力，研究漂浮方阵所受最大波浪力随行数变化的规律，固定列数为32，行数分别为2、4、6、8、10、12、14、16和18时，方阵所受最大波浪力随行数的变化见图21。同样，研究漂浮方阵所受最大波浪力随列数变化的规律，固定行数为32，列数分别为2、4、6、8、10、12、14、16和18时，方阵所受波浪力随列数的变化见图22。
图21 最大波浪力随行数变化
Fig.21 Wave loads as the number of rows
由图21可知，随着行数的增加，整体方阵所受波浪力的合力在21 kN附近震荡，整体方阵所受波浪力的合力几乎不再变化；由图22可知，整体方阵所受波浪力的合力随方阵的列数呈线性增加。
图22 最大波浪力随列数变化
Fig.22 Wave loads as the number of columns
3.3 整体漂浮方阵波浪载荷计算
由于漂浮方阵由上万个浮体组成，一个基本单元上布置67个面网格这方案下计算120×96的整个漂浮方阵则需要36万个面网格，软件网格总数限制为1.5万个面网格，最多能计算25×20的方阵，无法直接计算。
由上述研究的波浪载荷随行列的变化规律可知，定义自变量n，规定5×4方阵的自变量为1，则波浪载荷Fx和Fy与自变量n成线性关系，方阵5×4、10×8、15×12、20×16、25×20和120×96对应的自变量分别为1、2、3、4、5和24。采用以下方法来计算120×96方阵在单位波幅条件下所受的波浪载荷：根据5×4、10×8、15×12、20×16、25×20方阵在0°、45°和90°浪向所受最大波浪载荷来推算出的公式来计算120×96方阵所受的波浪载荷。
浪向0°时，波浪载荷计算见图23。浪向90°时，波浪载荷计算见图24。浪向45°时，波浪载荷计算见图25、图26和图27。
图23 0°浪向时波浪力Fx
Fig.23 Computed Fx in 0° wave direction
图24 90°浪向时波浪力Fy
Fig.24 Computed Fy in 90° wave direction
图25 45°浪向时波浪力Fx
Fig.25 Computed Fx in 45° wave direction
图26 45°浪向时波浪力Fy
Fig.26 Computed Fy in 45° wave direction
图27 45°浪向时波浪力矩Mz
Fig.27 Computed Mz in 45° wave direction
由图23～图27中的公式得到单位波幅下120×96方阵所受波浪载荷，结果见表3。
表3 单位波幅下方阵所受波浪载荷
Table 3 Wave loads on the floating array at incident wave with unit amplitude
浪向 Fx/kN Fy/kN Mz/(kN•m)0° 56.16 — —45° 43.49 45.01 7.41×102 90° — 72.75 —
可见浪向为90°时，漂浮方阵所受波浪力最大；浪向为0°时，漂浮方阵所受波浪力最小；浪向为45°时，漂浮方阵所受波浪力没有正浪时大，但会受到绕z轴的旋转力矩。
3.4 极限条件下方阵所受波浪载荷
目前进行湖泊的风浪特性预报有三种方法[16]：一是现场观测方法，二是经验公式计算，三是风浪生成数学模型方法。从能量的角度，极限条件下的湖泊波浪主要是由风能的传递引起的，由于没有实测数据，本文采用经验公式方法，用50年一遇的风作为输入来计算50年一遇的波浪。
在波浪问题中，对波浪的统计特征值有特殊的意义的两个波高分别为有义波高H1/3和最大波高H1/10，通常目测的波高接近有义波高，有义波高H1/3与最大波高H1/10的关系[17]：
采用《海港水文规范》推荐的方法[18]计算有义波高：
式中：g为重力加速度；U为风速；F为风区长度；h为水深。
对于本文，当重现期为50年一遇时，漂浮方阵所在湖泊的最大设计风速U=30 m/s，相当于11级风力[19]，在鄱阳湖环境条件下风区长度F=38000 m，平均水深h=6 m，由式(2)和式(3)可得有义波高H1/3=1.62 m和最大波高H1/10=2.06 m。
根据线性波浪理论可得50年一遇极限条件下方阵所受波浪载荷见表4。
表4 极限条件下方阵所受波浪载荷
Table 4 Wave loads on the floating array under the extreme conditions
浪向 Fx/kN Fy/kN Mz /(kN•m)0° 57.84 — —45° 44.79 46.36 7.63×102 90° — 74.93 —
4 漂浮方阵载荷数值预报方法小结
本文用数值方法分别计算漂浮方阵的风载荷、流载荷和波浪载荷。其中，风速是主要因素，来自气象预报等资料；考虑湖泊条件，流为风生流，其值一般较低，如果有洪水因素，则风生流公式不适用，而用别的方法得到了流速后，计算流载荷的方法依然有效；波浪载荷计算的主要参数之一是波高，其值由内湖风生波经验公式计算，对其他情况不适用。
针对漂浮式光伏电站方阵庞大而复杂的特点，为得到方阵整体风载荷，可参考图28所示流程进行风载荷计算，条件允许的情况下建议进行试验验证。首先，对单体光伏板+浮体模型进行数值计算与试验对比，验证采用简化模型进行计算的可靠性及CFD方法运用于此类问题研究的可行性。其次，对实体方阵进行简略CFD计算，得到整体风载荷在不同风向角下的分布规律，分析风载荷最大时对应的北风条件下的载荷分布规律，按2.5 D局部计算的概念进行数值计算及验证。最后，进行模型方阵的试验与计算对比，根据所确定数值方法进行高精度的2.5 D计算，得到针对简略计算的修正系数，进而获得较为准确的整体方阵风载荷。
流载荷与风载荷类似，计算方法可完全参考图28所示流程。
图28 风载荷计算流程图
Fig.28 Flow chart of wind load's computations
对漂浮方阵的波浪载荷进行数值预报可采用势流理论进行计算。针对计算单元数过多的问题，可研究波浪载荷随浮体数目变化的规律，进而得到漂浮方阵在不同浪向角下的波浪载荷。同样根据此规律，可获得极端条件下漂浮方阵所受的波浪载荷。
本文所示方法适用于本研究对应的淮南漂浮式光伏电站的浮体及太阳能板，当相关条件，如浮体形状、太阳能板角度、浮体间距等发生一定程度的变化时，推荐尝试以下办法计算风浪流载荷：
1) 假定风载荷和流载荷沿行、列的分布规律基本不变。
2) 对北风2.5 D建模计算，网格密度和计算方法参考本文相关段落。
3) 参考本文方法推算风载荷和流载荷。
4) 波浪载荷用小维数(较小的行数和列数)进行建模计算，再外推至整体方阵。
5 结论
本文综合风洞试验和数值计算，对水上漂浮式光伏电站整体风载荷、流载荷及波浪载荷进行了数值研究，为解决类似大规模计算问题提供了新的方案。采用CFD方法对170 m×170 m三维实体漂浮方阵在不同流向角下所受风载荷和流载荷进行计算和分析，获得了载荷的分布规律；为解决完整3 D计算复杂性带来的精度不足问题，提出了2.5 D计算方法并进行了验证；采取势流计算方法，参照鄱阳湖环境条件，对漂浮方阵进行了波浪载荷数值计算与分析，研究了在单位波幅规则波作用下漂浮方阵的波浪载荷随行和列的变化规律，并据此获得了整体漂浮方阵的波浪载荷。
本研究为大型漂浮方阵的环境载荷的数值预报提供了完整的解决方案。相关研究结论如下：
(1) CFD计算结果与风洞模型试验结果吻合良好，表明基于RANS方程的CFD方法可用于此类问题的风载荷和流载荷计算。
(2) 采用简化主浮体模型计算所得风载荷稍大于原始模型，不同风向下受力规律一致。
(3) 在各风向角下，北风为最大载荷风向，西风最小；光伏板为主要受风构件。
(4) 对北风这一典型工况，采用周期边界条件对单独一列组块进行计算，与方阵整体计算中间一列所得风载荷结果一致。
(5) 不同流速的计算结果表明在一定流速范围内，无量纲的流载荷系数与流速基本无关。
(6) 在由北向南行进的规则波中，最北端的 4行浮体单元承受了几乎全部波浪载荷；整体方阵所受波浪载荷与方阵列数呈线性关系。
[1]沈浩.光伏发电现状与发展综述[J].现代制造,2016(30): 71―73.Shen Hao.Overview of the status and development of photovoltaic power generation[J].Maschinen Markt,2016(30): 71―73.(in Chinese)
[2]陈东坡.我国水上光伏电站的新机遇、新发展和新挑战[J].电子产品世界, 2017(5): 3―5.Chen Dongpo.New opportunities, developments and challenges of Chinese photovoltaic power plants[J].Electronic Engineering ＆ Product World for Engineering Managers ＆ Designers, 2017(5): 3―5.(in Chinese)
[3]王方毓.水上太阳能光伏电站的技术特点及应用[J].工程技术研究, 2017(10): 76―77.Wang Fangyu.Technology and application of solar photovoltaic power plants[J].Engineering and Technological Research, 2017(10): 76―77.(in Chinese)
[4]Cabrera-Tobar A, Bullich-Massague E, Aragues-Penalba M, et al.Topologies for large scale photovoltaic power plants[J].Renewable ＆ Sustainable Energy Reviews,2016, 59: 309―319.
[5]Branker K, Pathak M J M, Pearce J M.A review of solar photovoltaic levelized cost of electricity[J].Renewable and Sustainable Energy Reviews, 2011, 15(9): 4470―4482.
[6]王建勃, 朱锐, 何惧.光伏电站阵列风荷载衰减特性数值模拟[J].太阳能, 2013(15): 20―22.Wang Jianbo, Zhu Rui, He Ju.Numerical simulation of wind load attenuation characteristics of PV array[J].Solar Energy, 2013 (15): 20―22.(in Chinese)
[7]黄张裕, 阎虹旭.太阳能光伏板风荷载体型系数群体遮挡效应数值模拟研究[J].特种结构, 2015(3): 2015,32(3): 18―22.Huang Zhangyu, Yan Hongxu.Numerical simulation study on the group shielding effect of wind power load coefficient of solar photovoltaic panels[J].Special Structures, 2015(3): 2015, 32(3):18―22.(in Chinese)
[8]阮辉, 廖伟丽, 王康生,等.光伏阵列表面风荷载数值研究[J].太阳能学报, 2015, 36(4): 871―877.Ruan Hui, Liao Weili, Wang Kangsheng, et al.Numerical research on surface wind load of PV array[J].Acta Energiae Solaris Sinica, 2015, 36(4): 871―877.(in Chinese)
[9]Shademan M, Hangan H.Wind loading on solar panels at different inclination angles[C]// Conference of American Society of Wind Engineers, 2009.
[10]Bitsuamlak G, Dangew A, Erwin J.Evaluation of wind loads on solar panel modules using CFD[C].Proceedings of the 5th International Symposium on Computational Wind Engineering (CWE2010), North Carolina, USA:Chapel Hill, 2010.
[11]Anderson J J D.Computational fluid dynamics: The basic with application[M].McGraw-Hill Book Company Europe : McGraw Hill, 1995.
[12]Shih T H, Lou W W, Shabir A, et al.A new eddy viscosity model for high Reynolds number turbulent flow model development and validation[J].Computers ＆ Fluids,1995, 24(3): 227―238.
[13]陈作钢, 李金成, 代燚, 等.多功能风洞及CFD优化设计[J].实验流体力学, 2012, 26(4): 71―76.Chen Zuogang, Li Jincheng, Dai Yi, et al.Versatile wind tunnel and CFD-based optimal design[J].Journal of Experiments in Fluid Mechanics, 2012, 26 (4): 71―76.(in Chinese)
[14]朱伟亮, 杨庆山.基于LES模型的近地脉动风场数值模拟[J].工程力学, 2010, 27(9): 17―21.Zhu Weiliang, Yang Qingshan.Large eddy simulation of near ground turbulent wind field[J].Engineering Mechanics, 2010, 27(9): 17―21.(in Chinese)
[15]戴遗山, 段文洋, 船舶在波浪中运动的势流理论[M].北京: 国防工业出版社, 2008.Dai Yishan, Duan Wenyang.Potential flow theory of ship motions in waves[M].Beijing: National Defense Industry Press, 2008.(in Chinese)
[16]程时长, 李良文.鄱阳湖的风情及风浪特性浅析[J].江西: 江西水利科技, 1993, 4: 007.Cheng Shichang, Li Liangwen.Analysis of characteristics of the wind and wave in Poyang Lake[J].Jiangxi: Jiangxi Hydraulic Science ＆ Technology, 1993, 4: 007.(in Chinese)
[17]张洪生, 辜俊波, 文武健.淀山湖风浪场的数值模拟[J].湖泊科学, 2011, 23(5): 783―788.Zhang Hongsheng, Gu Junbo, Wen Wujian.Numerical simulation of the wind wave fields in Dianshan Lake[J].Journal of Lake Sciences, 2011, 23(5): 783―788.(in Chinese)
[18]Faltinsen O M.Sea loads on ships and offshore structures[M].United Kingdom: Cambridge University Press, 1990.17―22.
[19]JTS145-2―2013, 海港水文规范[S].2013.JTS145-2―2013, Code for Sea port hydrology[S].2013.(in Chinese)
[20]徐祝龄.气象学[M].北京: 高等教育出版社, 1992:32―33.Xu Zhuling.Meteorology[M].Beijing: Higher Education Press, 1992: 32―33.(in Chinese)