优化问题中电池初始值异常：Ebat初始值不符定义的求助

问题描述 #

构建了一个基于Pyomo的电池充放电优化模型，目标是最小化外购能源成本。模型中设定初始时刻（t=1）电池能量Ebat为容量的一半（8.3/2=4.15kWh），但求解后初始Ebat值为0.415kWh（刚好是电池最小SOC对应的能量），与预期不符。

模型代码 #

import pyomo.environ as pyo
import numpy as np  # 补充原函数遗漏的导入

def self_consumption(pv, demand, prices):

    pv.index = np.arange(1, len(pv) + 1)
    demand.index = np.arange(1, len(demand) + 1)
    prices.index = np.arange(1, len(prices) + 1)

    # 定义模型
    model = pyo.ConcreteModel()

    # 定义时间步集合
    model.timesteps = pyo.Set(initialize = pyo.RangeSet(len(pv)), ordered = True)

    # 定义模型参数
    model.b_efficiency  = pyo.Param(initialize = ETTA)
    model.b_cap = pyo.Param(initialize = BATTERY_CAPACITY)
    model.b_min_soc = pyo.Param(initialize = BATTERY_SOC_MIN)
    model.b_max_soc = pyo.Param(initialize = BATTERY_SOC_MAX)
    model.b_charging_rate = pyo.Param(initialize = BATTERY_CHARGE_RATE)
    model.Ppv = pyo.Param(model.timesteps, initialize = pv.to_dict()['value'], within = pyo.Any)
    model.Pdemand = pyo.Param(model.timesteps, initialize = demand.to_dict()['value'], within = pyo.Any)
    model.day_ahead_prices = pyo.Param(model.timesteps, initialize = prices.to_dict(), within = pyo.Any)

    # 定义决策变量
    model.Pbat_ch = pyo.Var(model.timesteps, within = pyo.NonNegativeReals, bounds = (0, model.b_charging_rate))
    model.Pbat_dis = pyo.Var(model.timesteps, within = pyo.NonNegativeReals, bounds = (0, model.b_charging_rate))
    model.Ebat = pyo.Var(model.timesteps, within = pyo.NonNegativeReals, bounds = (model.b_min_soc * model.b_cap, model.b_max_soc * model.b_cap))
    model.Pgrid = pyo.Var(model.timesteps)
    model.is_charging = pyo.Var(model.timesteps, within = pyo.Binary)
    model.AbsPgrid = pyo.Var(model.timesteps, within=pyo.NonNegativeReals)

    # 定义约束
    def BatEnergyRule(model, t):
        if t == 1:
            return model.Ebat[t] == model.b_cap/2 # 初始电池能量设为容量一半
        else:  
            return model.Ebat[t] == model.Ebat[t-1] + (model.b_efficiency * model.Pbat_ch[t] - model.Pbat_dis[t]/model.b_efficiency)

    model.cons1 = pyo.Constraint(model.timesteps, rule = BatEnergyRule)

    def PowerBalanceRule(model, t):
        return model.Pgrid[t] == model.Pdemand[t] - model.Ppv[t] + model.Pbat_ch[t] - model.Pbat_dis[t]

    model.cons2 = pyo.Constraint(model.timesteps, rule = PowerBalanceRule)

    # 充放电互斥约束
    def charge_discharge_rule(model, t):
        return model.Pbat_ch[t] <= model.is_charging[t] * model.b_charging_rate

    model.charge_constraint = pyo.Constraint(model.timesteps, rule=charge_discharge_rule)

    def discharge_charge_rule(model, t):
        return model.Pbat_dis[t] <= (1 - model.is_charging[t]) * model.b_charging_rate

    model.discharge_constraint = pyo.Constraint(model.timesteps, rule=discharge_charge_rule)

    # 电网功率绝对值约束
    def AbsPgridRule1(model, t):
        return model.AbsPgrid[t] >= model.Pgrid[t]

    model.cons6 = pyo.Constraint(model.timesteps, rule=AbsPgridRule1)

    def AbsPgridRule2(model, t):
        return model.AbsPgrid[t] >= -model.Pgrid[t]

    model.cons7 = pyo.Constraint(model.timesteps, rule=AbsPgridRule2)

    # 目标函数：最小化外购能源成本
    def ObjRule(model):
        return sum(model.day_ahead_prices[t] * model.AbsPgrid[t] for t in model.timesteps)

    model.obj = pyo.Objective(rule=ObjRule, sense=pyo.minimize)
    
    # 选择求解器
    opt = pyo.SolverFactory('glpk')

    # 求解优化问题
    result = opt.solve(model, tee = True)

    solution = {
        "Pbat_ch": {t: model.Pbat_ch[t].value for t in model.timesteps},
        "Pbat_dis": {t: model.Pbat_dis[t].value for t in model.timesteps},
        "Ebat": {t: model.Ebat[t].value for t in model.timesteps},
        "Pgrid": {t: model.Pgrid[t].value for t in model.timesteps},
        }
    
    return solution

运行代码 #

import numpy as np
import pandas as pd

# 电池参数
BATTERY_CAPACITY = 8.3 # 电池容量（kWh）
BATTERY_CHARGE_RATE = 2.6  # 充放电功率（kW）
BATTERY_SOC_MAX = 1 # 最大SOC比例
BATTERY_SOC_MIN = 0.05 # 最小SOC比例
ETTA = 0.96 # 电池充放电效率

# 生成模拟数据
pv = pd.DataFrame({'value': np.random.uniform(0, 5, 24)})
demand = pd.DataFrame({'value': np.random.uniform(0, 5, 24)})
prices = pd.Series(np.random.uniform(0.1, 0.5, 24), name='prices')

# 运行优化
optimization_results = self_consumption(pv, demand, prices)

原因分析 #

核心问题是初始约束中使用Param对象的除法表达式model.b_cap/2时，Pyomo未正确将其解析为固定数值（4.15），导致初始约束未生效。求解器在无法满足失效约束的情况下，选择了满足变量边界条件的可行解——即电池最小SOC对应的能量值（0.05*8.3=0.415kWh）。

另外，原模型函数中遗漏了numpy的导入语句，虽然运行代码中已全局导入，但函数内直接使用np存在潜在报错风险，需补充。

解决方案 #

方案1：直接使用固定数值设置初始约束 #

修改BatEnergyRule函数中t=1的分支，直接用全局变量的数值定义初始能量，避免Param表达式解析问题：

def BatEnergyRule(model, t):
    if t == 1:
        return model.Ebat[t] == BATTERY_CAPACITY / 2  # 直接使用固定数值4.15
    else:  
        return model.Ebat[t] == model.Ebat[t-1] + (model.b_efficiency * model.Pbat_ch[t] - model.Pbat_dis[t]/model.b_efficiency)

方案2：将初始能量定义为独立Param #

在模型的Param定义部分添加初始能量参数，再用于约束：

# 在Param定义区域添加
model.b_initial_energy = pyo.Param(initialize=BATTERY_CAPACITY / 2)

# 修改BatEnergyRule
def BatEnergyRule(model, t):
    if t == 1:
        return model.Ebat[t] == model.b_initial_energy
    else:  
        return model.Ebat[t] == model.Ebat[t-1] + (model.b_efficiency * model.Pbat_ch[t] - model.Pbat_dis[t]/model.b_efficiency)

额外优化：补充函数内的numpy导入 #

在self_consumption函数开头添加：

import numpy as np

验证方法 #

修改后重新运行模型，查看optimization_results["Ebat"][1]的值，应为预期的4.15kWh。

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