单时段机组组合优化的粒子群算法实现（MATLAB）

使用粒子群算法(PSO)解决单时段机组组合问题的MATLAB实现。机组组合问题是在满足电力需求的前提下，确定哪些发电机组应该运行以及它们的出力水平，以最小化总成本。

% 单时段机组组合优化的粒子群算法实现
clear; clc; close all;

%% 机组参数：最小出力、最大出力、成本系数(a, b, c)
units = struct();
units(1).p_min = 10;  units(1).p_max = 50;  units(1).cost_coeff = [0.01, 0.5, 10];
units(2).p_min = 20;  units(2).p_max = 80;  units(2).cost_coeff = [0.02, 0.3, 20];
units(3).p_min = 30;  units(3).p_max = 100; units(3).cost_coeff = [0.01, 0.4, 30];
units(4).p_min = 40;  units(4).p_max = 120; units(4).cost_coeff = [0.03, 0.2, 40];
units(5).p_min = 50;  units(5).p_max = 150; units(5).cost_coeff = [0.01, 0.6, 50];

% 系统参数
load_demand = 300;  % 负荷需求
n_units = length(units);  % 机组数量

%% PSO参数
n_particles = 50;     % 粒子数量
max_iter = 100;       % 最大迭代次数
w = 0.7;              % 惯性权重
c1 = 1.5;             % 个体学习因子
c2 = 1.5;             % 群体学习因子

%% 初始化粒子群
particles = zeros(n_particles, n_units);
for i = 1:n_particles
    for j = 1:n_units
        % 随机生成机组出力，在最小和最大出力之间
        particles(i, j) = units(j).p_min + rand() * (units(j).p_max - units(j).p_min);
    end
end

% 初始化速度
velocities = zeros(n_particles, n_units);

%% 计算机组成本
function cost = calculate_cost(particle, units)
    cost = 0;
    for i = 1:length(particle)
        a = units(i).cost_coeff(1);
        b = units(i).cost_coeff(2);
        c = units(i).cost_coeff(3);
        p = particle(i);
        % 成本函数: a*p^2 + b*p + c
        cost = cost + a * p^2 + b * p + c;
    end
end

%% 计算约束违反程度
function violation = calculate_constraint_violation(particle, units, load_demand)
    % 计算总出力
    total_output = sum(particle);
    
    % 计算出力平衡约束违反程度
    balance_violation = abs(total_output - load_demand);
    
    % 计算机组出力限制约束违反程度
    limit_violation = 0;
    for i = 1:length(particle)
        p_min = units(i).p_min;
        p_max = units(i).p_max;
        if particle(i) < p_min
            limit_violation = limit_violation + (p_min - particle(i));
        elseif particle(i) > p_max
            limit_violation = limit_violation + (particle(i) - p_max);
        end
    end
    
    violation = balance_violation + limit_violation;
end

%% 适应度函数
function fit = fitness(particle, units, load_demand)
    penalty_factor = 1000;  % 罚函数系数
    cost = calculate_cost(particle, units);
    violation = calculate_constraint_violation(particle, units, load_demand);
    fit = cost + penalty_factor * violation;
end

%% PSO主算法
% 初始化个体最优位置和适应度
personal_best_positions = particles;
personal_best_fitness = zeros(n_particles, 1);
for i = 1:n_particles
    personal_best_fitness(i) = fitness(particles(i, :), units, load_demand);
end

% 初始化全局最优位置和适应度
[global_best_fitness, global_best_index] = min(personal_best_fitness);
global_best_position = personal_best_positions(global_best_index, :);

% 记录每次迭代的最佳适应度
fitness_history = zeros(max_iter, 1);

% PSO迭代
for iter = 1:max_iter
    for i = 1:n_particles
        % 更新速度
        r1 = rand(1, n_units);
        r2 = rand(1, n_units);
        velocities(i, :) = w * velocities(i, :) + ...
                           c1 * r1 .* (personal_best_positions(i, :) - particles(i, :)) + ...
                           c2 * r2 .* (global_best_position - particles(i, :));
        
        % 更新位置
        particles(i, :) = particles(i, :) + velocities(i, :);
        
        % 计算新位置的适应度
        current_fitness = fitness(particles(i, :), units, load_demand);
        
        % 更新个体最优
        if current_fitness < personal_best_fitness(i)
            personal_best_fitness(i) = current_fitness;
            personal_best_positions(i, :) = particles(i, :);
            
            % 更新全局最优
            if current_fitness < global_best_fitness
                global_best_fitness = current_fitness;
                global_best_position = particles(i, :);
            end
        end
    end
    
    % 记录本次迭代的最佳适应度
    fitness_history(iter) = global_best_fitness;
    
    % 打印当前迭代信息
    if mod(iter, 10) == 0
        fprintf('迭代次数: %d/%d, 最佳适应度: %.2f\n', iter, max_iter, global_best_fitness);
    end
end

%% 输出结果
fprintf('\n优化结果:\n');
fprintf('总成本: %.2f\n', global_best_fitness);
fprintf('各机组出力:\n');
total_output = 0;
for i = 1:n_units
    p_min = units(i).p_min;
    p_max = units(i).p_max;
    output = min(max(global_best_position(i), p_min), p_max);  % 确保在限制范围内
    total_output = total_output + output;
    fprintf('机组 %d: %.2f MW (限制: %d-%d MW)\n', i, output, p_min, p_max);
end
fprintf('总出力: %.2f MW, 负荷需求: %d MW\n', total_output, load_demand);

%% 绘制收敛曲线
figure;
plot(1:max_iter, fitness_history, 'LineWidth', 2);
title('PSO收敛曲线');
xlabel('迭代次数');
ylabel('最佳适应度(总成本)');
grid on;

算法

1. 问题建模：

   • 每个粒子代表一个可能的解决方案，即各机组的出力水平
   • 目标是最小化总发电成本，同时满足负荷需求和机组出力限制

2. 成本函数：

   • 每个机组的发电成本建模为二次函数：\(C_i(P_i) = a_i P_i^2 + b_i P_i + c_i\)
   • 总成本是所有机组成本之和

3. 约束处理：

   • 使用罚函数法处理约束条件
   • 约束包括：总出力等于负荷需求、各机组出力在最小和最大出力范围内

4. PSO参数：

   • 惯性权重w=0.7
   • 个体学习因子c1=1.5
   • 群体学习因子c2=1.5
   • 粒子数量50，最大迭代次数100

结果

算法运行后会输出：

1. 最优解的总成本
2. 各机组的最优出力分配
3. 收敛曲线，显示算法随迭代次数的优化过程

参考代码 PSO单时段机组寻优 www.youwenfan.com/contentcnj/54823.html

注意

1. 此代码实现了基本的PSO算法，对于复杂的机组组合问题，可能需要更精细的约束处理机制。
2. 罚函数系数的选择会影响优化结果，可能需要根据具体问题调整。
3. 对于大规模问题，可能需要增加粒子数量和迭代次数以获得更好的解。
4. 实际应用中，可能需要考虑机组的启停成本、爬坡率限制等更多约束条件。