MATLB|基于matpower优化调度的风力模型预测

 💥💥💥💞💞💞欢迎来到本博客❤️❤️❤️💥💥💥

🏆博主优势：🌞🌞🌞博客内容尽量做到思维缜密，逻辑清晰，为了方便读者
📝目前更新：🌟🌟🌟电力系统相关知识，期刊论文，算法，机器学习和人工智能学习。
🚀支持：🎁🎁🎁如果觉得博主的文章还不错或者您用得到的话，可以关注一下博主，如果三连收藏支持就更好啦！这就是给予我最大的支持！

👨‍🎓博主课外兴趣：中西方哲学，送予读者：

👨‍💻做科研，涉及到一个深在的思想系统，需要科研者逻辑缜密，踏实认真，但是不能只是努力，很多时候借力比努力更重要，然后还要有仰望星空的创新点和启发点。当哲学课上老师问你什么是科学，什么是电的时候，不要觉得这些问题搞笑，哲学就是追究终极问题，寻找那些不言自明只有小孩子会问的但是你却回答不出来的问题。建议读者按目录次序逐一浏览，免得骤然跌入幽暗的迷宫找不到来时的路，它不足为你揭示全部问题的答案，但若能让人胸中升起一朵朵疑云，也未尝不会酿成晚霞斑斓的别一番景致，万一它居然给你带来了一场精神世界的苦雨，那就借机洗刷一下原来存放在那儿的“真理”上的尘埃吧。

     或许，雨过云收，神驰的天地更清朗.......🔎🔎🔎

📋📋📋本文目录如下：⛳️⛳️⛳️

目录

1 概述

2 讲解

3 参数表示 

4 运行结果展现

5 Matlab代码实现 

---

1 概述

风光新能源出力受夭气、气候影响巨大,短期、中长期电力供应均处于难预测状态。受风向、风力、光照强度等影响,风光出力具有随机性、波动性、间歇性等特点,新能源日内出力波折;寒潮或连晴高温等极端气候下,电力需求陡增,但新能源最小出力可能处于较低水平，平衡支 撑能力不足。

大规模新能源并网,新能源电力、电量占比不断提高,消纳与利用率矛盾突出:与常规机组相互挤占通道,重要断面满功率运行压力较大;风电反调峰特性突出,加剧电网调峰困难,消纳与利用率相互制约,保证消纳水平的基础上,利用率存在下降的趋势。在大型风电场附近如果有500kV特高压输电线路,这个风电场风电利用率可以基本可以达到接近100% ,可以通过特高压线路输送到远方;如果没有特高压线路,只能就近消纳。还有新能源只能发有功,由于功角因素,有时潮流不是按照我们想要的方向流入,如部分流到其它地区啦,需要在-些地方加装大型无功装置,形成一 个“截面” ,不让新能源流过去。

风力预测配置文件基于计量经济学模型生成，该模型由捕获季节性效应的 OLS 组件和捕获残差中动态的时间序列组件组成。预测所需的不确定性由时间序列模型的白噪声残差的var-cov矩阵实现，使其能够反映不同风场之间的相关性。本文还提供风力实现配置文件，这些配置文件由实际实现数据或计量经济学模型生成。

2 讲解

使用来自 16 个站点 NPCC 模型的站点 2、6 和 15 的模型创建一个 WY 风模型对象：

widx = [2;6;15];
wm = wy_wind_model('wind_model_npcc', widx);

生成 24 个周期范围的MOST 转换概率，每个周期有 4 个风力情景。

tp = wm.transition_probs(24, 4);

加载原始 NPCC 历史风速数据（以 m/s 为单位的原始风速）并将其转换为 ，因为这是此模型所基于的。log(raw_wind_speed+1)

s = load('wind_data_npcc');
log_wind_data = log10(s.wind_data + 1);

从原始 NPCC 历史数据中提取 3 个感兴趣地点的 48 小时风速实现值，从 2005 年 6 月 1 日中午开始（数据从 2004 年 1 月 1 日凌晨 1 点开始）。

pidx1 = wy_wind_date2pidx(24, [2004 1 1 1 0 0], [2005 6 1 12 0 0]);
wsr = wm.realizations(pidx1, 48, log_wind_data);

从模型中为 3 个感兴趣地点生成 24 小时风速实现，从同一小时开始。 

wsr = wm.realizations(pidx1, 24);

使用默认（多涡轮机）功率曲线将风速实现转换为风力发电实现，以装机容量的分数表示。

wpr = wm.speed2power(wsr);

