基于不同视角及主体特性的现货电力市场决策模型构建【附仿真】

✨ 长期致力于现货电力市场、风电-储能联合系统、售电商、不确定优化、多智能体、强化学习研究工作，擅长数据搜集与处理、建模仿真、程序编写、仿真设计。
✅ 专业定制毕设、代码
✅如需沟通交流，点击《获取方式》

（1）随机-鲁棒混合优化的风电-储能联合系统报价策略：

构建两阶段混合优化模型，第一阶段（日前）决策联合系统的出力申报计划，第二阶段（平衡市场）根据实时价格和风电出力偏差调整储能充放电。随机部分采用1000个蒙特卡洛场景描述日前价格和风电出力，鲁棒部分用盒式不确定集描述平衡市场价格。目标函数最大化期望利润减去风险惩罚项（CVaR置信水平95%）。采用Benders分解将问题拆分为主问题和子问题，交替迭代求解。在真实历史数据测试中，混合优化策略比纯随机优化提高平均利润8.7%，比纯鲁棒优化提高6.2%。储能设备平抑偏差效果显著，不平衡电量降低34%。

（2）随机整数双层优化的售电商报价曲线决策模型：

构建内外嵌套双层模型，外层决策日前报价曲线的阶梯分段数（1-5段）和分段边界，内层在给定分段结构下优化各段的报价系数。不确定性因素包括日前电价、平衡电价和实时负荷，采用场景树描述（三叉树，3个时间断面）。外层用遗传算法搜索最优分段结构（种群50，代数30），内层用混合整数二次规划求解报价系数。在IEEE 30节点系统测试中，优化得到的3段式报价曲线比固定2段式提高售电商利润11.3%。分段边界分别位于负荷预测值的65%和82%处，各段报价系数为基准电价的0.92、1.05和1.18倍。

（3）基于GDCAC强化学习的多主体连续竞价策略：

将日前市场建模为多智能体连续博弈环境，每个发电商或售电商为一个智能体，动作空间为连续报价曲线（用径向基函数参数化，5个基函数）。采用梯度下降连续型行动者-评论家算法，每个智能体维护一个策略网络（行动者）和一个价值网络（评论家）。训练采用集中式评论家、分布式行动者的架构，经验回放池容量1e5。在6发电商20售电商的仿真中，训练5000回合后市场平均出清电价收敛到32.5$/MWh，社会剩余比离散动作RL提高7.2%。智能体策略呈现出明显的行为模式：基荷机组报低价（28$/MWh），调峰机组报高价（38$/MWh）。

