wordpress站酷主题,如何做设计师个人网站,个人全屏网站模板,网页制作素材网站推荐2024华为杯研究生数学建模竞赛A题精品成品论文已出#xff01; A题 风电场有功功率优化分配
一、问题分析
A题是一道工程建模与优化类问题#xff0c;其目的是根据题目所给的附件数据资料分析风机主轴及塔架疲劳损伤程度#xff0c;以及建立优化模型求解最优有功功率分配…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材料常数C9.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_(i1),e_(i2)。通过计算它们之间的幅值变化可以识别闭合应力循环。
幅值计算
定义幅值变化ΔS1和ΔS2为 根据这两个幅值变化判断是否识别到闭合循环。
闭合循环识别条件
闭合循环的条件是 如果满足该条件则认为e_i和e_(i1)之间存在一个完整的应力循环幅值为ΔS1,均值为 此时将e_i和e_(i1)从极值序列中移除并继续处理新的极值序列。
无闭合循环的处理
如果不满足闭合循环的条件即ΔS1ΔS2,则移动到下一个三点组合e_(i1),e_(i2),e_(i3)进行分析。
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_thrust thrust_data(end, 2:end); % 塔架推力的累积疲劳损伤值% 材料参数
sigma_b 50000000; % 材料拉伸断裂时的最大载荷值
m 10; % Wohler指数
C 9.77e70; % S-N曲线常数
design_life 42565440.4361; % 风机设计寿命,
NPeaks 100;% 计算每台风机的指标
num_wind_turbines size(torque_wind_turbines, 2); % 风机数量calculated_equivalent_load_torque zeros(1, num_wind_turbines);
calculated_fatigue_damage_torque zeros(1, num_wind_turbines);
calculated_equivalent_load_thrust zeros(1, num_wind_turbines);
calculated_fatigue_damage_thrust zeros(1, num_wind_turbines);cpu_times zeros(1, num_wind_turbines); % 记录每台风机的计算时间for i 1:num_wind_turbines% 开始计时tic;% 提取每台风机的主轴扭矩数据torque_data_turbine torque_wind_turbines(:, i);% 提取每台风机的塔架推力数据thrust_data_turbine thrust_wind_turbines(:, i);% 计算主轴扭矩的等效疲劳载荷与累积疲劳损伤值[ranges_torque, means_torque, counts_torque] three_point_rainflow(torque_data_turbine, NPeaks); % 三点式雨流计数法L_i_torque arrayfun((r, m) goodman_correction(r, m, sigma_b), ranges_torque, means_torque);[calculated_equivalent_load_torque(i), equivalent_load_torque_vector] calculate_equivalent_load(L_i_torque, counts_torque, m, design_life);[calculated_fatigue_damage_torque(i), fatigue_damage_torque_vector] calculate_fatigue_damage(L_i_torque, counts_torque, C, m);% 计算塔架推力的等效疲劳载荷与累积疲劳损伤值[ranges_thrust, means_thrust, counts_thrust] three_point_rainflow(thrust_data_turbine, NPeaks); % 三点式雨流计数法L_i_thrust arrayfun((r, m) goodman_correction(r, m, sigma_b), ranges_thrust, means_thrust);[calculated_equivalent_load_thrust(i), equivalent_load_thrust_vector] calculate_equivalent_load(L_i_thrust, counts_thrust, m, design_life);[calculated_fatigue_damage_thrust(i), fatigue_damage_thrust_vector] calculate_fatigue_damage(L_i_thrust, counts_thrust, C, m);% 停止计时记录时间cpu_times(i) toc;% 输出中间变量图像if mod(i, 25) 0 % 每隔25个风机选取一个figure;subplot(4, 1, 1);plot(time, torque_data_turbine);title(sprintf(风机 %d 的主轴扭矩, i));xlabel(时间 (s));ylabel(主轴扭矩 (Nm));subplot(4, 1, 2);plot(time, thrust_data_turbine);title(sprintf(风机 %d 的塔架推力, i));xlabel(时间 (s));ylabel(塔架推力 (N));% 分开绘制等效疲劳载荷和累积疲劳损伤值的累计图subplot(4, 1, 3);plot(cumsum(equivalent_load_torque_vector), b, DisplayName, 等效疲劳载荷累计值);title(sprintf(风机 %d 的等效疲劳载荷累计, i));xlabel(索引);ylabel(累积等效疲劳载荷);legend(show);grid on;subplot(4, 1, 4);plot(cumsum(fatigue_damage_torque_vector), r, DisplayName, 累积疲劳损伤值);title(sprintf(风机 %d 的累积疲劳损伤值, i));xlabel(索引);ylabel(累积疲劳损伤);legend(show);grid on;end
end% 结果对比与图形化
% 绘制每台风机的CPU时间消耗
figure;
plot(1:num_wind_turbines, cpu_times, b-o, LineWidth, 2);
title(每台风机计算的CPU时间);
xlabel(风机编号);
ylabel(CPU时间 (秒));
grid on;% 主轴扭矩的等效疲劳载荷与参考值对比
figure;
subplot(2,1,1);
plot(1:num_wind_turbines, calculated_equivalent_load_torque, b, DisplayName, 计算值);
hold on;
plot(1:num_wind_turbines, ref_equivalent_load_torque, r--, DisplayName, 参考值);
title(主轴扭矩 - 等效疲劳载荷);
xlabel(风机编号);
ylabel(等效疲劳载荷);
legend(show);
grid on;% 主轴扭矩的累积疲劳损伤值与参考值对比
subplot(2,1,2);
plot(1:num_wind_turbines, calculated_fatigue_damage_torque, b, DisplayName, 计算值);
hold on;
plot(1:num_wind_turbines, ref_fatigue_damage_torque, r--, DisplayName, 参考值);
title(主轴扭矩 - 累积疲劳损伤值);
xlabel(风机编号);
ylabel(累积疲劳损伤值);
legend(show);
grid on;% 塔架推力的等效疲劳载荷与参考值对比
figure;
subplot(2,1,1);
plot(1:num_wind_turbines, calculated_equivalent_load_thrust, b, DisplayName, 计算值);
hold on;
plot(1:num_wind_turbines, ref_equivalent_load_thrust, r--, DisplayName, 参考值);
title(塔架推力 - 等效疲劳载荷);
xlabel(风机编号);
ylabel(等效疲劳载荷);
legend(show);
grid on;% 塔架推力的累积疲劳损伤值与参考值对比
subplot(2,1,2);
plot(1:num_wind_turbines, calculated_fatigue_damage_thrust, b, DisplayName, 计算值);
hold on;
plot(1:num_wind_turbines, ref_fatigue_damage_thrust, r--, DisplayName, 参考值);
title(塔架推力 - 累积疲劳损伤值);
xlabel(风机编号);
ylabel(累积疲劳损伤值);
legend(show);
grid on;% 计算平均相对误差
% 主轴扭矩平均相对误差
mean_relative_error_equivalent_load_torque mean(abs((calculated_equivalent_load_torque - ref_equivalent_load_torque) ./ ref_equivalent_load_torque));
mean_relative_error_fatigue_damage_torque mean(abs((calculated_fatigue_damage_torque - ref_fatigue_damage_torque) ./ ref_fatigue_damage_torque));% 塔架推力平均相对误差
mean_relative_error_equivalent_load_thrust mean(abs((calculated_equivalent_load_thrust - ref_equivalent_load_thrust) ./ ref_equivalent_load_thrust));
mean_relative_error_fatigue_damage_thrust mean(abs((calculated_fatigue_damage_thrust - ref_fatigue_damage_thrust) ./ ref_fatigue_damage_thrust));% 输出平均相对误差
fprintf(主轴扭矩的等效疲劳载荷平均相对误差: %.2f%%\n, mean_relative_error_equivalent_load_torque * 100);
fprintf(主轴扭矩的累积疲劳损伤值平均相对误差: %.2f%%\n, mean_relative_error_fatigue_damage_torque * 100);
fprintf(塔架推力的等效疲劳载荷平均相对误差: %.2f%%\n, mean_relative_error_equivalent_load_thrust * 100);
fprintf(塔架推力的累积疲劳损伤值平均相对误差: %.2f%%\n, mean_relative_error_fatigue_damage_thrust * 100);function [ranges, means, counts] three_point_rainflow(load_data, NPeaks)% 输入 load_data 是原始的载荷-时间历程数据% NPeaks 是波峰的最大数量% 输出 ranges 是识别出的循环幅值% means 是识别出的循环均值% counts 是每个循环的次数等于循环起点与终点的索引差% 使用MATLAB内置的findpeaks函数识别波峰和波谷[peaks, peak_indices] findpeaks(load_data, NPeaks, NPeaks); % 限制波峰识别数量[valleys, valley_indices] findpeaks(-load_data, NPeaks, NPeaks); % 限制波谷识别数量取负valleys -valleys; % 取反恢复波谷的原始值% 合并波峰和波谷extremes [peaks; valleys]; % 将波峰和波谷的值组合indices [peak_indices; valley_indices]; % 将波峰和波谷的索引组合% 按索引排序确保极值按时间顺序排列[indices, sort_order] sort(indices);extremes extremes(sort_order); % 按相应索引顺序排序极值% 初始化变量ranges [];means [];counts [];% 应用三点式雨流计数法识别循环while length(extremes) 3i 1;loop_detected false; % 标志是否识别出循环while i length(extremes) - 2% 计算 ΔS1 和 ΔS2deltaS1 abs(extremes(i1) - extremes(i));deltaS2 abs(extremes(i2) - extremes(i1));if deltaS1 deltaS2% 识别出一个闭合循环记录循环幅值和均值S_ai deltaS1; % 循环幅值S_mi (extremes(i) extremes(i1)) / 2; % 循环均值% 存储幅值、均值和计数起点和终点的索引差ranges(end1) S_ai;means(end1) S_mi;counts(end1) abs(indices(i1) - indices(i))2; % 记录索引差% 移除闭合循环中的两个点i和i1形成新的序列extremes(i:i1) [];indices(i:i1) [];% 标志设置为识别出循环loop_detected true;% 重置 i 为1继续从头处理剩余数据i 1;else% 如果没有闭合循环移动到下一个点i i 1;endendif ~loop_detectedtotal_diff 0;total_mean 0;total_count 0;for j 1:(length(extremes) - 1)S_ai abs(extremes(j1) - extremes(j)); % 计算幅值S_mi (extremes(j1) extremes(j)) / 2; % 计算均值count abs(indices(j1) - indices(j)); % 记录索引差% 累加幅值、均值和计数total_diff total_diff S_ai;total_mean total_mean S_mi;total_count total_count count;end% 将剩余部分作为一个完整的循环处理取均值ranges(end1) total_diff / (length(extremes) - 1); % 平均幅值means(end1) total_mean / (length(extremes) - 1); % 平均均值counts(end1) total_count; % 计数为索引差的总和break; % 结束循环endend
end% Goodman修正公式
function L_i goodman_correction(S_ai, S_mi, sigma_b)L_i S_ai / (1 - S_mi / sigma_b); % 根据Goodman修正公式计算L_i
end% 计算等效疲劳载荷
function [L_N, L_i_vector] calculate_equivalent_load(L_i, counts, m, design_life)L_i_vector (L_i.^m) .* counts; % 求和前的向量L_N (sum(L_i_vector) / design_life)^(1/m); % 根据公式计算等效疲劳载荷
end% 计算累积疲劳损伤
function [D, D_vector] calculate_fatigue_damage(L_i, counts_torque, C, m)N_F round(exp(-m * log(L_i) log(C))); % 根据S-N曲线公式计算疲劳寿命 N_FD_vector 1 ./ N_F; % 求和前的向量D sum(D_vector); % 累积损伤值
end 三、问题二模型建立与求解 3.1 问题二求解思路 随着全球对可再生能源需求的不断增长风力发电已经成为满足能源需求的主要可再生能源之一。风力发电机在运行时受风力、气动载荷以及结构载荷的影响长时间运转可能导致设备的疲劳损伤和性能下降。为了准确评估风机的运行状态、延长其使用寿命必须对风机的关键力学参数进行实时监测和估算。塔架推力和主轴扭矩是两个重要的参数直接关系到风机的结构稳定性和功率输出性能。 通常情况下风速与风机的发电功率之间存在正相关关系。当风机在理想条件下运行时风能被风机叶片最大化地捕获并转化为电能。此时塔架推力和主轴扭矩达到最小值。然而当风速过高或发电功率不足时多余的风能作用于风机系统会导致塔架推力和主轴扭矩的增加进而加剧风机的疲劳损伤。因此准确估算风机在不同风速条件下的塔架推力和主轴扭矩对于风机的设计、运行管理和寿命评估具有重要意义。 本章旨在建立一个基于风速和功率参考值的数学模型来估算风机的塔架推力和主轴扭矩。我们结合空气动力学理论和能量守恒原理推导了相应的力学公式并通过100台风机的监测数据进行验证。本文将详细介绍模型的构建、求解过程以及与实际参考值的误差分析。 3.2 受力分析 风力发电机受力分析涉及风力对叶片的推力、气动扭矩和主轴的扭矩传递。通过调整叶片迎角、叶片半径以及控制转速可以优化风力发电机的输出功率并保持系统的稳定性。 1主轴扭矩分析
在风力发电机运行过程中叶片旋转带动主轴产生扭矩。这一扭矩由风速和叶片的空气动力学特性决定。图片1中展示了主轴旋转的方向主轴上的扭矩是由叶片旋转带来的空气动力学作用力产生的直接作用于主轴并被传递给发电机。
2塔架推力分析
塔架推力是风力对叶片产生的作用力这一作用力通过叶片传递到整个塔架。塔架推力在图片1中展示为叶片前方的推力方向它是沿着风向作用的。塔架必须承受风产生的推力并保证整个结构的稳定性。
3迎角与气动力矩分析
图片1中还标注了叶片的迎角和气动扭矩。迎角是指叶片相对于来流风的角度它直接影响气动性能和产生的力矩。通过控制叶片的迎角可以优化气动效率从而调整主轴扭矩和发电功率。
叶片旋转时不仅在旋转平面内产生推力还会因空气动力效应产生气动扭矩。气动扭矩决定了叶片如何受到风力影响并最终通过传动系统转化为电能。
4叶片几何参数的影响
叶片半径是另一个重要参数它决定了风轮扫掠面积的大小。图片1中标明了叶片半径较大的叶片半径意味着更大的扫掠面积和更高的捕风能力从而产生更大的推力和扭矩。
3.3 模型求解与分析
1风能与发电功率的关系
风力发电机的工作原理是将风能转化为机械能然后通过发电机将机械能转化为电能。因此风机的功率输出直接与风速及其所捕获的风能相关。风的动能公式为 其中ρ是空气密度,通常取1.225lg/m^3A是风轮的扫掠面积等于Aπr^2其中r为风机叶片的半径V是风速。风轮从风中捕获的能量与其功率系数C_p有关风机提取的有效功率可表示为 其中C_p是风机的功率系数取决于风机的设计和运行条件。功率系数的理论最大值为 0.593 (Betz定律)但实际运行中功率系数通常在 0.35 到 0.5 之间。 发电机的输出功率P_ref与风轮提取的功率P_w通过发电机效率η相关联。发电功率参考值P_ref的表达式为 在本模型中已知风速V和功率参考值P_ref,我们通过风轮的功率系数C_p以及发电机效率η来反 推出塔架推力和主轴扭矩。
2主轴扭矩的估算
主轴扭矩是风机从风中捕获的机械能通过叶片传递至主轴的扭矩。主轴扭矩可以通过风机的发电功率与转子转速之间的关系推导出来。风轮上的机械功率P_mech与主轴扭矩T_s和角速度ω_r的关系为 其中
⋅Pmech是主轴上的机械功率
⋅T_s是主轴扭矩
⋅ ω_r是风轮的角速度。
因此主轴扭矩可以通过发电功率P_ref和角速度ω_r计算得到 角速度ω_r可通过风速V和风机叶片的半径r以及尖速比λ计算得到具体关系为 其中λ为风轮的尖速比是风机叶尖速度与风速之比通常取值在6到8之间。因此主轴扭矩的 最终计算公式为 该公式表明主轴扭矩与功率参考值、叶片半径和风速成正比。
3塔架推力的估算
塔架推力是指风机叶片在旋转过程中风推动叶片产生的水平推力。塔架推力的大小决定了风机塔架所承受的载荷对风机结构的稳定性具有重要影响。塔架推力可以通过气动力学推导得出其公式为 其中·C_t是推力系数通常取决于风机的设计和运行状况取值范围一般为 0.7 到 0.9;⋅ρ为空气密度Aπr^2是风轮的扫掠面积·V为风速。
推力系数C_t与功率系数C_p之间存在一定的关联一般来说当风机功率系数较大时推力系数也较大。因此通过风速和功率参考值可以估算出塔架推力的大小。
4误差分析
为了评估模型的准确性本文将模型估算得到的主轴扭矩和塔架推力与风机实际参考值进行对比。通过计算估算值与参考值之间的误差平方和来衡量模型的精度。误差平方和的计算公式为 其中
Tsi和Fti分别为第i台风机的估算主轴扭矩和塔架推力
·Trefi和Frefi分别为实际参考主轴扭矩和塔架推力。 通过分析误差平方和我们可以确定模型在不同风速和功率工况下的表现并找出可能导致误差的原因。
3.4 模型求解分析
在该问题中优化的目标是找到能够最小化主轴扭矩误差和塔架推力误差的最优参数组合。对于每一个风机利用双目标优化算法计算最优参数
1. 被优化的参数
功率系数C_p:表示风能转化为机械能的效率其范围通常为[0.3,0.5]。
空气密度ρ: 空气密度会因风场位置和环境条件而变化影响推力和扭矩优化范围为[1.1,1.3]kg/m。
叶片半径R: 叶片的长短直接影响风轮的扫风面积进而影响推力和扭矩的大小优化范围为[60,65]米。
2. 优化方法
为了更好地优化风机的功率系数、叶片半径和空气密度三个可调节参数我们采用了三种多目标优化方法来寻找这些参数的最优解。这三种方法分别是 NSGA-II非支配排序遗传算法、加权和法以及 Pareto 前沿构造法。以下对每种优化方法的特点和应用进行详细描述。
1NSGA-II非支配排序遗传算法
NSGA-II 是一种经典的基于遗传算法的多目标优化方法能够在不损失多个目标函数的情况下生成一组 Pareto 最优解集。在风机系统的优化中我们同时最小化两个目标函数主轴扭矩误差和塔架推力误差。
NSGA-II 的工作流程
初始化种群随机生成若干个初始解种群中的每个个体代表一组参数组合包括功率系数、叶片半径和空气密度。
目标函数评估对于种群中的每个个体计算两个目标函数值即主轴扭矩误差和塔架推力误差。误差通过估算值与参考值的平方和计算。
选择、交叉和变异操作通过基于适应度的选择策略选择出优质个体进行交叉和变异操作生成下一代解。
非支配排序和拥挤距离计算对新生成的个体进行非支配排序找出表现优异的解同时通过拥挤距离对解的分布进行调整确保解在 Pareto 前沿上的均匀分布。
迭代重复选择、交叉和变异操作逐渐逼近 Pareto 最优解集直到达到设定的终止条件如计算时间达到上限。
NSGA-II 的优势在于它能够同时优化多个目标并提供一个 Pareto 最优解集供决策者根据需要在解集中选择最优解。
2加权和法 (Weighted Sum Method)
加权和法是一种将多个目标函数合并为一个目标函数的多目标优化方法。该方法通过引入权重参数 w_1和w_2,将主轴扭矩误差和塔架推力误差的影响合并成一个单一目标函数。优化问题转化为 其中
· ETs表示主轴扭矩误差的平方和。·EFt表示塔架推力误差的平方和。· w1和w2是相应的权重用于调整两个目标之间的相对重要性。
加权和法的工作流程
设置权重首先设定权重w1和w2的初始值。权重的调整可以反映不同优化目标的优先级。例如若更重视最小化主轴扭矩误差则增大w1的值。
单目标优化将合并后的目标函数作为单一目标使用如 fminsearch 这样的单目标优化算法进行求解。通过迭代找到能够最小化加权和误差的参数组合。
调整权重并重复若优化结果不符合预期则可调整权重 w1和w2重复上述过程直到得到最优解。
加权和法的优点是实现了从多目标到单一目标的转化简化了优化过程。然而该方法对权重的选择较为敏感若权重设置不当可能无法找到所有的 Pareto 最优解。
3Pareto 前沿构造法Pareto-front Method
Pareto 前沿构造法是一种直接寻找 Pareto 最优解集的方法。在多目标优化问题中Pareto 最优解集代表了所有在不同目标之间无法互相超越的解。这意味着任何一个 Pareto 最优解在某一目标上变得更优的同时会导致另一个目标变得更差。
Pareto 前沿构造法的工作流程
初始化参数范围设定待优化参数的范围包括 CpC_pCp、叶片半径 RRR 和空气密度 ρ\rhoρ。
目标函数评估对于每一组参数组合计算两个目标函数主轴扭矩误差和塔架推力误差。
Pareto 排序对所有候选解进行 Pareto 排序找出那些在两个目标上均无法被其他解支配的解集。
构造 Pareto 前沿从所有候选解中构造 Pareto 前沿生成一组 Pareto 最优解。这些解在不同目标函数之间表现出合理的权衡关系供决策者选择。
优势Pareto 前沿构造法能够直接生成一组 Pareto 最优解允许在不同目标之间进行灵活的权衡选择。这在多目标优化问题中非常有用。
3. 优化结果
下图给出了针对100个风机三种方法计算得到全部时刻估算值与参照值之差的平方和统计图。表格给出了针对每种算法100个风机平方和的总和表格根据表格中的全部时刻估算值与参照值之差的平方和加权和法表现较好因此将加权和法获得的全部时刻估算值填入附件6中。 加权和法 NSGA-II Pareto-front 所有风机主轴扭矩总误差平方和 2.3317e22 2.3317e22 2.3317e22 所有风机塔架推力总误差平方和 1.6234e15 1.0655e16 1.8554e16
3.5 结果分析
通过对100台风机的数据分析可以看出本文建立的数学模型能够有效估算风机的主轴扭矩和塔架推力。模型的误差平方和较小特别是在额定风速附近估算结果与实际参考值非常接近。这表明基于风速和功率参考值的模型在大多数工况下是可靠的。
然而在风速超过额定值时风机的控制策略通常会限制叶片的转速以保护设备这可能导致主轴扭矩的变化不再遵循理论计算结果。因此未来的模型改进可以考虑引入更多风机的控制策略参数以提高模型的适应性。
此外风机运行时受到的载荷不仅与风速和功率有关还与风切变、湍流强度等气象因素相关。这些因素可能会影响塔架推力和主轴扭矩的实际值。在未来的研究中可以通过结合更多的气象数据来改进模型从而进一步提高估算精度。 3.6 参考代码
clc,clear,close all;
% MATLAB多转子风机控制代码
% 调用三个优化函数
weighted_sum_optimization();
nsga2_optimization();
pareto_front_optimization();% NSGA-II 优化
function nsga2_optimization()% 加载数据Pref readmatrix(sheet1.csv);V readmatrix(sheet2.csv);Tshaft_ref readmatrix(sheet3.csv);Ft_ref readmatrix(sheet4.csv);omega_r readmatrix(sheet6.csv);omega_f readmatrix(sheet7.csv);% 固定参数generator_efficiency 0.9; % 发电机效率gear_ratio 100; % 齿轮箱速比% 参数范围定义只优化Cp、blade_radius和air_densitylb [0.3, 60, 1.1]; % 下界Cp, blade_radius, air_densityub [0.5, 65, 1.3]; % 上界% 设定时间限制为5分钟options optimoptions(gamultiobj, MaxTime, 300, PopulationSize, 100, Display, iter);torque_errors zeros(100, 1); % 用于存储每个风机的主轴扭矩总误差平方和thrust_errors zeros(100, 1); % 用于存储每个风机的塔架推力总误差平方和for turbine 1:100% 目标函数使用固定的generator_efficiency和gear_ratioobjective_func_nsga2 (x) calculate_torque_and_thrust_fixed(x, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, generator_efficiency, gear_ratio);% 使用 NSGA-II 优化[optimal_params, ~] gamultiobj(objective_func_nsga2, 3, [], [], [], [], lb, ub, options);% 计算并保存结果[torque_errors(turbine), thrust_errors(turbine)] save_optimization_results(optimal_params, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, generator_efficiency, gear_ratio);end% 绘制100个风机主轴扭矩和塔架推力的总误差平方和figure;subplot(2, 1, 1);bar(torque_errors);title(NSGA-II: 主轴扭矩总误差平方和);xlabel(风机编号);ylabel(总误差平方和);subplot(2, 1, 2);bar(thrust_errors);title(NSGA-II: 塔架推力总误差平方和);xlabel(风机编号);ylabel(总误差平方和);% 计算所有风机的误差平方和总和total_torque_error sum(torque_errors);total_thrust_error sum(thrust_errors);fprintf(NSGA-II: 所有风机主轴扭矩总误差平方和: %.4e\n, total_torque_error);fprintf(NSGA-II: 所有风机塔架推力总误差平方和: %.4e\n, total_thrust_error);
end% 加权和法优化
function weighted_sum_optimization()% 加载数据Pref readmatrix(sheet1.csv);V readmatrix(sheet2.csv);Tshaft_ref readmatrix(sheet3.csv);Ft_ref readmatrix(sheet4.csv);omega_r readmatrix(sheet6.csv);omega_f readmatrix(sheet7.csv);% 固定参数generator_efficiency 0.9; % 发电机效率gear_ratio 100; % 齿轮箱速比% 权重w1 0.5;w2 0.5;torque_errors zeros(100, 1); % 用于存储每个风机的主轴扭矩总误差平方和thrust_errors zeros(100, 1); % 用于存储每个风机的塔架推力总误差平方和for turbine 1:100% 目标函数只优化Cp, blade_radius, air_densityobjective_func_weighted_sum (x) calculate_weighted_error_fixed(x, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, generator_efficiency, gear_ratio, w1, w2);% 初始参数init_params [0.4, 63, 1.225]; % 初始值: Cp, blade_radius, air_density% 设置优化选项options optimset(MaxIter, 300, Display, iter);% 使用 fminsearch 优化optimal_params fminsearch(objective_func_weighted_sum, init_params, options);% 计算并保存结果[torque_errors(turbine), thrust_errors(turbine)] save_optimization_results(optimal_params, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, generator_efficiency, gear_ratio);end% 绘制100个风机主轴扭矩和塔架推力的总误差平方和figure;subplot(2, 1, 1);bar(torque_errors);title(加权和法: 主轴扭矩总误差平方和);xlabel(风机编号);ylabel(总误差平方和);subplot(2, 1, 2);bar(thrust_errors);title(加权和法: 塔架推力总误差平方和);xlabel(风机编号);ylabel(总误差平方和);% 计算所有风机的误差平方和总和total_torque_error sum(torque_errors);total_thrust_error sum(thrust_errors);fprintf(加权和法: 所有风机主轴扭矩总误差平方和: %.4e\n, total_torque_error);fprintf(加权和法: 所有风机塔架推力总误差平方和: %.4e\n, total_thrust_error);
end% Pareto-front 优化
function pareto_front_optimization()% 加载数据Pref readmatrix(sheet1.csv);V readmatrix(sheet2.csv);Tshaft_ref readmatrix(sheet3.csv);Ft_ref readmatrix(sheet4.csv);omega_r readmatrix(sheet6.csv);omega_f readmatrix(sheet7.csv);% 固定参数generator_efficiency 0.9; % 发电机效率gear_ratio 100; % 齿轮箱速比% 参数范围定义只优化Cp, blade_radius和air_densitylb [0.3, 60, 1.1]; % 下界ub [0.5, 65, 1.3]; % 上界% 设置时间限制为5分钟options optimoptions(paretosearch, MaxTime, 300, Display, iter);torque_errors zeros(100, 1); % 用于存储每个风机的主轴扭矩总误差平方和thrust_errors zeros(100, 1); % 用于存储每个风机的塔架推力总误差平方和for turbine 1:100% 目标函数objective_func_pareto (x) calculate_torque_and_thrust_fixed(x, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, generator_efficiency, gear_ratio);% 使用 Pareto-front 优化[optimal_params, ~] paretosearch(objective_func_pareto, 3, [], [], [], [], lb, ub, [], options);% 计算并保存结果[torque_errors(turbine), thrust_errors(turbine)] save_optimization_results(optimal_params, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, generator_efficiency, gear_ratio);end% 绘制100个风机主轴扭矩和塔架推力的总误差平方和figure;subplot(2, 1, 1);bar(torque_errors);title(Pareto-front: 主轴扭矩总误差平方和);xlabel(风机编号);ylabel(总误差平方和);subplot(2, 1, 2);bar(thrust_errors);title(Pareto-front: 塔架推力总误差平方和);xlabel(风机编号);ylabel(总误差平方和);% 计算所有风机的误差平方和总和total_torque_error sum(torque_errors);total_thrust_error sum(thrust_errors);fprintf(Pareto-front: 所有风机主轴扭矩总误差平方和: %.4e\n, total_torque_error);fprintf(Pareto-front: 所有风机塔架推力总误差平方和: %.4e\n, total_thrust_error);
end% 计算与保存结果的函数
function [torque_error_sum_of_squares, thrust_error_sum_of_squares] save_optimization_results(optimal_params, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, generator_efficiency, gear_ratio)% 计算结果并与参考值对比num_time_steps size(Pref, 1);Cp optimal_params(1);blade_radius optimal_params(2);air_density optimal_params(3);Tshaft_est zeros(num_time_steps, 1);Ft_est zeros(num_time_steps, 1);for t 1:num_time_stepsP_ref Pref(t, turbine);wind_speed V(t, turbine);omega_r_turbine omega_r(t, turbine);omega_f_turbine omega_f(t, turbine);% 主轴扭矩计算Tshaft_est(t) (P_ref / (generator_efficiency * omega_r_turbine)) * gear_ratio;Tshaft_est(t) min(max(Tshaft_est(t), 0), 10^9);% 塔架推力计算swept_area pi * blade_radius^2;Ft_est(t) 0.5 * air_density * swept_area * Cp * (wind_speed^2);Ft_est(t) min(max(Ft_est(t), 0), 10^7);end% 计算总误差平方和torque_error_sum_of_squares sum((Tshaft_est - Tshaft_ref(:, turbine)).^2);thrust_error_sum_of_squares sum((Ft_est - Ft_ref(:, turbine)).^2);% 保存数据到附件6的表格中% result_data [Tshaft_est, Ft_est];% writematrix(result_data, 附件6.xlsx, Sheet, sprintf(WindTurbine_%d, turbine));
end% 主轴扭矩和塔架推力估算的函数
function [torque_error, thrust_error] calculate_torque_and_thrust_fixed(params, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, generator_efficiency, gear_ratio)% 提取参数Cp params(1); % 功率系数blade_radius params(2); % 叶片半径air_density params(3); % 空气密度% 初始化估算值num_time_steps size(Pref, 1);Tshaft_est zeros(num_time_steps, 1);Ft_est zeros(num_time_steps, 1);for t 1:num_time_stepsP_ref Pref(t, turbine);wind_speed V(t, turbine);omega_r_turbine omega_r(t, turbine);omega_f_turbine omega_f(t, turbine);% 主轴扭矩计算Tshaft_est(t) (P_ref / (generator_efficiency * omega_r_turbine)) * gear_ratio;Tshaft_est(t) min(max(Tshaft_est(t), 0), 10^9); % 限制估算值% 塔架推力计算swept_area pi * blade_radius^2;Ft_est(t) 0.5 * air_density * swept_area * Cp * (wind_speed^2);Ft_est(t) min(max(Ft_est(t), 0), 10^7); % 限制推力估算值end% 计算主轴扭矩和塔架推力的误差平方和torque_error sum((Tshaft_est - Tshaft_ref(:, turbine)).^2);thrust_error sum((Ft_est - Ft_ref(:, turbine)).^2);
endfunction weighted_error calculate_weighted_error_fixed(params, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, w1, w2, generator_efficiency, gear_ratio)% 提取可调参数Cp params(1); % 功率系数blade_radius params(2); % 风机叶片半径air_density params(3); % 空气密度% 初始化估算值num_time_steps size(Pref, 1);Tshaft_est zeros(num_time_steps, 1); % 估算的主轴扭矩Ft_est zeros(num_time_steps, 1); % 估算的塔架推力% 循环计算每个时间步的主轴扭矩和塔架推力for t 1:num_time_stepsP_ref Pref(t, turbine); % 功率调度指令wind_speed V(t, turbine); % 输入风速omega_r_turbine omega_r(t, turbine); % 低速轴转速omega_f_turbine omega_f(t, turbine); % 高速轴转速% 主轴扭矩估算 if wind_speed 0% 基于功率、发电机效率和齿轮箱速比估算主轴扭矩Tshaft_est(t) (P_ref / (generator_efficiency * omega_r_turbine)) * gear_ratio;Tshaft_est(t) min(max(Tshaft_est(t), 0), 10^9); % 限制扭矩估算值在合理范围内elseTshaft_est(t) 0;end% 塔架推力估算 swept_area pi * blade_radius^2; % 计算叶片扫过的面积Ft_est(t) 0.5 * air_density * swept_area * Cp * (wind_speed^2); % 基于风速计算推力Ft_est(t) min(max(Ft_est(t), 0), 10^7); % 限制推力估算值在合理范围内end% 计算主轴扭矩和塔架推力的误差平方和torque_error sum((Tshaft_est - Tshaft_ref(:, turbine)).^2);thrust_error sum((Ft_est - Ft_ref(:, turbine)).^2);% 加权总误差加权和法weighted_error w1 * torque_error w2 * thrust_error;
endclc,clear,close all;
% MATLAB多转子风机控制代码
% 调用三个优化函数
weighted_sum_optimization();% 加权和法优化
function weighted_sum_optimization()% 加载数据Pref readmatrix(sheet1.csv);V readmatrix(sheet2.csv);Tshaft_ref readmatrix(sheet3.csv);Ft_ref readmatrix(sheet4.csv);omega_r readmatrix(sheet6.csv);omega_f readmatrix(sheet7.csv);% 权重w1 0.5;w2 0.5;torque_errors zeros(100, 1); % 用于存储每个风机的主轴扭矩总误差平方和thrust_errors zeros(100, 1); % 用于存储每个风机的塔架推力总误差平方和result1zeros(2000,100);result2zeros(2000,100);for turbine 1:100% 目标函数objective_func_weighted_sum (x) calculate_weighted_error(x, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, w1, w2);% 初始参数init_params [0.4, 63, 1.225, 0.9, 100];% 设置优化选项options optimset(MaxIter, 300, Display, iter);% 使用 fminsearch 优化optimal_params fminsearch(objective_func_weighted_sum, init_params, options);% 计算并保存结果[result1(:,turbine), result2(:,turbine), torque_errors(turbine), thrust_errors(turbine)] save_optimization_results(optimal_params, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f);end% 保存数据到附件6的表格中%result_data [Tshaft_est, Ft_est];writematrix(result1, 附件6.xlsx, Sheet, sprintf(WindTurbine_%d, 1));writematrix(result2, 附件6.xlsx, Sheet, sprintf(WindTurbine_%d, 2));% 绘制100个风机主轴扭矩和塔架推力的总误差平方和figure;subplot(2, 1, 1);bar(torque_errors);title(加权和法: 主轴扭矩总误差平方和);xlabel(风机编号);ylabel(总误差平方和);subplot(2, 1, 2);bar(thrust_errors);title(加权和法: 塔架推力总误差平方和);xlabel(风机编号);ylabel(总误差平方和);% 计算所有风机的误差平方和总和total_torque_error sum(torque_errors);total_thrust_error sum(thrust_errors);fprintf(加权和法: 所有风机主轴扭矩总误差平方和: %.4e\n, total_torque_error);fprintf(加权和法: 所有风机塔架推力总误差平方和: %.4e\n, total_thrust_error);
end% 计算与保存结果的函数
function [Tshaft_est, Ft_est,torque_error_sum_of_squares, thrust_error_sum_of_squares] save_optimization_results(optimal_params, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f)% 计算结果并与参考值对% 计算结果并与参考值对比num_time_steps size(Pref, 1);Cp optimal_params(1);blade_radius optimal_params(2);air_density optimal_params(3);generator_efficiency optimal_params(4);gear_ratio optimal_params(5);Tshaft_est zeros(num_time_steps, 1);Ft_est zeros(num_time_steps, 1);for t 1:num_time_stepsP_ref Pref(t, turbine);wind_speed V(t, turbine);omega_r_turbine omega_r(t, turbine);omega_f_turbine omega_f(t, turbine);% 主轴扭矩计算Tshaft_est(t) (P_ref / (generator_efficiency * omega_r_turbine)) * gear_ratio;Tshaft_est(t) min(max(Tshaft_est(t), 0), 10^7);% 塔架推力计算swept_area pi * blade_radius^2;Ft_est(t) 0.5 * air_density * swept_area * Cp * (wind_speed^2);Ft_est(t) min(max(Ft_est(t), 0), 10^7);end% 计算总误差平方和torque_error_sum_of_squares sum((Tshaft_est - Tshaft_ref(:, turbine)).^2);thrust_error_sum_of_squares sum((Ft_est - Ft_ref(:, turbine)).^2);% 打印结果fprintf(风机 %d 主轴扭矩总误差平方和: %.4e\n, turbine, torque_error_sum_of_squares);fprintf(风机 %d 塔架推力总误差平方和: %.4e\n, turbine, thrust_error_sum_of_squares);endfunction weighted_error calculate_weighted_error(params, turbine, Pref, V, Tshaft_ref, Ft_ref, omega_r, omega_f, w1, w2)% 提取输入参数Cp params(1); % 功率系数blade_radius params(2); % 风机叶片半径air_density params(3); % 空气密度generator_efficiency params(4); % 发电机效率gear_ratio params(5); % 齿轮箱速比% 初始化估算值num_time_steps size(Pref, 1);Tshaft_est zeros(num_time_steps, 1); % 估算的主轴扭矩Ft_est zeros(num_time_steps, 1); % 估算的塔架推力% 循环计算每个时间步的主轴扭矩和塔架推力for t 1:num_time_stepsP_ref Pref(t, turbine); % 功率调度指令wind_speed V(t, turbine); % 输入风速omega_r_turbine omega_r(t, turbine); % 低速轴转速omega_f_turbine omega_f(t, turbine); % 高速轴转速% 主轴扭矩估算 if wind_speed 0% 基于功率、发电机效率和齿轮箱速比估算主轴扭矩Tshaft_est(t) (P_ref / (generator_efficiency * omega_r_turbine)) * gear_ratio;Tshaft_est(t) min(max(Tshaft_est(t), 0), 10^9); % 限制扭矩估算值在合理范围内elseTshaft_est(t) 0;end% 塔架推力估算 swept_area pi * blade_radius^2; % 计算叶片扫过的面积Ft_est(t) 0.5 * air_density * swept_area * Cp * (wind_speed^2); % 基于风速计算推力Ft_est(t) min(max(Ft_est(t), 0), 10^7); % 限制推力估算值在合理范围内end% 计算主轴扭矩和塔架推力的误差torque_error sum((Tshaft_est - Tshaft_ref(:, turbine)).^2);thrust_error sum((Ft_est - Ft_ref(:, turbine)).^2);% 加权总误差加权和法weighted_error w1 * torque_error w2 * thrust_error;
end
...
参考文献
汤国生.基于SWOT分析法的大学生数学建模创新实践基地建设探索——以江苏科技大学为例[J].高教学刊,2023,9(11):53-56.DOI:10.19980/j.CN23-1593/G4.2023.11.013.谢雨婧,韩惠丽.北师大版初中数学教材中数学建模的多维度分析[J].教学与管理,2023(09):73-76.刘志梅.数学建模与高职数学教学的深度融合[J].佳木斯职业学院学报,2023,39(03):152-154.许亚桃,吴立宝.基于Delphi-AHP高中数学建模教学评价指标体系的研究[J].内江师范学院学报,2023,38(02):113-119.DOI:10.13603/j.cnki.51-1621/z.2023.02.018.杨本朝,石雅男,段乾恒,李光松,于刚.大学生数学建模竞赛开展全周期教学实践探究[J].大学教育,2023(04):44-46.黄健,徐斌艳.国际视野下数学建模教与学研究的发展趋势——基于第14届国际数学教育大会的分析[J].数学教育学报,2023,32(01):93-98.
...
