优化求解器Cplex（Python）在微网（源网荷储的简单场景）能源系统优化应用

在Python中通过Cplex求解器对简单微网多状态建模求解（入门版），欢迎学习讨论。

笙枫

1346人浏览 · 2024-10-18 17:14:18

笙枫 · 2024-10-18 17:14:18 发布

因为在公司开发项目（新加坡），对于商用软件的使用较为敏感。所在项目需要使用到Cplex求解器，先利用python中可以直接使用的免费cplex库（限1000/1000）进行了简单的学习。

本文内容主要参考了此博客：Cplex +Python下的虚拟电厂智能调度实例（附代码）_cplex python 例子-CSDN博客

ps：因为AI的发展，现在的应用编程基本对代码的要求较低，初学可以利用一些AI软件进行辅助。

1. 借助conda软件（因为公司没买，这里我用的是miniconda）安装Python3.8，如果是vscode可以利用jupyter内核选择快速创建（点击右侧内核选择，然后点击中间选择其他内核，选择3.8版本）

2. 我这里conda库不能直接conda install，所以我是选定好环境后在notebook利用pip安装并重启内核。

pip install cplex docplex

以下是本次优化任务（初版来源于开头提到的参考博客）

虚拟电厂单用户调度（含储能）

光储荷经济性调度，针对单个用户、单个光伏、单个储能设备，短周期内的经济最优目标

变量

充电变量： $(x_t) \in [0, \alpha]$
放电变量： $(y_t) \in [-\alpha, 0]$
储能设备能量状态： $(s_t) \in [0, 1]$

以上三个变量为 continuous_var_list 连续变量列表

参数

负荷：load
功率：power
电价：price
充电效率：charge_efficiency
放电效率：discharge_efficiency
储能单次充放电限制：nominal_power

目标函数

最小化经济性最优目标：

$\text{minimize } Cost = \sum_{t=1}^T (load(t) - power(t) + x(t) + y(t)) \times price(t) \div costbase$

其中 costbase 为不使用储能的情况下的花费：
$\sum_{t=1}^T (load(t) - power(t)) \times price(t)$

约束条件

储能能量状态：
- 储能初始状态为 0： $s_0 = 0$
- 储能状态更新： $st=st−1+c1⋅xt+ytc2s_t = s_{t-1} + c_1 \cdot x_t + \frac{y_t}{c_2}$

from docplex.mp.model import Model
import random

#设置一个随机种子
random.seed(1234)

# 创建模型，在这里我们对连续10个状态进行动态优化,电价有4个时段
time = 4
time1 = 5

# 生成随机的用户负荷和用户功率
user_loads = [round(random.uniform(0, 1), 2) for _ in range(time1)]
user_powers = [round(random.uniform(0, 1), 2) for _ in range(time1)]

elec_price = [0.54, 0.22, 0.32, 0.24]

#设置储能充放电效率和单次充放限制
charge_eff = 0.91
discharge_eff = 0.95

nominal_power = 0.8

model = Model(name='Electricity')
x = model.continuous_var_list(time, lb=0,ub=nominal_power,name='x')
y = model.continuous_var_list(time, lb = -nominal_power,ub = 0, name='y')
soc = model.continuous_var_list(time1,lb=0,ub=nominal_power,name='soc')

#无储能下的花费
cost_base = sum(((user_loads[i] - user_powers[i]) * elec_price[i]) for i in range(time))
#总花费
total_cost = model.sum(((user_loads[i] - user_powers[i] + x[i] + y[i]) * elec_price[i] / cost_base) for i in range(time))

model.minimize(total_cost)
model.add_constraint(soc[0] == 0)

for i in range(time):
    model.add_constraint(soc[i+1] == soc[i] + x[i] * charge_eff + y[i]/discharge_eff)
model.solve()

print('The optimized solution is:')
print("charging(+)/discharging(-) power for each time step:")
for i in range(time):
    print("time{}:solution:".format(i),x[i].solution_value+y[i].solution_value)

