一种风电场代表年风速确定方法及系统与流程

1.本发明涉及风力发电技术领域，特别是涉及一种风电场代表年风速确定方法及系统。


背景技术：

2.本部分的陈述仅仅是提供了与本发明相关的背景技术信息，不必然构成在先技术。
3.建设风电场最基本的条件是要有能量丰富、风向稳定的风能资源，风电场选址时应尽量选择风能资源丰富的地点。所以对风电场风能资源的评估是整个风电场建设、运行的重要环节。
4.在风资源评估中，计算代表年的目的是通过风场附近气象站近多年的年平均风速的变化规律及风场与气象站实测风速相关关系，分析计算出一套能够反映风电场今后若干年运行期间平均水平的代表年的风速系列、计算出风电场今后运行若干年的平均发电量及上网电价，其中风速与风功率密度两个指标是风电场的可行性分析中重要的指标。
5.风电场风资源评估一般遵循以下步骤：
6.(1)在风电场中选取对风电场代表性良好的点设立测风塔；
7.(2)待测风期满一周年后，核实数据完整率，收集附近气象站同期气象资料及30年以上年平均风速；
8.(3)将风场短期测风数据订正为代表年风速数据；其方法如下：
9.1)作测风塔与对应年份气象站各风向象限风速相关曲线。
10.某一风向象限内风速相关曲线的具体做法是：建一直角坐标系，横坐标轴为气象站风速，纵坐标轴为测风塔风速。取测风塔在该象限内的某一风速值(某一风速值在一个风向象限内一般有许多个，分别出现在不同时刻)为纵坐标，找出气象站对应时刻的风速值(这些风速值不一定相同，风向也不一定与风电场测风塔风向相对应)，求其平均值作为横坐标即可定出相关曲线的一个点。对测风塔在该象限内的其余每一个风速重复上述过程，就可得到这一象限内的风速相关曲线。对其余各象限重复上述过程，可获得16个测风塔与气象站的风速相关曲线。
11.2)对每个风速相关曲线，在横坐标轴上表明气象站多年年平均风速，以及与测风塔观测同期的气象站的年平均风速，然后在纵坐标轴上找到对应的风电场测风塔的两个风速值，并求出这两个风速值得代数差值(共有16个代数差值)。
12.3)测风塔数据的各个风向象限内的每个风速都加上对应的风速代数差值，即可获得订正后的测风塔风速资料。
13.(4)将订正后的风速资料代入wasp软件中，计算其发电量，从而评估风电场经济效益，核实其可行性。
14.上述步骤(3)中的方法适用于绝大多数风电场，但某些风电场距离气象站较远或者受到地形影响，测风塔和气象站数据的相关性并不好，此时则无法应用上述规范中方法
进行订正；另外，若测风年为大风年，采用上述方法进行订正时可能出现风速为负的情况。


技术实现要素：

