基于Scipy罚函数的风电机组启停优化及PyWake仿真实现咨询

1 优化逻辑梳理 #

你的需求属于典型的0-1整数规划问题，优化目标为：
最大化 全场总发电量(AEP) - 约束罚项
由于Scipy优化器默认执行最小化求解，我们将目标转换为最小化 -(AEP - 罚项) 即可，其中罚项根据你的实际约束自定义（比如开机数量上限、电网出力限制等，违反约束时罚项为正，降低目标得分）。

2 完整实现代码 #

import numpy as np
from scipy.optimize import milp, Bounds
from py_wake.examples.data.hornsrev1 import V80 
from py_wake.examples.data.hornsrev1 import Hornsrev1Site
from py_wake import BastankhahGaussian
import function

# 初始化风场和模型参数
site = Hornsrev1Site()
x, y = site.initial_position.T
windTurbines = V80()
wf_model = BastankhahGaussian(site, windTurbines)
wt_count = len(x) # 风机总数量

# 自定义罚函数，可根据实际约束修改
def calculate_penalty(c):
    penalty = 0
    # 示例约束：最多允许开机60台，超出部分每台罚100（罚值要远大于单台风机的AEP贡献，保证约束生效）
    max_on_wt = 60
    on_count = np.sum(c)
    if on_count > max_on_wt:
        penalty += 100 * (on_count - max_on_wt)
    # 可添加更多约束的罚项计算逻辑
    return penalty

# 优化目标函数
def objective(c):
    # 转换为整数类型（避免求解器返回浮点值）
    c = np.round(c).astype(int)
    # 过滤停机风机坐标
    x_new, y_new = function.funC(x, y, c)
    # 仿真计算AEP
    sim_res = wf_model(x_new, y_new, h=None, type=0, wd=None, ws=None)
    total_aep = sim_res.aep().sum().item()
    # 计算罚项
    penalty = calculate_penalty(c)
    # 返回负的净收益，适配最小化求解逻辑
    return -(total_aep - penalty)

# ----------------------
# 方案1：使用MILP求解（适合中小规模风场，风机数<100）
# ----------------------
# 设置变量边界：每个c的取值为0或1，整数类型
bounds = Bounds(lb=np.zeros(wt_count), ub=np.ones(wt_count))
integrality = np.ones(wt_count, dtype=bool) # 所有变量都是整数

# 执行优化
res = milp(objective, bounds=bounds, integrality=integrality, options={"disp": True})

# 输出结果
optimal_c = np.round(res.x).astype(int)
print("最优启停序列c：", optimal_c)
print("最优净发电量：", -res.fun)

# ----------------------
# 方案2：使用差分进化求解（适合大规模风场，鲁棒性更强）
# ----------------------
'''
from scipy.optimize import differential_evolution
# 设置每个变量的取值范围
bounds_de = [(0, 1)] * wt_count
res_de = differential_evolution(objective, bounds=bounds_de, integrality=np.ones(wt_count, dtype=bool), disp=True)
optimal_c = np.round(res_de.x).astype(int)
print("最优启停序列c：", optimal_c)
print("最优净发电量：", -res_de.fun)
'''

3 关键配置说明 #

• 罚函数调整：你可以根据实际业务约束修改calculate_penalty的逻辑，比如新增风机检修约束、区域出力限制等，罚值需设置为大于单台风机的AEP贡献，确保违反约束的解会被优先排除
• 求解器选择：风机数量少于80台时优先用MILP求解，速度更快；风机数量更多时切换为差分进化方案，避免求解时间过长
• 结果验证：得到最优c序列后可代入你原有的仿真逻辑，直接验证发电量输出是否符合预期

内容的提问来源于stack exchange，提问作者sadra