Pyomo 电池性能的分段线性表达式

Question

我尝试使用 pyomo 在能量优化问题中实现电池性能的分段线性表达式。我正在使用 pyomo.Piecewise 库。

我对电池充电过程的虚构性能函数具有以下形式：

y = (0.6*x - 0.2*x**2 - 0.01)

（y=充电，x=充电）

下图中 .

它显然有负 y 截距，因此在低 x 值时有负 y 值。 我的问题来了。在我下面的 mwe 中，电池充电功率 (x) 强制为大约 0.0125，充电功率 (y) 为 0，以防电池不充电。由于负 y 截距：y(0.0125)=0.

• 有没有办法在 pyomo.Piecewise 中实现这样一个方程（负截距）并利用没有二进制变量的简化线性表示？
• 是否可以只在 y>0 的范围内定义方程？但是如何实现充电功率也可以在不操作的时候为0呢？

我找不到任何方法来实现这样的方程式。万一忽略了 y 截距，实施起来很顺利。

我希望我能把我的问题说清楚。 谢谢

import pyomo.environ as pyo
import random

##% Generate some random data for PV and Load
pv = [random.randint(0, 5) for _ in range(48)]
pv_dict = (dict(enumerate(pv,1)))

load_el = [random.randint(0, 8) for _ in range(48)]
load_el_dict = (dict(enumerate(load_el,1)))

#%% 
## Define model
model = pyo.ConcreteModel()

# Define timeperiod set
model.T = pyo.RangeSet(len(pv_dict))

# Define model parameters
model.pv = pyo.Param(model.T, initialize=pv_dict)
model.load_el= pyo.Param(model.T, initialize=load_el_dict)
model.grid_cost_buy = pyo.Param(model.T, initialize=0.4)

model.battery_eoCH = pyo.Param(initialize=1.0)
model.battery_eoDCH = pyo.Param(initialize=0.1)
model.battery_capacity = pyo.Param(initialize=5)
model.battery_power_max = pyo.Param(initialize=100)

# Define the variables             
model.battery_soc = pyo.Var(model.T, bounds=(model.battery_eoDCH, model.battery_eoCH))     # battery soct with end of ch/DCH levels 
model.grid_power_import = pyo.Var(model.T, domain=pyo.NonNegativeReals)        # grid import power 
model.grid_power_export = pyo.Var(model.T, domain=pyo.NonNegativeReals)        # grid export power 
model.battery_power_DCH = pyo.Var(model.T, domain=pyo.NonNegativeReals)        # battery discharging power
# PWA variables
model.battery_power_CH_in = pyo.Var(model.T, domain=pyo.NonNegativeReals)
model.battery_power_CH_out = pyo.Var(model.T, domain=pyo.NonNegativeReals)
model.battery_power_CH_in_norm = pyo.Var(model.T, domain=pyo.NonNegativeReals, bounds=(0,1))
model.battery_power_CH_out_norm = pyo.Var(model.T, domain=pyo.NonNegativeReals)

#%% Linearization of charge efficiency

# Define function for PWA
def f(model,t,x):
    # Normalized charge performance function y=P charge out and x=P charge in
    y = (0.6*x - 0.2*x**2 - 0.01)
    return y
# Define breakpoints
breakpoints = [0.0, 0.3, 0.6, 1.0]
# Create breakpoints dict with same index as variables
PW_PTS = {}
for idx in model.battery_power_CH_in_norm.index_set():
    PW_PTS[idx] = breakpoints
    
# Define the PWA function
model.battery_power_CH_PWA_func = pyo.Piecewise(model.T, 
                                                model.battery_power_CH_out_norm, #y
                                                model.battery_power_CH_in_norm, #x
                                                pw_pts=PW_PTS, 
                                                f_rule=f, 
                                                pw_constr_type='UB', 
                                                pw_repn='CC',
                                                force_pw=False)

# Change normalized values to absolute values
def battery_power_ch_in_rule(m,t):
    return (m.battery_power_CH_in[t] == model.battery_power_CH_in_norm[t] * m.battery_power_max)
model.battery_power_ch_in_rule_c = pyo.Constraint(model.T, rule=battery_power_ch_in_rule)

def battery_power_ch_out_rule(m,t):
    return (m.battery_power_CH_out[t] == model.battery_power_CH_out_norm[t] * m.battery_power_max)
model.battery_power_ch_out_rule_c = pyo.Constraint(model.T, rule=battery_power_ch_out_rule)

#%% Further battery constraints
# Battery SoC constraint
def battery_soc_rule(m, t):
    if t == m.T.first():
        return m.battery_soc[t] == ((m.battery_power_CH_out[t] - m.battery_power_DCH[t]) / model.battery_capacity)
    
    return m.battery_soc[t] == m.battery_soc[t-1] + ((m.battery_power_CH_out[t] - m.battery_power_DCH[t]) / model.battery_capacity)
model.battery_soc_c = pyo.Constraint(model.T, rule=battery_soc_rule)

# Define balanced electricity bus rule
def balanced_bus_rule(m, t):
    return (0 == (m.pv[t] - m.load_el[t] 
                 + m.battery_power_DCH[t] - m.battery_power_CH_in[t]
                 + m.grid_power_import[t] - m.grid_power_export[t]))
model.bus_c = pyo.Constraint(model.T, rule=balanced_bus_rule)

## Define the cost function
def obj_rule(m):
    return sum(m.grid_power_import[t]*m.grid_cost_buy[t] for t in m.T)
model.obj = pyo.Objective(rule=obj_rule, sense=1)


## Solve the problem
solver = pyo.SolverFactory('gurobi')
results = solver.solve(model)
print('Total operation costs:',pyo.value(model.obj))