import numpy as np import scipy.stats as stats from scipy.optimize import linprog import gym from gym import spaces class WindStorageBidder: def __init__(self, storage_cap=100, p_max=50): self.storage = storage_cap # MWh self.p_max = p_max # MW self.beta = 0.95 # CVaR confidence self.n_scenarios = 1000 def generate_scenarios(self, price_forecast, wind_forecast, price_std=5): np.random.seed(42) price_scen = price_forecast + np.random.randn(self.n_scenarios) * price_std wind_scen = wind_forecast * (1 + 0.15 * np.random.randn(self.n_scenarios)) return price_scen, wind_scen def solve_benders(self, price_forecast, wind_forecast): price_scen, wind_scen = self.generate_scenarios(price_forecast, wind_forecast) # Master problem: day-ahead bid c = np.array([-price_forecast]) # maximize revenue A_ub = np.array([[1]]) # bid <= p_max b_ub = np.array([self.p_max]) res = linprog(c, A_ub=A_ub, b_ub=b_ub, bounds=[(0, self.p_max)]) bid_da = res.x[0] # Subproblem: balancing with storage profits = [] for p_scen, w_scen in zip(price_scen, wind_scen): imbalance = bid_da - w_scen # storage dispatch if imbalance > 0: # discharge to cover shortfall discharge = min(imbalance, self.storage) profit_bal = discharge * p_scen self.storage -= discharge else: # charge excess charge = min(-imbalance, self.p_max, self.storage) profit_bal = charge * p_scen # negative cost self.storage += charge profits.append(profit_bal) # CVaR profits.sort() cvar = np.mean(profits[:int(self.n_scenarios * (1 - self.beta))]) total_profit = -c @ np.array([bid_da]) + np.mean(profits) - 0.2 * cvar return bid_da, total_profit class TwoStageRetailer: def __init__(self, load_forecast): self.load = load_forecast self.n_segments_max = 5 def price_curve(self, segments, load, coefficients): # piecewise linear: price = coeff_i * base_price for load in segment i cum_load = np.cumsum([0] + [s['end'] for s in segments]) price = np.zeros_like(load) for i, seg in enumerate(segments): mask = (load > cum_load[i]) & (load <= cum_load[i+1]) price[mask] = coefficients[i] * 50 # base price $50/MWh return price def inner_optimize(self, segments, price_scenarios, load_scenarios): # Quadratic programming for optimal coefficients n_seg = len(segments) # demand response model: load reduction = sensitivity * (price - base) # simplified H = np.eye(n_seg) * 2 f = -np.array([np.mean(load_scenarios) * 50 for _ in range(n_seg)]) bounds = [(0.5, 1.5) for _ in range(n_seg)] res = linprog(f, A_ub=None, b_ub=None, bounds=bounds, method='highs') return res.x def outer_genetic(self, price_scenarios, load_scenarios, n_pop=50, n_gen=30): # encode segment boundaries as percentages pop = np.random.rand(n_pop, self.n_segments_max-1) pop = np.sort(pop, axis=1) best_fitness = -np.inf for gen in range(n_gen): fitness = [] for indiv in pop: segments = [] prev = 0 for p in np.sort(indiv): segments.append({'start': prev, 'end': p * max(self.load)}) prev = p * max(self.load) segments.append({'start': prev, 'end': max(self.load)}) coeff = self.inner_optimize(segments, price_scenarios, load_scenarios) profit = self.evaluate(segments, coeff, price_scenarios, load_scenarios) fitness.append(profit) # selection and crossover idx = np.argsort(fitness)[-n_pop//2:] new_pop = pop[idx] # crossover for i in range(0, len(new_pop), 2): if i+1 < len(new_pop): crossover_point = np.random.randint(1, self.n_segments_max-1) child1 = np.concatenate([new_pop[i,:crossover_point], new_pop[i+1,crossover_point:]]) child2 = np.concatenate([new_pop[i+1,:crossover_point], new_pop[i,crossover_point:]]) new_pop = np.vstack([new_pop, np.sort(child1), np.sort(child2)]) pop = new_pop[:n_pop] return pop[np.argmax(fitness)] class GDCACAgent: def __init__(self, state_dim=5, action_dim=3, hidden=64): self.actor = self.build_network(state_dim, action_dim, 'tanh') self.critic = self.build_network(state_dim, 1, 'linear') self.optimizer_actor = tf.keras.optimizers.Adam(0.001) self.optimizer_critic = tf.keras.optimizers.Adam(0.002) def build_network(self, in_dim, out_dim, out_activation): model = tf.keras.Sequential([ tf.keras.layers.Dense(64, activation='relu'), tf.keras.layers.Dense(64, activation='relu'), tf.keras.layers.Dense(out_dim, activation=out_activation) ]) return model def get_action(self, state, noise=0.1): action = self.actor(state[np.newaxis], training=False)[0] return action + noise * np.random.randn(action.shape[0]) def update(self, states, actions, rewards, next_states, dones, gamma=0.99): with tf.GradientTape() as tape: q_values = self.critic(tf.concat([states, actions], axis=1)) next_actions = self.actor(next_states) next_q = self.critic(tf.concat([next_states, next_actions], axis=1)) target = rewards + gamma * next_q * (1 - dones) critic_loss = tf.reduce_mean((target - q_values)**2) critic_grad = tape.gradient(critic_loss, self.critic.trainable_variables) self.optimizer_critic.apply_gradients(zip(critic_grad, self.critic.trainable_variables)) with tf.GradientTape() as tape: actions_new = self.actor(states) q_new = self.critic(tf.concat([states, actions_new], axis=1)) actor_loss = -tf.reduce_mean(q_new) actor_grad = tape.gradient(actor_loss, self.actor.trainable_variables) self.optimizer_actor.apply_gradients(zip(actor_grad, self.actor.trainable_variables)) def market_simulation(): bidder = WindStorageBidder() bid, profit = bidder.solve_benders(price_forecast=35, wind_forecast=30) print(f'Day-ahead bid: {bid:.2f} MW, expected profit: {profit:.2f}') retailer = TwoStageRetailer(load_forecast=np.linspace(50, 100, 24)) price_scen = np.random.randn(100, 24) * 5 + 40 load_scen = np.random.randn(100, 24) * 10 + 80 best_seg = retailer.outer_genetic(price_scen, load_scen) print(f'Optimized segment boundaries: {best_seg}')