（3）主从博弈

先贴一个主从博弈的框架，以风光机组发电为主导者，以风光的不确定性为跟随者，进行主从博弈。后边再填内容

%% 清除
clc;clear

%% 参数
num_hours=720;
p_max_wg=100;
p_max_pv=100;
%% 主程序
%初始化

% 生成随机的风速和太阳辐射数据（实际应用中应使用实际数据）
wind_speed = rand(1, num_hours);
solar_irradiance = rand(1, num_hours);

% 风能和光伏发电出力的计算
p_wg = p_max_wg * wind_speed; % 风能发电出力
p_pv = p_max_pv * solar_irradiance; % 光伏发电出力

% Stackelberg博弈的主导者和跟随者的目标函数
objective_leader = @(p_wg_sched, p_pv_sched) c_om_wg * p_wg_sched + c_om_pv * p_pv_sched + c_buy * (p_wg_sched - p_wg) + c_buy * (p_pv_sched - p_pv);
objective_follower = @(p_wg_sched, p_pv_sched) c_waste * (p_wg - p_wg_sched) + c_waste * (p_pv - p_pv_sched);

% 优化设置（实际应用中可能需要进行调整）
options = optimoptions('fmincon', 'Algorithm', 'sqp', 'Display', 'iter');

% 主导者的优化问题
x0 = [p_wg(1), p_pv(1)]; % 初始解
lb = [0, 0]; % 变量下界
ub = [p_max_wg, p_max_pv]; % 变量上界
[p_opt, cost_opt] = fmincon(@(x) objective_leader(x(1), x(2)) + objective_follower(x(1), x(2)), x0, [], [], [], [], lb, ub, [], options);

%% 输出最优解和最小成本
disp(['最优风能发电调度: ', num2str(p_opt(1)), ' kW']);
disp(['最优光伏发电调度: ', num2str(p_opt(2)), ' kW']);
disp(['最小成本: ', num2str(cost_opt), ' 元']);