2024华为杯研究生数学建模竞赛A题精品成品论文已出!
A题 风电场有功功率优化分配
一、问题分析
A题是一道工程建模与优化类问题,其目的是根据题目所给的附件数据资料分析风机主轴及塔架疲劳损伤程度,以及建立优化模型求解最优有功功率分配策略。整个问题的模型主要可以分为三大块:计算疲劳损伤程度的实时模型、计算应力/扭矩值的模型、求解最优有功功率分配策略的优化模型。相对应的,也需要参赛团队有对应方面的一些知识积累和技能掌握:第一,须具备扎实的数学建模能力,能够基于风力机的物理特性建立合理的数学模型。对于问题一和二,特别需要对风能转换、力学受力分析、疲劳损伤理论(如S-N曲线、Palmgren-Miner理论)等领域有深入理解。这涉及到风速与功率、扭矩与推力、疲劳损伤的关系建模,并将这些物理关系转换为数学表达式。第二,解决问题三和四需要对多目标优化算法有深入理解,能够设计优化问题并找到最优解。同时,求解大规模实时问题需要考虑算法的计算效率和精度。第三,处理风电场的数据时,必须具备数据分析和信号处理能力,特别是对于噪声和延迟的处理。问题四中的测量噪声和通信延迟的处理需要掌握滤波技术(如卡尔曼滤波)、信号分析方法以及鲁棒控制的理论。
问题一:要求我们利用主轴扭矩和塔架推力荷载数据分别实现对风机的主轴和塔架累积疲劳损伤程度的低复杂度计算。对于题目中给出的数据,首先应使用雨流计数法统计不同幅值载荷相应的循环次数(参考[2]),然后利用附件B中的步骤,计算等效疲劳载荷和累计疲劳损伤值,这里应当对基于雨流计数法计算所得到的累积疲劳损伤值和等效疲劳载荷与给出的数据进行对比,误差应该相对较小。为了实现对雨流统计法的简化,可以考虑[8],[9]中的方法,例如参考文献[8]中提到点的选取遵循单向选取、循环进行的原则这种处理方法使得雨流计数判定只向一个方向进行避免回退取点来进行计数判定的情况的出现降低了程序实现的复杂性。此问求解较为复杂,需要计算时间小于1.00s且求解时间要尽可能地短。也可以考虑利用能量进行计算和改进。
问题二:根据风机所处位置的风速条件和功率参考值,估算当前风机所承受的应力/扭矩。结合受力分析和能量守恒,关键是理解风速、功率与主轴扭矩和塔顶推力之间的物理关系,中间参数应包括风机实际输出功率,桨距角低速轴转速,高速轴转速等。[11]涵盖了风力机的空气动力学分析,包括主轴扭矩和塔顶推力的计算方法。[10]详细介绍了风力发电的基本原理,包括风速、功率系数和风机输出功率的计算。风机的工作过程涉及风能向机械能的转化,其主轴和塔架在此过程中承受不同的力。因此,首先需要通过风速和功率的关系来计算出风机的受力情况,进而推导出主轴扭矩和塔顶推力。风机从风中提取能量,发电功率与风速的三次方成正比。通过风轮扫掠的风能转化为风机的输出功率,能量转换效率由功率系数决定,它反映了风机叶片在不同风速条件下的工作效率。通过这个关系,可以确定风速对风机输出功率的影响。主轴扭矩是通过风机发电功率和风机转速的关系计算的。主轴扭矩与发电功率和转速成反比,而转速则与风速有关。通常,风机的转速由叶尖速比(TSR)控制,它表示叶片旋转速度与风速的比例。基于此,可以通过风速和功率关系,结合叶尖速比计算主轴的扭矩。风机的塔顶推力是风对风轮施加的水平力,它与风速的平方成正比。塔顶推力直接作用于风机的塔架,风轮的受风面积、风速、空气密度以及推力系数都是影响推力的关键因素。推力系数通常根据风机设计确定,并与叶片的角度和气动特性有关。通过这些参数,可以估算出塔顶推力的大小。问题二需要通过模型计算全部时刻的应力/扭矩计算,并展示计算结果与实际数据的对比结果。
问题三:通过优化算法对风电场的有功功率进行实时分配,以降低风机的疲劳损伤并保证功率平衡。该问题可以通过构建多目标优化模型来解决,目标是最小化风机主轴和塔架的累积疲劳损伤,同时满足有功功率总和等于电网调度指令 Pt,并遵守功率参考值的约束。模型的目标函数可以设定为最小化所有风机的总体疲劳损伤程度。单个风机的总疲劳损伤包含主轴疲劳损伤和塔架疲劳损伤两部分. 模型需要满足功率平衡约束(所有风机的有功功率参考值之和应等于电网调度指令), 功率参考值约束:每台风机的功率参考值不得超过风机的额定功率(5MW), 平滑度约束:每台风机的功率分配值与平均分配值的差异不超过1MW. 该问题可通过连续优化方法进行求解,以确保优化过程的实时性。为了加快求解速度,可以使用梯度下降法或拉格朗日乘子法进行快速迭代,确保每秒进行一次功率分配计算。
问题四:在风电场的有功功率优化调度过程中,考虑通信延迟和测量噪声对优化结果的影响。优化模型需要在不确定性条件下保持稳健性,因此在设计时必须增强模型的鲁棒性,以应对噪声和延迟带来的影响。首先是测量噪声和通信延迟的处理,可以通过构建鲁棒优化模型对噪声进行处理,采用预测控制方法(如模型预测控制,MPC[15])来解决延迟问题。同时,噪声在优化中可以通过设置误差容忍度来建模,以确保噪声干扰不会显著影响优化结果[16]。
二、问题一模型建立与求解
2.1 问题一求解思路
风机的主轴和塔架承受着随机波动的载荷,这种载荷的循环特征不易通过简单的方法统计,通常使用雨流计数法来分析载荷的循环特性。然而,标准雨流计数法计算复杂且无法满足实时计算需求。因此,我们提出了基于三点式雨流计数法的低复杂度算法,能够在较短的时间内有效计算累积疲劳损伤。
2.2 问题一模型建立
1. 累积疲劳损伤计算
累积疲劳损伤的计算基于Palmgren-Miner线性累积损伤理论。该理论表明材料在不同应力幅值下的损伤是可累积的,其损伤程度𝐷可表示为:
其中,ni是应力幅值Si下的循环次数,Ni是该应力幅值下的疲劳寿命。累积疲劳损伤D达到1时,材料将失效。
2. S-N曲线模型
S-N曲线用于描述材料在不同应力幅值下的疲劳寿命。Wohler指数m和常数C决定了材料的疲劳 特性,其表达式为:
通过双对数变形,疲劳寿命N可计算为:
在模型中,m取值10,材料常数C=9.77×10^70,通过此公式可以计算材料在特定应力幅值下的疲劳寿命。
3. Goodman修正公式
在实际应用中,由于应力循环的平均应力S_m不为零,累积疲劳损伤需进行修正。Goodman修正公 式用于计算平均应力S_m作用下的等效应力幅值S_e:
其中,S_a是应力幅值,σ_b是材料的最大强度。通过此修正公式,可以得到更准确的疲劳损伤计算结果。
4. 改进的三点式雨流计数法
为了识别载荷数据中的应力循环,我们使用了改进的三点式雨流计数法。该方法通过分析载荷数据中的波峰和波谷,识别应力循环,并统计不同幅值的循环次数。
(1)波峰与波谷识别
假设载荷时间历程为一个离散时间序列X={
x1,x2,...,xn} ’findpeaks 函数的目标是从序列X中识别出局部极大值(波峰)和局部极小值(波谷)。
波峰识别的核心是找到序列中满足以下条件的点x_i:
其中i为波峰所在的位置。波谷可以通过对序列取反来识别
(2)索引排序
识别出的波峰和波谷可以表示为两个集合:
为了保证极值序列的时间顺序,需要将波峰和波谷的索引I_P和I_V合并并进行升序排列,生成极值 序列E={ e_1,e_2,...,e_l}及其对应的索引I_E,使得:
并且极值序列E满足时间上的顺序e_1,e_2,...,e_l.
(3)三点式雨流计数法的规则
在三点式雨流计数法中,每次考察极值序列E中的三个连续点e_i,e_(i+1),e_(i+2)。通过计算它们之间的幅值变化,可以识别闭合应力循环。
幅值计算
定义幅值变化ΔS1和ΔS2为:
根据这两个幅值变化,判断是否识别到闭合循环。
闭合循环识别条件
闭合循环的条件是:
如果满足该条件,则认为e_i和e_(i+1)之间存在一个完整的应力循环,幅值为ΔS1,均值为:
此时,将e_i和e_(i+1)从极值序列中移除,并继续处理新的极值序列。
无闭合循环的处理:
如果不满足闭合循环的条件,即ΔS1>ΔS2,则移动到下一个三点组合e_(i+1),e_(i+2),e_(i+3)进行分析。
2.3 问题一模型求解与分析
为了展示模型的计算结果,选择5-10个有代表性的风机,绘制主轴扭矩和塔架推力的时序数据、等效疲劳载荷及累积疲劳损伤的累积图,并与参考数据进行对比。最后,输出所有100台风机在100秒时长内每秒的累积疲劳损伤值。
-
为了保证模型的实时性,我们对标准雨流计数法进行了简化,同时减少了复杂的迭代步骤。整个计算过程可以在1秒内完成,适合在CPU环境下实时运行。
综上所述,基于三点式雨流计数法、Goodman修正及S-N曲线模型的低复杂度算法能够高效、准确地计算风机主轴及塔架的累积疲劳损伤值,建模方法100台风机的所有数据样本均有效,并满足实时计算需求。
-
2.4 参考代码
-
clc, clear, close all % ================== 读取 .xls 文件 ================== filename = 'data1.xls'; % 文件名 % 读取第2个工作表(主轴扭矩数据) torque_data = xlsread(filename, 2); % 主轴扭矩 time = torque_data(1:end-2, 1); % 第一列为时间 torque_wind_turbines = torque_data(1:end-2, 2:end); % 每台风机的扭矩载荷时序数据 % 读取第3个工作表(塔架推力数据) thrust_data = xlsread(filename, 3); % 塔架推力 thrust_wind_turbines = thrust_data(1:end-2, 2:end); % 每台风机的塔架推力载荷时序数据 % 提取参考的等效疲劳载荷与累积疲劳损伤值(最后两行) ref_equivalent_load_torque = torque_data(end-1, 2:end); % 主轴扭矩的等效疲劳载荷 ref_fatigue_damage_torque = torque_data(end, 2:end); % 主轴扭矩的累积疲劳损伤值 ref_equivalent_load_thrust = thrust_data(end-1, 2:end); % 塔架推力的等效疲劳载荷 ref_fatigue_damage_t