for i in range(time1):
    print('time{}:soc{}'.format(i,soc[i].solution_value))

print('Total cost:{}'.format(total_cost.solution_value))

4. 当我们把时间和状态扩展到24小时，每10分钟采样一次并指定策略，将电价设置为实时波动，可以将以上代码进行如下改进：

from docplex.mp.model import Model
import random
import matplotlib.pyplot as plt

# 设置一个随机种子
random.seed(1234)

# 创建模型，在这里我们对一天24小时进行动态优化，每10分钟采样一次
sampling_interval = 10  # 每10分钟采样一次
minutes_per_day = 24 * 60  # 一天共有 1440 分钟
time = minutes_per_day // sampling_interval  # 采样点数量

# 生成随机的用户负荷和用户功率
user_loads = [round(random.uniform(0, 1), 2) for _ in range(time)]
user_powers = [round(random.uniform(0, 1), 2) for _ in range(time)]

# 生成随机电价，基于 0.5 生成 ±0.2 的波动
elec_price = [round(0.5 + random.uniform(-0.2, 0.2), 2) for _ in range(time)]

# 设置储能充放电效率和单次充放限制
charge_eff = 0.91
discharge_eff = 0.95
nominal_power = 0.8

# 创建模型
model = Model(name='Electricity')

# 创建变量
x = model.continuous_var_list(time, lb=0, ub=nominal_power, name='x')
y = model.continuous_var_list(time, lb=-nominal_power, ub=0, name='y')
soc = model.continuous_var_list(time, lb=0, ub=nominal_power, name='soc')

# 无储能下的花费
cost_base = sum(((user_loads[i] - user_powers[i]) * elec_price[i]) for i in range(time))

# 总花费
total_cost = model.sum(((user_loads[i] - user_powers[i] + x[i] + y[i]) * elec_price[i] / cost_base) for i in range(time))

# 设置目标函数
model.minimize(total_cost)

# 添加约束条件
model.add_constraint(soc[0] == 0)
for i in range(time - 1):
    model.add_constraint(soc[i + 1] == soc[i] + x[i] * charge_eff + y[i] / discharge_eff)

# 求解模型
solution = model.solve()

