答卷编号：
论文题目：风电功率预测问题
指导教师：王化东
参赛学校：白城师范学院
报名序号：500
证书邮寄地址：
吉林省白城市中兴西大路57号白城师范学院西校区数学系马超
答卷编号：
风电功率预测问题
摘 要
本文通过建立数据预测模型对风电机发电功率尽可能地进行准确的预测，以减少大幅度风电功率波动对电网的功率平衡和频率调节带来的不利影响。
本文建立了三个数学预测模型，如下: 1, 灰色预测法
在提取部分历史数据的情况下，构建生成列，通过引用(1,1)GM 的灰色微分方程
(0)
(1)
()()x
k a z k b +=和此微分方程的白化方程
(1)
(1)
dx
ax
b dt
+=，并列出(0)
(1)
()()x
k az
k b
+=的时间响应序列，得到简单预测模型(0)
(1)
(1)
(1)(1)()x k x k x k ∧∧∧+=+-；
2.最小二乘法
同样提取部分数据，构建生成列，将数据的变化趋势明显化，即通过构建生成列完成，然后通过最小二乘法，将生成列里面的数据通过m atlab 拟合成函数曲线，并得到曲线图，在此之前通过观察生成列的数据分布图，可初步断定数据基本按直线增长，所以为一次曲线。得到拟合曲线的表达式，即预测模型。
3．指数平滑法
构建生成列，在E xcel 表格中画出生成列曲线和生成列平滑值曲线，经过两次平滑值求解，即通过求一次指数平滑值和二次指数平滑值，得到预测模型为 t t t T y a b T ∧
+=+
1,2,...T =；(0)
(1)
16(16)(161)T x
T y x
T ∧∧
++=-+-
(1)
(2)
2t t
t
a S S =-，(1)
(2)
()1t t
t
b S S α
α
=
--。
关键词：灰色预测法 最小二乘法 指数平滑法 微分方程白化方程拟合曲线一次指数平滑值和二次指数平滑值 预测模型 生成列 响应序列
一、问题重述
大规模风电场接入电网运行时，大幅度地风电功率波动会对电网的功率平衡和频率调节带来不利影响。
根据电力调度部门安排运行方式的不同需求，运用日前预测和实时预测两种方法进行风电功率预测，就能够根据风电功率变化预先安排调度计划，保证电网的功率平衡和运行安全。我们的任务是：
（一）采用不少于三种预测方法构造数学模型，进行风电实时预测及误差分析，对给定
数据进行风电功率实时预测并检验预测结果是否满足附件1中的关于预测精度的相关要求。并分析方法的准确性及推荐最优方法。
（二）试分析风电机组的汇聚对于预测结果误差的影响。根据问题一的预测结果构建数
学模型，总结带有普遍性的规律，进而对风电机组汇聚给风电功率预测误差带来的影响做出预期。
（三）再问题一基础上，构建数学模型提供更高精度的实时预测方法，用预测结果说明
有效性，分析预测精度进一步改善的因素，说明风电功率预测精确度能否无限提高。
二、符号说明
(1)
t
S ：一次指数平滑； (2)
t
S ：二次指数平滑；
t
：当前时点数；
T :为由当前时期t 到预测期的时期数；
t a ：截距； t b ：斜率；
α：加权系数；
(0)
X ：原始数据数列；
(1)
X
：通过累加得到的生成数列；
(0)
()x i ，1,2...,i n =：风电机的单位时点的输出功率；
(1)
()x
i ，1,2...,i n =：单个风电机输出功率累加生成数列的第i 项；
(1)
Z
：(1)X 的紧邻均值()MEAN 生成序列；
(1)
()z
k ：(1)
Z 数列的第k 项；
(1,1)GM :G →Grey (灰色)， M →M odel ，（1→1阶方程，1→1个变量）；
a ：称为发展系数；
b ：灰色作用量；
(1)
(1)x
k ∧+：
通过代入灰色微分方程的时间序列中k 值求得的所建生成数列中的第1k +项预测值；
(1)
()x
k ∧：通过代入灰色微分方程的时间序列中k 值求得的所建生成数列中的第k 项预测
值；
(0)
(1)x
k ∧+：通过代入灰色微分方程的时间序列中k 值求得的所建原始数据数列中的第
1k +项预测值；
三、模型假设
(1)指数平滑法中，首先假设数据的总体变化起伏幅度不大，并按一定的函数曲线走向进行变化，有较明显的规律性；
(2) 最小二乘法中假设数据按一次函数模型变化，并且线性相关性较为明显； (3)灰色预测法
四、问题分析
通过三种预测方法，我们将对风电记得实时，和日前输出功率进行较为精确的预测，我们考虑到了通过历史数据进行预测，预测量的数量较多，也是必要的，模型在经过大量计算后，才可以进行实际方面的应用；风电机组的汇聚对于预测的结果结果误差的影响，我们可通过多组数据对比，进行结论，我们计算过了58p 的预测值，通过实际数据的对
比，分析误差与小数量风电机或单个风电机的误差波动性有明显的指数上升趋势，数量越多，误差的产生率越大，数量越少误差变化越为平稳。
问题三中，要求我们建立更精确的模型进行预测。首先我要考虑到，风电机汇聚的影响系数，然后是数据采集与应用方面的做法，数据采集要预测时点近期，采集数量要适中，由于物理方面的知识可知，发电机输出功率的变化成波动曲线变化，较容易收到外界的电压，电离子流动的影响，所以，第一选好时间点，选近期数据；充分考虑外界环境，与自身产生的自发影响，我们可因此建立观察模型，在通过预测模型对风电功率进行预测。真正做到最为精确的预测，还需要时间的不断增长，发现事物的主潜在变化规律，才能做好预测工作。
五、模型的建立与求解
（1）灰色预测法
下面我们将用灰色预测法建立灰色预测模型。
首先我们利用的是数据库中所给出的风电机在每一时点的输出功率数值，在某一区
段内按时间先后的顺序排列成的数列，一台风电发电机的数据排列成一条数列，数据局限于数据库中的数据；我们用
(0)
(0)
(0)
(0)
((1),(2),...,()),X
x
x
x
n =
来表示功率数据的数列；
接下来我们要通过累加的方法得出生成列
(1)
(1)
(1)
(1)
((1),(2),...,()),X
x
x
x
n =
(1)
(0)
1
()(),1,2,...,k
i x
k x
i k n
==
=∑；
令(1)Z 为(1)X 的紧邻均值()MEAN 生成序列
(1)
(1)
(1)
(1)
((2),(3),...,())Z z
z
z
n =；
(1)
(1)
(1)
()0.5()0.5(1)z
k x k x
k =+-，2,3,...,k n =；
则(1,1)GM 的定义型，即(1,1)GM 的灰微分方程模型为
(0)
(1)
()()x
k az
k b +=；
此灰微分方程的白化方程，也称影子方程为：
(1)
(1)
dx
ax
b dt
+=；
(0)
(1)
()()x
k az
k b +=的时间响应序列为
(1)
(1)
(1)[(0)]ak
b b x k x
e
a
a
∧-+=-
+
，1,2,...,k n =
取(1)(0)(0)(1)x x =，则
(1)
(0)
(1)[(1)]ak
b b x k x
e
a
a
∧-+=-
+
，1,2,...,k n =；
还原值：
(0)
(1)
(1)
(1)(1)()x k x k x
k ∧∧∧+=+-；
上式便为预测方程。
通过代入(0)(1)()()x k az k b +=临近5月31日两组非检验区的(0)1()x k ，(1)1()z k 和
(0)
2()x
k ，(1)
2()z
k ，求得a ，b
的数值；
在计算过程中，规定5月31日0时0分的次序在每一个生成列里为1n +；规定从
5月30号20时0分开始，建立一个累加生成列，实际数据在表中均列出；P：(1)A X（5月30日20时0分到23时45分），
A
P：(1)B X(5月30日20时0分到23时45分)，
B
P：(1)C X(5月30日20时0分到23时45分)，
C
P：(1)D X(5月30日20时0分到23时45分)，
D
P:(1)4X(5月30日20时0分到23时45分);
4
（2）最小二乘法
下面将应用生成数列和最小二乘法并通过m atlab 实现拟合函数的创建。
首先要提取数据库中，各台指定检验风电机的各时点的输出功率，按时间顺序构建数列
(0)
(0)
(0)
(0)
((1),(2),...,()),X
x
x
x
n =
通过累加的方法得出生成数列
(1)
(1)
(1)
(1)
((1),(2),...,()),X
x
x
x
n =
(1)
(0)
1
()(),1,2,...,k
i x
k x
i k n
==
=∑
在在时间区域内，取部分输出功率数据，构建累加生成列：
四台风电发电机数据区域分别为非随机抽取日期非随机时点（我们欲预测5月31日0时0分至23时0分的输出功率，所以我们建立近期生成列来预测下阶段四个小时内的各时点的输出功率）：
A P ：(1)
A
X
（5月30日20时0分到23时45分）， B P ：(1)
B
X (5月30日20时0分到23时45分)，
C
P ：(1)
C X (5月30日20时0分到23时45分)，
D P ：(1)
D
X (5月30日20时0分到23时45分)，
4P :(1)
4
X (5月30日20时0分到23时45分);
通过观察数据表格，发现任意其中一台风电发电机的输出功率在某个阶段，数值较为平稳，也可通过E xcel 对数据进行列表，画出数据分布曲线图，如图（1），来辨别生成数列的走向满足一次曲线还是二次曲线，最后应用最小二乘法，通过m atlab 的实现，求出拟合曲线的表达式； m atlab 程序如下： 见附录1
（3）指数平滑法
下面我们将应用指数平滑法建立预测模型。 首先通过前面对数据的直观分析，数据并不进行一定的特别明显的增减变化趋势，于是，我们创建生成列，以构建有较为明显的增减趋势的数列，从而方便建模预测； 1. 生成列
(1)
(1)
(1)
(1)
((1),(2),...,()),X
x
x
x
n =
(1)
(0)
1
()(),1,2,...,k
i x
k x
i k n
==
=∑；
1.分别列出：
(1)
A
X，(1)B X，(1)C X，(1)D X；
2分别.用E xcel进行统计列表，画出数据分布曲线图；
时点编号
输出功率
一次指数平
滑
二次指数平
滑
1 347.0625
2 650.7188 347.0625
3 1002.563 529.25628 347.0625
4 1408.
5 813.240072 474.598146
5 1711.5 1170.396089 711.6474942
6 1962.375 1495.058496 1032.77151
7 2232.656 1775.448458 1356.3724
8 2438.906 2049.773223 1649.725641
9 2609.531 2283.253129 1929.758949
10 2900.156 2479.020092 2177.204875
11 3167.906 2731.701877 2388.475527
12 3371.156 2993.424591 2628.733972
13 3678.844 3220.063676 2884.017405
14 3957.938 3495.331811 3119.249795
15 4383.094 3772.895344 3382.507206
16 4578.282 4139.014538 3655.778903
以上为A P 的指数平滑曲线，B P ，C P ，D P ，4P ，58P ，我们做同样的处理，省略图表的构建；
3.通过图表可看出数据呈明显线性增长趋势，预测应使用一次和二次指数平滑，意在修正滞后偏差；
设生成列中时间序列为12,,...,...t y y y ，则一次指数平滑公式为：
(1)
(1)
1(1)t
t t S y S αα-=+-；
二次指数平滑公式为：
(2)
(1)
(2)
1
(1)t
t
t S S S αα-=+-；
式中(1)t S 为第t 周期的一次指数平滑值，(2)t S 为二次指数平滑值；α为加权系数，0<α<1.
预测模型为：
t t t T y a b T ∧
+=+ 1,2,...T =；
(1)
(2)
2t t
t
a S S =-，(1)
(2)
()1t t
t
b S S α
α
=
--；
3.通过图标中的平滑值可求得16a ，16b （在图表中α取值为0.3）
(1)
(2)
161616
2a S S =-
(1)
(2)
161616()1b S S α
α
=
--
A P :16a =4622.250173
B P :16a =4002.759027
C P :16a =4423.46649
16b =207.101 16b =192.1004 16b =205.795521
D P :16a =3933.504385
4P ：16a =18493.17587
16b =171.626
16b =841.417
161616T y a b T
∧
+=+
(0)
(1)
16(16)(161)T x
T y x
T ∧∧
++=-+-
求解
图e
六、模型的检验与评价
通过检验，预测值误差均不超过，题所给范围。
参考文献
（1）《点集拓扑学》余玄冰编译北京师范大学出版社出版年限1983
（2）《点集拓扑讲义》徐森林胡白胜金亚东薛春华丛书名高等学校教材点集拓扑学出版日期2007年7月
（3）百度文库《数学预测法》
(4)豆丁网《数学预测法》
（5）《运筹学》主编胡运权清华大学出版社第三版
附录：
附录
x=[1,3,5,7,9,11,13,15];
y=[347.0625, 650.7188, 1002.563 , 1408.5 , 1711.5, 1962.375 , 2232.656 , 2438.906] ;
p=polyfit (x,y,1)
x1=1:2:15;
y1=polyval(p,x1);
plot(x,y,'*r',x1,y1,'-b')
051015
153.1847243.8077
A A y x =+
x=[1,3,5,7,9,11,13,15];
y=[309, 599.625, 914.0625, 1099.0313 , 1283.0626, 1458.5626, 1771.4064 , 1977.9377] ; p=polyfit (x,y ,1) x1=1:2:15;
y1=polyval(p,x1); plot(x,y ,'*r',x1,y1,'-b')
115.2322254.7288B B y x =+
x=[1,3,5,7,9,11,13,15];
y=[407.3438, 858.8438, 1064.1563, 1222.1251 , 1427.8126, 1594.6876, 1958.8126, 2250.1876] ; p=polyfit (x,y ,1) x1=1:2:15;
y1=polyval(p,x1); plot(x,y ,'*r',x1,y1,'-b')
120.2204386.2327
C C y x =+
x=[1,3,5,7,9,11,13,15];
y=[345.9375, 671.9063, 1014.7501, 1152.8439 , 1304.0627, 1489.5957, 1840.502, 2031.6583] ; p=polyfit (x,y ,1) x1=1:2:15;
y1=polyval(p,x1); plot(x,y ,'*r',x1,y1,'-b')
5
10
15
114.3975316.2271D D y x =+
x=[1,3,5,7,9,11,13,15]; y=[1490.344, 2944.094, 4241.532, 5212.5009, 6141.4385,7006.2198, 8391.3758,
9374.6884] ; p=polyfit (x,y ,1) x1=1:2:15;
y1=polyval(p,x1); plot(x,y ,'*r',x1,y1,'-b')