从模型中为 3 个感兴趣地点生成 24 小时风速预报，从同一小时开始，使用 4 个条柱来表示预报分布。 

ws0 = log_wind_data(pidx1-1, widx);
wsf = wm.forecasts(pidx1, ws0, 24, 4);

使用 中定义的第 4 条功率曲线（海上）将风速预测转换为风力发电实现，以装机容量的分数表示。'wind_power_curve_EIC.txt'

s2p = wy_wind_power_curve(4, 'wind_power_curve_EIC.txt');
wm = wy_wind_model('wind_model_npcc', [2;6;15], s2p);
wpf = wm.speed2power(wsf);

3 参数表示 

np_all — 原始风数据中的周期数
nw_all — 原始数据中的风站数量（wind_data_npcc.mat 为 16）
np - 时期数（例如，对于计划范围）
nw - 风站点的数量
npd — 每天的周期数（通常为 24，对于每小时数据）
nb - 用于 MOST 风模型输入的箱数
widx — () 感兴趣的风场指数向量 nw x 1
pidx — 原始历史风数据的标量周期指数
pidx1 - 第一周期的标量周期指数，相对于原始历史风数据的第一周期
dt — 1 x 6 标准 Matlab 日期向量 [yr, mo, day, hr, min, sec]，指定原始历史风数据中的特定时期，yr = 4 位数年份，mo = 月份（1 -12)，day = 日期 (1-31)，hr = 小时 (0-23)，min = 分钟 (0-59)，sec = 秒 (0-59)
wind_data - (np_all x nw_all) 风速矩阵（以 m/s 为单位），对应于时期、地点、特定的 和起始日期/时间 (np_allnw_allnpddt0)
log_wind_data — () 矩阵等于 np_all x nw_alllog10(wind_data + 1)
bins — bin 规范，提供为：
箱数 (), ornb
(1 x nb-1) bin 边界向量（标准差系数），其中假设初始和最终，但不包括-Inf+Inf
模型 - 具有字段的结构：
type — 风速模型的类型
0 = 基于 raw_wind_speed 以 m/s 为单位
1 = 基于日志（raw_wind_speed + 1）
npd - 每天的周期数（每小时模型为 24）
dt0 - 对应于创建模型的第一期数据的 Matlab 日期向量，用于从模型创建的预测和实现是 w.r.t.这个日期pidx1
ar1 — () 单个站点的 AR1 系数向量 nw_all x 1
ols - () 单个站点的 OLS 估计参数矩阵：
[C CY1 SY1 CY2 SY2 CD1 SD1 CD2 SD2nw_all x 9]
var_wnr — () 各个站点的协方差矩阵 nw_all x nw_all
ar1_total — 总风的标量 AR1 系数
ols_total — 总风量 1 x 9 的 OLS 估计参数向量
var_wnr_total — 总风的标量方差
ws——风速量，单位取决于模型的类型

4 运行结果展现

 

 

 

部分代码： 

%% MATLB|基于matpower优化调度的风力模型预测

%% 加载数据
s = load('wind_data_npcc');
wind_data = s.wind_data;
log_wind_data = log10(wind_data + 1);

widx = [1;3;5;7];   %% 风场
np = 24;            %% 24小时
nb = 5;           
pidx0 = wy_wind_date2pidx(24, [2004 1 1 1 0 0], [2004 7 27 8 0 0]); %% 5000
nw = length(widx);  %% 风场景数
nrh = 24;           %% x轴长度

%% 创建风的模型
wm = wy_wind_model('wind_model_npcc', widx);

%% 初始化输出
tp = cell(nrh, 1);
wsr = cell(nrh, 1);
wpr = cell(nrh, 1);
wsf = cell(nrh, 1);
wpf = cell(nrh, 1);

%% 可视化讲解（风模型参数输入）
for t = 1:nrh
    %% 创建转移概率矩阵
    tp{t} = wm.transition_probs(np, nb);

    %% x轴时间的起始期指数
    pidx1 = pidx0 - 1 + t;

    %% 产生已实现的风力发电
    wsr{t} = wm.realizations(pidx1, np);
    wpr{t} = wm.speed2power(wsr{t});

    %% 生成预测风电
    ws0 = log_wind_data(pidx1-1, widx);
    wsf{t} = wm.forecasts(pidx1, ws0, np, nb);
    wpf{t} = wm.speed2power(wsf{t});
end

%% 风电预测可视化
for i = 1:nw
    figure(i)
    plot(wsf{1}(:, :, i));
    xlabel('（时间/小时）')
    ylabel('功率/KW')
end

5 Matlab代码实现 

博客主页：电气辅导帮