15.为了解决上述问题，本发明提出了一种风电场代表年风速确定方法及系统，通过分析长期气象站与测风塔的关联性以确定代表年风速，解决气象站与测风塔相关性较差时代表年风速的确定问题。
16.为了实现上述目的，本发明采用如下技术方案：
17.第一方面，本发明提供一种风电场代表年风速确定方法，包括：
18.获取测风塔的测风年风速和气象站的测风风速，所述气象站的测风风速包括历史年风速和测风年风速；
19.其中，历史年风速为与气象站的多年平均风速相同年份的风速，测风年风速为与测风年同期的风速；
20.对测风塔的测风年风速、气象站历史年风速与气象站测风年风速分别进行边缘分布曲线拟合；
21.根据测风塔的测风年风速与气象站历史年风速和气象站测风年风速得到秩相关系数；
22.根据秩相关系数、测风塔的拟合结果和气象站的拟合结果构建测风年连接函数和历史年连接函数，以分别构建测风塔在测风年与代表年的风速分布图；
23.根据测风塔在测风年与代表年的风速差得到风速订正值，以此修订代表年风速。
24.作为可选择的实施方式，所述秩相关系数τ为：
[0025][0026]
其中，sign()是符号函数，当(x
i-xj)(y
i-yj)》0时，sign＝1；(x
i-xj)(y
i-yj)《0时，sign＝-1；(x
i-xj)(y
i-yj)＝0时，sign＝0；xi、xj为第i时刻和第j时刻测风塔测风年逐时风速，yi、yj为第i时刻和第j时刻气象站历史年逐时风速或测风年逐时风速。
[0027]
作为可选择的实施方式，构建测风塔和气象站连接函数的过程中，包括：根据秩相关系数τ得到copula函数的相关性参数θ，具体为：
[0028]
作为可选择的实施方式，构建测风塔和气象站连接函数的过程中，包括：根据相关性参数、测风塔测风年风速的边缘分布、气象站测风年风速的边缘分布计算测风年的copula函数；根据相关性参数、测风塔测风年风速的边缘分布、气象站历史年风速的边缘分布计算历史年的copula函数。
[0029]
作为可选择的实施方式，测风年连接函数和历史年连接函数为：
[0030]
c(u1,u2,
…
,un)＝exp(-((-lnu1)
θ
+(-lnu2)
θ
+
…
+(-lnun)
θ
)
1/θ
),θ≥1
[0031]
其中，u1为测风塔测风年风速的边缘分布，u2……
un为每个气象站测风年风速或历史年风速的边缘分布，θ为相关性参数。
[0032]
作为可选择的实施方式，构建测风塔在代表年与测风年的风速分布图的过程中，包括：根据气象站的历史年风速和测风年风速分别得到气象站的历史年平均风速和测风年平均风速，根据气象站的历史年平均风速采用历史年连接函数构建测风塔在代表年的风速
分布图，根据气象站的测风年平均风速采用测风年连接函数构建测风塔在测风年的风速分布图。
[0033]
作为可选择的实施方式，测风塔在测风年的风速分布为：
[0034][0035]
其中，u1为测风塔测风年风速的边缘分布，u2……
un为每个气象站测风年风速或历史年风速的边缘分布，u2……
un代表每个气象站历史年平均风速或气象站测风年平均风速的边缘概率。
[0036]
第二方面，本发明提供一种风电场代表年风速确定系统，包括：
[0037]
风速获取模块，被配置为获取测风塔的测风年风速和气象站的测风风速，所述气象站的测风风速包括气象站的历史年风速和测风年风速；
[0038]
边缘分布模块，被配置为对测风塔的测风年风速、气象站历史年风速与气象站测风年风速分别进行边缘分布曲线拟合；
[0039]
相关计算模块，被配置为根据测风塔的测风年风速与气象站历史年风速和气象站测风年风速得到秩相关系数；
[0040]
风速分布模块，被配置为根据秩相关系数、测风塔的拟合结果和气象站的拟合结果构建测风年连接函数和历史年连接函数，以分别构建测风塔在测风年与代表年的风速分布图；
[0041]
风速修订模块，被配置为根据测风塔在测风年与代表年的风速差得到风速订正值，以此修订代表年风速。
[0042]
第三方面，本发明提供一种电子设备，包括存储器和处理器以及存储在存储器上并在处理器上运行的计算机指令，所述计算机指令被处理器运行时，完成第一方面所述的方法。
[0043]
第四方面，本发明提供一种计算机可读存储介质，用于存储计算机指令，所述计算机指令被处理器执行时，完成第一方面所述的方法。
[0044]
与现有技术相比，本发明的有益效果为：
[0045]
常规方法处理只采用相关系数表征气象站与测风塔风速之间关系，单个系数值不能真正反映关联性，且气象站只能选取距离风电场最近的气象站。而本发明提出一种风电场代表年风速确定方法及系统，采用copula函数评估气象站与测风塔风速之间关联性，从而计算风电场代表年风速；气象站不再局限于一个气象站，可以在任意气象站风速条件下计算测风塔风速的条件分布，可以根据多个气象站的数据共同确定风电场代表年风速。
[0046]
本发明提出一种风电场代表年风速确定方法及系统，通过秩相关系数、气象站与测风塔风速的联合分布，以及在任意气象站风速条件下测风塔风速的条件分布更深层次剖析气象站与测风塔风速的关联性，结果真实可信，在实际应用中更具有效性。
[0047]
本发明附加方面的优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本发明的实践了解到。
附图说明
[0048]
构成本发明的一部分的说明书附图用来提供对本发明的进一步理解，本发明的示
意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。
[0049]
图1为本发明实施例1提供的风电场代表年风速确定方法流程图；
[0050]
图2为本发明实施例1提供的测风塔测风年风速的边缘分布图；
[0051]
图3为本发明实施例1提供的测风塔代表年与测风年的风速条件分布图；
[0052]
图4为本发明实施例1提供的测风年气象站与测风塔的逐时风速散点图；
[0053]
图5(a)-5(b)为本发明实施例1中的方法与现有方法结果对比图。
具体实施方式
[0054]
下面结合附图与实施例对本发明做进一步说明。
[0055]
应该指出，以下详细说明都是示例性的，旨在对本发明提供进一步的说明。除非另有指明，本文使用的所有技术和科学术语具有与本发明所属技术领域的普通技术人员通常理解的相同含义。
[0056]
需要注意的是，这里所使用的术语仅是为了描述具体实施方式，而非意图限制根据本发明的示例性实施方式。如在这里所使用的，除非上下文另外明确指出，否则单数形式也意图包括复数形式，此外，还应当理解的是，术语“包括”和“具有”以及他们的任何变形，意图在于覆盖不排他的包含，例如，包含了一系列步骤或单元的过程、方法、系统、产品或设备不必限于清楚地列出的那些步骤或单元，而是可包括没有清楚地列出的或对于这些过程、方法、产品或设备固有的其它步骤或单元。
[0057]
在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。
[0058]
实施例1
[0059]
如图1所示，本实施例提供一种风电场代表年风速确定方法，通过分析长期气象站与测风塔的关联性以确定代表年风速，解决气象站与测风塔相关性较差时，代表年风速的获取与修订，为风电场的开发提供科学决策依据，包括以下步骤：
[0060]
获取测风塔的测风年风速和气象站的测风风速，所述气象站的测风风速包括气象站的历史年风速和测风年风速；
[0061]
对测风塔的测风年风速和气象站的测风风速分别进行边缘分布曲线拟合；
[0062]
根据测风塔的测风年风速与气象站的测风风速得到秩相关系数；
[0063]
根据秩相关系数、测风塔的拟合结果和气象站的拟合结果构建测风年连接函数和历史年连接函数，以分别构建测风塔在测风年与代表年的风速分布图；
[0064]
根据测风塔在测风年与代表年的风速差得到风速订正值，以此修订代表年风速。
[0065]
下面对本实施例的风电场代表年风速确定方法进行详细阐述。
[0066]
s1：在风电场区域设立测风塔，并获取测风塔的测风年风速；
[0067]
在本实施例中，通过在风电场内选取有代表性的地点设立测风塔，测风塔包括风速仪、风向仪、温度传感器和压力传感器等设备，在测风年进行测风并获取测风年风速。
[0068]
s2：获取周边或中尺度气象站的测风风速，气象站的测风风速包括历史年风速和测风年风速；
[0069]
s3：判断气象站的测风风速与测风塔的测风年风速的相关性，若气象站的测风风速与测风塔的测风年风速相关性超过阈值，则根据其相关性确定代表年风速；可以理解的，该过程采用现有方法即可，本实施例不做重点介绍。
[0070]
s4：若气象站的测风风速与测风塔的测风年风速相关性低于阈值，或气象站为多个时，采用下述方法计算代表年风速。
[0071]
s401：计算测风塔的测风年风速、气象站历史年测风风速和气象站测风年的边缘分布，进行曲线拟合，如图2所示；
[0072]
所述边缘分布包括正态分布、伽马分布、指数分布、威布尔分布、对数正态分布和logistic分布等；采用标准差mse评估边缘分布的拟合度，标准差mse值越小，说明边缘分布拟合越好。
[0073]
标准差mse表示为：
[0074][0075]
其中，e为期望值；n为样本总数；xc(i)为边缘分布的第i个计算值；x0(i)为第i个观测值，观测值由经验概率表示，随机变量的经验不超过概率评估，采用gringorten公式：
[0076][0077]
其中，k为升序排列中的第k个观测值，n为样本大小。
[0078]
s402：根据测风塔的测风年风速与气象站历史年或测风年测风风速计算秩相关系数τ；
[0079][0080]
其中，测风年为8760个小时，xi，xj代表测风塔测风年逐时风速，yi，yj代表气象站历史年或测风年逐时风速，i，j代表不同的时间；sign()是符号函数，当(x
i-xj)(y
i-yj)》0时，sign＝1；(x
i-xj)(y
i-yj)《0时，sign＝-1；(x
i-xj)(y
i-yj)＝0时，sign＝0。
[0081]
s403：根据步骤s401得到的边缘分布和步骤s402得到的秩相关系数τ，计算测风塔和气象站在测风年与历史年的copula函数；
[0082]
copula函数是定义域为[0,1]均匀分布的多维联合分布函数，可以将多个随机变量的边际分布连接起来构造联合分布。优点在于：随机变量之间的相关性可以通过kendall’s表征；在不用假设随机变量边缘分布类型条件下可建立变量之间的联合分布。copula函数包括椭圆型、二次型、archimedean型等，本实施例采用archimedean型中的对称gumbel-hougaardcopula函数。
[0083]
在本实施例中，秩相关系数τ与copula函数的相关性参数θ存在函数关系：
[0084][0085]
根据步骤s402得到秩相关系数τ计算相关性参数θ，根据相关性参数θ、测风塔测风年风速的边缘分布、气象站测风年风速的边缘分布计算测风年的copula函数；根据相关性参数θ、测风塔测风年风速的边缘分布、气象站历史年风速的边缘分布计算历史年的copula函数；
[0086]
c(u1,u2,
…
,un)＝exp(-((-lnu1)
θ
+(-lnu2)
θ
+
…
+(-lnun)
θ
)
1/θ
),θ≥1(5)
[0087]
其中，u1为测风塔测风年风速的边缘分布，u2……
un为每个气象站测风年风速或历
史年风速的边缘分布，θ为相关性参数。
[0088]
给定u2＝u2,
…
,un＝un情况下，u1的条件概率分布为：
[0089][0090]
式中，f(u1|u2＝u2,
…
,un＝un)为测风塔风速条件分布，u2……
un为每个气象站测风年或历史年风速的边缘分布，u2……
un代表每个气象站给定的特定数值的边缘概率，即历史年平均风速或气象站测风年平均风速的边缘概率。
[0091]
根据气象站的历史年风速和气象站的测风年风速分别得到气象站的历史年平均风速和气象站的测风年平均风速；根据气象站的历史年平均风速和测风年平均风速，采用公式(6)得到两条条件概率曲线，如图3所示，即为测风塔在代表年与测风年的风速分布图。
[0092]
s404：根据两条条件概率曲线的差值得到代表年风速订正值，以此修订代表年风速。
[0093]
在本实施例中，根据修订后的代表年风速可以计算风电场风资源情况，计算发电量，解决气象站与测风塔数据相关性不好时的代表年风速计算过程，且可用于多气象站订正，所得数据准确。
[0094]
本实施例的上述步骤s1-s3均采用常规方法即可，下面对步骤s4进行详细说明，但是本实施例的保护范围不局限于所述实施例。
[0095]
在风电场采用一台70m高的测风塔进行测风工作，分别在70m、50m、10m高度各安装了一个风速仪，并在70m、10m高度各安装了一个风向仪，10m高度安装了温度计和压力计，经检验，其现场实测数据有效完整率为92.6％，满足规范要求。
[0096]
获取距离风电场最近的气象站最近30年的年平均风速，气象站与测风塔风速统计值如表1所示；
[0097]
表1气象站与测风塔风速统计表
[0098][0099]
可以看出，气象站近30年的平均风速为2.62m/s，气象站测风年风速为2.93m/s，所以测风年为大风年，根据现场的测风资料与气象站同期的风速、风向小时记录进行相关分析，合理推算风电场代表年的风能要素。按照《风电场风能资源评估方法》中的方法，计算得到风电场70m处风速与气象站风速各个扇区相关性如表2所示；
[0100]
表2各扇区相关分析结果表
[0101][0102][0103]
可以看出，nne、ne、sse、nw四个扇区的相关系数较好，均在0.8以上，可以直接采用相关方程进行订正，而se、wsw、w三个扇区的相关系数较差，均在0.5以下，若用此方法计算，得出的订正风速误差较大；又因为测风年为大风年，比平风年高出约0.3m/s，16个扇区的一次项系数均大于1，特别是ssw扇区一次项系数为3.26，相当于每个ssw风向的风速均减去3.26*0.3≈1m/s，也就是说若测风年风速在0.9m/s以下，得出的代表年风速为负，这明显是不合乎规律的，而本实施例中的技术方法则可以避免这种情况的发生。
[0104]
气象站的年平均风速以及测风塔70m处逐时风速的边缘分布采用mse标准评估结果如表3，结果表明二者均可用威布尔分布进行描述。
[0105]
表3不同分布类型对应的mse值
[0106][0107]
根据秩相关系数τ，得到相关性参数θ＝1.40；根据最大似然法，得到相关性参数θ＝1.51，二者相差不大，本实施例中采用秩相关系数τ的计算结果。
[0108]
气象站的历史年平均风速为2.6m/s，测风年平均风速为2.9m/s，根据威布尔分布参数得到平均风速为2.6m/s时其边缘分布为0.482，风速为2.9m/s时其边缘分布为0.694；所谓条件copula函数计算，就是计算当风速分别为2.6m/s和2.9m/s时测风塔的风速分布。
[0109]
其中，气象站的历史年平均风速为代表年风速，即年平均风速为2.6m/s时的风速即为代表年风速，年平均风速为2.9m/s时的风速即为测风年风速，则相同概率下的差即为代表年风速订正值。
[0110]
两种方法计算得到的代表年月平均风速及风功率密度如图5(a)-5(b)所示，该风场由《风电场风能资源评估方法》计算得到风功率密度为276w/m2，采用copula方法计算的风功率密度为251w/m2，二者相差不大，说明了此方法的适用性。
[0111]
为了分析两种方法的计算误差，以某风电场为例进行说明，共收集到测风塔周围两个气象站的数据。由于缺少长系列观测资料，对各月数据进行订正。测风塔、两个气象站的月平均风速如表4；
[0112]
表4月平均风速表单位：m/s
[0113][0114]
测风塔年平均风速为6.62m/s，风功率密度为340.83w/m2；气象站1年平均风速为2.44m/s；气象站2年平均风速为2.61m/s。若为理想状态，则计算得出的每个月的平均风速均为6.62m/s左右，因此，将其作为真实值与两种方法得到的计算值进行对比分析。分别计算各月测风塔与气象站的秩相关系数如表5；
[0115]
表5各月份kendallτ秩相关系数值
[0116][0117][0118]
可知，各月的秩相关系数均不一致，其中四月份两个气象站相关系数基本相等，而8月份却相差较大，气象站1仅为0.0191，即气象站与测风塔之间相关性极差。
[0119]
在求解两个气象站的边缘分布时，综合考虑逐时风速分布与逐月风速分布情况两种，以秩相关系数作为分配系数，得到其边缘分布如表6。
[0120]
表6各月气象站平均风速频率表
[0121][0122]
同样，不同月份测风塔的威布尔参数也不一致，如表7。
[0123]
表7各月测风塔威布尔参数值
[0130][0131][0132]
计算五种计算方法相对于理想值的rmse值，如表10和表11。
[0133]
表10五种方法计算结果与理想值风速rmse分析
[0134]
[0135][0136]
表11五种方法计算结果与理想值风功率密度rmse分析
[0137][0138]
[0139]
由表10和11可知，在五种方法中rmse无论是风速还是风功率密度，均是copula方法两个气象站均考虑时rmse值最小，也就是其最接近真实值，而规范方法气象站1计算结果rmse均为最大，也就是最偏离真实值。
[0140]
针对气象站与测风塔相关性较差以及可能出现负风速的缺点，本实施例提出一种风电场代表年风速确定方法，该方法适用于一个或多个气象站，不考虑测风塔与气象站的相关性好坏且保证全部风速均为正值。通过实例计算可以看出，本实施例方法结果可靠，求解过程较为简单，具有很好的实用性。
[0141]
实施例2
[0142]
本实施例提供一种风电场代表年风速确定系统，包括：
[0143]
风速获取模块，被配置为获取测风塔的测风年风速和气象站的测风风速，所述气象站的测风风速包括气象站的历史年风速和测风年风速；
[0144]
边缘分布模块，被配置为对测风塔的测风年风速、气象站历史年风速与气象站测风年风速分别进行边缘分布曲线拟合；
[0145]
相关计算模块，被配置为根据测风塔的测风年风速与气象站历史年风速和气象站测风年风速得到秩相关系数；
[0146]
风速分布模块，被配置为根据秩相关系数、测风塔的拟合结果和气象站的拟合结果构建测风年连接函数和历史年连接函数，以分别构建测风塔在测风年与代表年的风速分布图；
[0147]
风速修订模块，被配置为根据测风塔在测风年与代表年的风速差得到风速订正值，以此修订代表年风速。
[0148]
此处需要说明的是，上述模块对应于实施例1中所述的步骤，上述模块与对应的步骤所实现的示例和应用场景相同，但不限于上述实施例1所公开的内容。需要说明的是，上述模块作为系统的一部分可以在诸如一组计算机可执行指令的计算机系统中执行。
[0149]
在更多实施例中，还提供：
[0150]
一种电子设备，包括存储器和处理器以及存储在存储器上并在处理器上运行的计算机指令，所述计算机指令被处理器运行时，完成实施例1中所述的方法。为了简洁，在此不再赘述。
[0151]
应理解，本实施例中，处理器可以是中央处理单元cpu，处理器还可以是其他通用处理器、数字信号处理器dsp、专用集成电路asic，现成可编程门阵列fpga或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件等。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器等。
[0152]
存储器可以包括只读存储器和随机存取存储器，并向处理器提供指令和数据、存储器的一部分还可以包括非易失性随机存储器。例如，存储器还可以存储设备类型的信息。
[0153]
一种计算机可读存储介质，用于存储计算机指令，所述计算机指令被处理器执行时，完成实施例1中所述的方法。
[0154]
实施例1中的方法可以直接体现为硬件处理器执行完成，或者用处理器中的硬件及软件模块组合执行完成。软件模块可以位于随机存储器、闪存、只读存储器、可编程只读存储器或者电可擦写可编程存储器、寄存器等本领域成熟的存储介质中。该存储介质位于存储器，处理器读取存储器中的信息，结合其硬件完成上述方法的步骤。为避免重复，这里
不再详细描述。
[0155]
本领域普通技术人员可以意识到，结合本实施例描述的各示例的单元即算法步骤，能够以电子硬件或者计算机软件和电子硬件的结合来实现。这些功能究竟以硬件还是软件方式来执行，取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能，但是这种实现不应认为超出本技术的范围。
[0156]
上述虽然结合附图对本发明的具体实施方式进行了描述，但并非对本发明保护范围的限制，所属领域技术人员应该明白，在本发明的技术方案的基础上，本领域技术人员不需要付出创造性劳动即可做出的各种修改或变形仍在本发明的保护范围以内。