- Introduction
A microgrid system incorporating wind power, photovoltaic (PV) generation, energy storage, and micro gas turbine is considered in this study. The optimization objective is to minimize the total operational cost, which includes grid power purchase costs, maintenance costs for PV panels and wind turbines, battery charging/discharging maintenance costs, gas turbine operating costs, and pollution gas treatment costs. The optimization model is formulated with the following constraints:
- Power balance constraint
- Gas turbine ramp rate constraint
- Grid exchange power constraint
- Energy storage device constraint
- Controllable micro power source output constraint
- System Model
2.1 Photovoltaic System
The open-circuit voltage of a PV cell is related to solar spectral irradiance and is independent of the cell area. Under a solar spectral irradiance of 100MW/cm², the open-circuit voltage of monocrystalline silicon PV cells ranges from 450-600mV, with a maximum of 690mV. When the incident spectral irradiance changes, the open-circuit voltage is proportional to the logarithm of the spectral irradiance. As ambient temperature increases, the open-circuit voltage decreases approximately 2-3mV per 1°C increase.
The relationship between PV cell parameters can be expressed through the equivalent circuit equations:
Where I₀ represents the reverse saturation current of the equivalent diode PN junction (typically constant), U_D is the diode terminal voltage, q is the electron charge, k is Boltzmann's constant, T is thermodynamic temperature, and A is the PN结 curve constant.
The V-I characteristics of PV cells depend on irradiance intensity and cell temperature. The mathematical model uses irradiance and cell temperature as parameters. When the actual temperature T and irradiance S differ from reference conditions, corrections must be applied to the model based on reference values of I_sc, U_oc, I_m, and U_m.
2.2 Wind Power Generation
The output power of a wind turbine can be expressed as:
Where P_wt is the wind turbine output power, ρ is air density, V is wind speed, R is the blade radius, and C_p is the wind energy utilization coefficient.
The relationship between wind turbine output and wind speed is:
Where P_rated is rated power, V_rated is rated wind speed, V_cutin is cut-in wind speed, and V_cutout is cut-out wind speed.
2.3 Micro Gas Turbine
A hydrogen-fueled micro gas turbine is selected as the energy conversion device for hydrogen-electricity coupling. This equipment uses hydrogen as fuel, with water as the only combustion product, achieving zero emissions.
- Objective Function
The objective function represents the microgrid economic operation cost, which includes:
- Power exchange cost with the main grid
- Maintenance costs for PV, wind turbines, and batteries
- Gas turbine operating costs
- Pollution gas treatment costs
- Gas turbine startup costs
- Solution Algorithm
The Particle Swarm Optimization (PSO) algorithm is employed, where each particle represents an individual in the N-dimensional search space. The particle's current position corresponds to a candidate solution, and the flight process represents the search processs. The velocity can be dynamically adjusted based on the particle's historical best position and the swarm's global best position.
Each particle has two attributes: velocity (representing movement speed) and position (representing movement direction). The best solution found by each particle individually is called the personal best, while the best among all personal bests becomes the global best. Through iterative updates of velocity and position, the algorithm converges to an optimal solution.
Algorithm Flow:
- Initialization: Set maximum iterations, number of decision variables, maximum particle velocity, and search space bounds. Initialize velocities and positions randomly within the specified ranges. Set the swarm size and assign random initial velocities to each particle.
- Personal Best and Global Best: Define a fitness function. The personal best is the optimal solution found by each particle. The global best is the best among all personal bests. Compare with historical global best for updates.
- Update Velocity and Position
- Termination Criteria:
- Maximum iterations reached
- Difference between generations below minimum threshold
Algorithm Enhancement: This study improves upon standard PSO by adaptively adjusting the inertia weight coefficient based on the global best position. The weight dynamically varies according to particle position rather than remaining fixed.
- Matlab Implementation
The following code demonstrates the hybrid PSO implementation for microgrid day-ahead scheduling:
%% Initialization
clc;
clear all;
close all;
%% System Parameters
flag_pso = 1;
population_size = 80;
battery_max_output = 400000; % Maximum battery output: 400kW
battery_min_output = -150000; % Minimum battery output: -150kW
pv_max_output = 100000; % Maximum PV output: 100kW
wind_max_output = 100000; % Maximum wind output: 100kW
mt_max_output = 50000; % Maximum gas turbine output: 50kW
mt_min_output = 5000; % Minimum continuous operation output: 5kW
grid_max_exchange = 150000; % Maximum grid exchange: 150kW
%% Time Configuration
step_duration = 3600; % Step duration in seconds
hours_simulation = 24; % Simulation period: 24 hours
total_simulation = hours_simulation * 1;
total_units = 5; % Five generation units
%% Initialize Power Bus
power_bus = zeros(total_units, hours_simulation);
%% Unit Indexing
idx_mt = 1;
idx_pv = 2;
idx_wind = 3;
idx_battery = 4;
idx_load = 5;
idx_grid = 6;
initial_vector(1:hours_simulation) = 0.5;
%% Micro Gas Turbine Parameters
cost_mt = 0.3571; % Fuel cost coefficient
flow_rate = 0.0085; % Gas consumption rate (m³/Wh)
time_step = 1; % Optimization time step (hours)
startup_cost = 0.0587; % Startup cost ($/m³)
hot_start = 30; % Hot startup time (seconds)
cold_start = 200; % Cold startup time (seconds)
cooling_duration = 520; % Cooling time (seconds)
min_up_time = 600; % Minimum up time (seconds)
min_down_time = 300; % Minimum down time (seconds)
time_on = 0;
min_on_time = 2;
time_off = 0;
min_off_time = 2;
status_mt = ones(1, hours_simulation); % Assume MT always running
ramp_up_limit = 75000; % Ramp-up rate limit (75kWh)
ramp_down_limit = 75000; % Ramp-down rate limit (75kWh)
%% Wind Turbine Parameters
cost_wind = 0.0001; % Wind power cost ($/Wh)
status_wind = ones(1, total_simulation);
%% Battery Parameters
state_charge = zeros(1, total_simulation);
state_charge(1, :) = 0.5; % Initialize SOC to 50%
status_battery = ones(1, total_simulation);
min_charge_power = 1000; % Minimum charge power
max_charge_power = 16000; % Maximum charge power
battery_efficiency = 0.9; % Battery efficiency
%% PV System Parameters
cost_pv = 0.2; % PV cost ($/kWh)
%% Load Parameters
load_min = 150 * 10^3; % Minimum load demand
load_max = 300 * 10^3; % Maximum load demand
%% Grid Parameters
grid_price = 1.05;
counter_dispatch = 0;
grid_power = zeros(1, total_simulation);
elapsed_time = 0;
%% Environmental Data
solar_irradiance = zeros(total_simulation, 1);
ambient_temperature = zeros(total_simulation, 1);
wind_speed_data = zeros(total_simulation, 1);
%% Define Daily Profiles
solar_irradiance = [0 0 0 0 0 0 0 0.2 0.3 0.4 0.6 0.8 0.8 0.8 0.8 0.7 0.6 0.4 0.3 0.1 0 0 0 0]';
ambient_temperature = [10 10 9 9 8 9 10 12 15 18 20 22 25 25 25 25 24 21 18 16 14 12 11 10]';
wind_speed_data = [10 11 12 12 12 12 12 13 15 16 17 17 18 19 19 20 21 22 23 24 25 26 27 28]';
%% Main Simulation Loop
for time_index = 1:total_simulation
elapsed_time = elapsed_time + 1;
% Load demand assignment based on time-of-use
if elapsed_time <= hours_simulation
if (elapsed_time == 1 || elapsed_time == 8 || elapsed_time == 9 || ...
elapsed_time == 16 || elapsed_time == 17 || elapsed_time == 24 || ...
elapsed_time == 23 || elapsed_time == 22)
power_bus(idx_load, time_index) = 30000;
elseif (elapsed_time == 11 || elapsed_time == 12 || elapsed_time == 15 || elapsed_time == 19)
power_bus(idx_load, time_index) = 60000;
elseif (elapsed_time == 10 || elapsed_time == 18)
power_bus(idx_load, time_index) = 90000;
elseif (elapsed_time == 21)
power_bus(idx_load, time_index) = 220000;
elseif (elapsed_time == 14 || elapsed_time == 20)
power_bus(idx_load, time_index) = 180000;
elseif (elapsed_time == 13)
power_bus(idx_load, time_index) = 280000;
else
power_bus(idx_load, time_index) = 8000;
end
else
% For simulation hours beyond 24
time_mod = mod(elapsed_time, hours_simulation);
if (time_mod == 1 || time_mod == 8 || time_mod == 9 || ...
time_mod == 16 || time_mod == 17 || time_mod == 24)
power_bus(idx_load, time_index) = 30000;
elseif (time_mod == 11 || time_mod == 12 || time_mod == 15 || ...
time_mod == 18 || time_mod == 19 || time_mod == 23)
power_bus(idx_load, time_index) = 60000;
elseif (time_mod == 10)
power_bus(idx_load, time_index) = 90000;
elseif (time_mod == 21 || time_mod == 22)
power_bus(idx_load, time_index) = 150000;
elseif (time_mod == 14 || time_mod == 20)
power_bus(idx_load, time_index) = 210000;
elseif (time_mod == 13)
power_bus(idx_load, time_index) = 300000;
else
power_bus(idx_load, time_index) = 3000;
end
end
% Calculate renewable energy output
power_bus(idx_wind, time_index) = compute_wind_power(...
wind_speed_data(time_index), status_wind(time_index), wind_max_output);
power_bus(idx_pv, time_index) = compute_pv_output(...
solar_irradiance(time_index), ambient_temperature(time_index), pv_max_output);
% Initialize optimization at end of each day
if mod(elapsed_time, hours_simulation) == 0
if flag_pso == 1
temp_soc = state_charge;
% Initialize PSO swarm
[decision_vars, grid_exchange, soc_values, mt_status] = initialize_population(...
power_bus(idx_wind, time_index-hours_simulation+1:time_index), ...
power_bus(idx_pv, time_index-hours_simulation+1:time_index), ...
power_bus(idx_load, time_index-hours_simulation+1:time_index), ...
state_charge(time_index-hours_simulation+1:time_index), ...
status_mt(time_index-hours_simulation+1:time_index), ...
battery_max_output, battery_min_output, mt_max_output, mt_min_output, ...
grid_max_exchange, population_size);
end
end
end
%% Wind Power Calculation Function
function power_output = compute_wind_power(wind_velocity, unit_status, rated_power)
cut_in = 3; % Cut-in speed (m/s)
rated_velocity = 12; % Rated speed (m/s)
cut_out = 25; % Cut-out speed (m/s)
if unit_status == 0 || wind_velocity < cut_in || wind_velocity > cut_out
power_output = 0;
elseif wind_velocity >= cut_in && wind_velocity < rated_velocity
power_output = rated_power * ((wind_velocity - cut_in) / (rated_velocity - cut_in))^3;
else
power_output = rated_power;
end
end
%% PV Output Calculation Function
function pv_output = compute_pv_output(irradiance, temperature, max_power)
reference_irradiance = 1000; % Reference irradiance (W/m²)
reference_temp = 25; % Reference temperature (°C)
% Temperature coefficient
temp_coeff = -0.004; % Power temperature coefficient
if irradiance <= 0
pv_output = 0;
else
% Current temperature correction
temp_variation = temperature - reference_temp;
output_factor = irradiance / reference_irradiance * (1 + temp_coeff * temp_variation);
pv_output = max_power * output_factor;
end
end
- Conclusion
This article presents a hybrid particle swarm optimization approach for solving the day-ahead optimal scheduling problem in a wind-solar-storage microgrid. The algorithm incorporates adaptive inertia weight adjustment to improve convergence performance. The optimization model considers multiple generation sources including wind turbines, photovoltaic panels, battery energy storage, and micro gas turbines, while minimizing total operational costs under various technical constraints.