# 检查是否找到最优解
if solution:
    print('The optimized solution is:')
    print("charging(+)/discharging(-) power for each time step:")
    charging_discharge_values = []
    soc_values = []
    
    for i in range(time):
        charge_discharge = x[i].solution_value + y[i].solution_value
        charging_discharge_values.append(charge_discharge)
        soc_values.append(soc[i].solution_value)
        if i % (60 // sampling_interval) == 0:  # 每小时输出一次状态
            print(f"Hour {i // (60 // sampling_interval)}: solution: {charge_discharge}")
    
    # 可视化结果
    hours = [i * sampling_interval / 60 for i in range(time)]
    
    plt.figure(figsize=(14, 10))
    
    # 充放电功率可视化
    plt.subplot(2, 1, 1)
    plt.plot(hours, charging_discharge_values, label='Charging/Discharging Power', color='b')
    plt.xlabel('Time (hours)')
    plt.ylabel('Power (kW)')
    plt.title('Charging (+) and Discharging (-) Power Over Time')
    plt.legend()

    # SOC 状态可视化
    plt.subplot(2, 1, 2)
    plt.plot(hours, soc_values, label='State of Charge (SOC)', color='g')
    plt.xlabel('Time (hours)')
    plt.ylabel('SOC')
    plt.title('State of Charge Over Time')
    plt.legend()

    plt.tight_layout()
    plt.show()

    print(f'Total cost: {total_cost.solution_value}')
else:
    print("No solution found.")

以上代码输出求解结果：
The optimized solution is:
charging(+)/discharging(-) power for each time step:
Hour 0: solution: 0.10839999999999994
Hour 1: solution: 0.2851999999999997
Hour 2: solution: 0.8
Hour 3: solution: 0.10840000000000005
Hour 4: solution: 0.8
Hour 5: solution: 0.8
Hour 6: solution: 0.8
Hour 7: solution: -0.7537304800462696
Hour 8: solution: 0.0
Hour 9: solution: -0.7537304800462696
Hour 10: solution: 0.0
Hour 11: solution: -0.6746096009253907
Hour 12: solution: 0.10839999999999994
Hour 13: solution: 0.17679999999999985
Hour 14: solution: 0.0
Hour 15: solution: 0.0
Hour 16: solution: 0.17679999999999985
Hour 17: solution: 0.10839999999999994
Hour 18: solution: -0.7537304800462696
Hour 19: solution: 0.0
Hour 20: solution: -0.6283400809716602
Hour 21: solution: 0.8
Hour 22: solution: 0.0
Hour 23: solution: 0.10839999999999994

我们可视化充放电状态和SOC状态：
在这里插入图片描述

此时我们能了解到单个系统内，一负载一发电光伏一储能并连接电网的简单模型优化，在此基础上长期优化还应当考虑以下约束：

增加复杂的约束条件：如电网容量、最小 SOC、电网平衡、冷却时间等。
改进目标函数：采用多目标优化，包括最小化总成本、最大化收益、平滑负荷等。
考虑不确定性和随机性：加入随机电价、负荷等因素，通过鲁棒优化提高模型的适用性。
引入实际电力市场的特性：如峰谷电价、惩罚机制等，确保模型能贴近实际运营情况。
设备寿命和维护：减少设备磨损，考虑储能充放电对设备寿命的影响。

比如我们这里可以考虑加入储能SOC下限（一般来说商用设备是在0.2~0.8间能够较好地延长使用寿命），并考虑到充放的冷却间隔（防止过热）

# 添加约束条件
min_soc = 0.2 * nominal_power  # 最小 SOC 限制
model.add_constraint(soc[0] == 0.5 * nominal_power)  # 初始 SOC 设置为 0.5

for i in range(time - 1):
    # SOC 状态更新约束
    model.add_constraint(soc[i + 1] == soc[i] + x[i] * charge_eff + y[i] / discharge_eff)
    # 最小 SOC 约束
    model.add_constraint(soc[i] >= min_soc)

# 放电后的冷却时间约束（通过引入额外的二进制变量来表示冷却状态）
for i in range(time - 2):
    # 创建二进制变量来表示冷却状态
    cooling_20 = model.binary_var(name=f'cooling_20_{i}')
    cooling_10 = model.binary_var(name=f'cooling_10_{i}')
    
    # 冷却 20 分钟约束，当放电量超过 50% 时触发冷却
    model.add_indicator(cooling_20, y[i] <= -0.5 * nominal_power, 1)
    model.add_constraint(x[i + 1] <= nominal_power * (1 - cooling_20))  # 冷却期间不可充电
    
    # 冷却 10 分钟约束，当放电量在 30% - 50% 之间时触发冷却
    model.add_indicator(cooling_10, y[i] <= -0.3 * nominal_power, 1)
    model.add_constraint(x[i + 1] <= nominal_power * (1 - cooling_10))  # 冷却期间不可充电

在这里插入图片描述
可以看到这里的约束添加在随机数据中表现并不理想，可能原因很多，比如随机数据的波动让优化求解面临局限；约束的添加不合理等等。

天启AI社区

GitCode 天启AI是一款由 GitCode 团队打造的智能助手，基于先进的LLM（大语言模型）与多智能体 Agent 技术构建，致力于为用户提供高效、智能、多模态的创作与开发支持。它不仅支持自然语言对话，还具备处理文件、生成 PPT、撰写分析报告、开发 Web 应用等多项能力，真正做到“一句话，让 Al帮你完成复杂任务”。

更多推荐