? ? ? ? ?本程序實現了風光柴儲微網中的粒子群優化（PSO）算法，用于優化微網的能源調度問題。具體來說，程序考慮了光伏發電、風力發電、柴油機發電（柴儲），并使用粒子群算法來優化這些能源的調度，以滿足負載需求。以下是對代碼的詳細說明：
1.?PSO算法參數設置
- 最大迭代次數?
maxgen：設定粒子群算法的最大迭代次數為2000次。 - 種群規模?
N：設定粒子群中的粒子數量為100。 - 粒子速度限制?
v_max,?v_min：粒子速度的上限和下限，分別為2和-2。 - 慣性權重?
w_max,?w_min：控制粒子群在搜索空間中探索的幅度。 - 學習因子?
c1,?c2：分別表示個體認知學習因子和社會學習因子。
2.?初始化
- 負載?
P_load：表示微網的負荷需求。 - 光伏發電?
solar：光伏電池板的發電量，按時間分布。 - 風力發電?
pv：風力發電機組的發電量，按時間分布。 - 發電總量?
P_pv：光伏和風力發電總和。
? ? ?粒子群的初始化過程中，粒子位置 x 和速度 v 被隨機生成。粒子的位置對應了微網中不同能源的發電量（包括風力、光伏、柴油發電機組等），這些位置的初始值是根據各能源的最大發電能力進行約束的。
3.?適應度計算
每個粒子的適應度是通過調用 fitness 函數來計算的。適應度函數衡量了粒子方案的優劣，旨在通過最小化成本、最大化系統效率或實現其他優化目標。
4.?主循環
在每次迭代中，粒子的位置和速度會根據以下公式進行更新：
- 速度更新：結合了慣性項、個體認知項和社會認知項。
- 位置更新：通過粒子的當前速度更新位置。
每個粒子的邊界都會根據各個能源的最大值和最小值進行限制，并采用反轉速度的方式進行邊界處理，以確保粒子不會越過邊界。
5.?粒子群更新
每個粒子在更新其位置和速度后，會重新評估其適應度。如果當前粒子的適應度比之前的最優適應度更好，則更新該粒子的最佳位置。全局最優解（即整個群體中最優的粒子）會在每次迭代后更新。
6.?結果可視化
程序會繪制多個圖表展示優化結果：
- 圖1：展示已消納的功率與總光伏發電功率的關系。
-
- 圖2：展示已消納的功率與總風電發電功率的關系。
-
- 圖3：展示光伏發電的已消納功率。
-
- 圖4：展示風電發電的已消納功率。
這些圖表幫助分析不同能源的調度效果，以及微網在負荷需求下的能源使用情況。
7 主程序代碼
clc;
clear;
%rng default;
%pso參數設置
format long;
maxgen=2000;%迭代次數
N=100;%種群規模
%c1=2;%自我學習因子
%c2=2;%群體學習因子
v_max=2; v_min=-2;%個體速度
w_max=0.9;w_min=0.1;
dmaxp=10;
dminp=0;
wminp=-10;
wmaxp=10;
s=1;
P_load=[52 50 50 51 56 63 70 75 76 80 78 74 72 72 76 80 85 88 90 87 78 71 65 56 ];%負載
%P_load=[152 150 150 151 156 163 170 175 176 180 178 174 172 172 176 180 185 188 190 187 178 171 165 156 ];
solar=[0 0 0 0 1 3 5 10 15 16 17 20 20 20 18 15 15 10 10 5 0 0 0 0 ]%光伏
pv=[75 65 60 55 50 57 60 66 50 58 60 60 57 50 61 58 55 65 56 60 72 66 75 76 ];%風力
P_pv=pv+solar;
%初始化種群
for i=1:N;for j=1:96;%4個維度v(i,j)=0.0;if j<=24;x(i,j)=rand()*pv(j);elseif j>24&&j<=48;x(i,j)=rand()*solar(j-24);elseif j>48&&j<=72;x(i,j)=dminp+rand()*(dmaxp-dminp);elseif j>72&&j<=96; x(i,j)=wminp+rand()*(wmaxp-wminp);endend
end%計算各個粒子的適應度，并初始化Pi和Pg
for i=1:Np(i)=fitness(x(i,:),s);y(i,:)=x(i,:);%每個粒子的個體尋優值
end
Pbest=fitness(x(1,:),s);
pg=x(1,:);%Pg為全局最優
for i=2:Nif fitness(x(i,:),s)<fitness(pg,s)Pbest=fitness(x(i,:),s);pg=x(i,:);%全局最優更新end
end
%進入主循環
for t=1:maxgenfor i=1:Nw=w_max-(w_max-w_min)*t/maxgen;%慣性權重更新c1=(0.5-2.5)*t/maxgen+2.5; %認知c2=(2.5-0.5)*t/maxgen+0.5; %社會認識
v(i,:)=w*v(i,:)+c1*rand()*(y(i,:)-x(i,:))+c2*rand()*(pg-x(i,:));
for m=1:96if(v(i,m)>v_max)v(i,m)=v_max;elseif(v(i,m)<-v_max)v(i,m)=-v_max;end
endx(i,:)=x(i,:)+v(i,:);%對粒子邊界處理for n=1:96if n<25if x(i,n)<0x(i,n)=0;v(i,n)=-v(i,n); elseif x(i,n)>pv(n)x(i,n)=pv(n);v(i,n)=-v(i,n); %反轉運動方向endelseif n>24&&n<49if x(i,n)<0x(i,n)=0;v(i,n)=-v(i,n); elseif x(i,n)>solar(n-24)x(i,n)=solar(n-24);v(i,n)=-v(i,n); endelseif n>49&&n<73if x(i,n)<dminpx(i,n)=dminp;v(i,n)=-v(i,n); elseif x(i,n)>dmaxpx(i,n)=dmaxp;v(i,n)=-v(i,n); endelse if x(i,n)<wminpx(i,n)=wminp;v(i,n)=-v(i,n); elseif x(i,n)>wmaxpx(i,n)=wmaxp;v(i,n)=-v(i,n); endendendend %對粒子進行評價，尋找最優值if fitness(x(i,:),t)<p(i)p(i)=fitness(x(i,:),t);y(i,:)=x(i,:);end if p(i)<PbestPbest=p(i);pg=y(i,:);s=t;end
endfor m=1:24pg1(m)=pg(m);
end
for m=25:48pg2(m-24)=pg(m);
end
for m=49:72pg3(m-48)=pg(m);
end
for m=73:96pg4(m-72)=pg(m);
end figure(1);
plot( pg1,'r-')
hold on
plot( pv,'-')
xlim([1 24])
grid
legend('已消納功率','總光伏功率');
xlabel('時間');
ylabel('功率');
title('風機發電出力')figure(2);
plot( pg2,'g-')
hold on
plot( solar,'-')
xlim([1 24])
grid
legend('已消納功率','總風電功率');
xlabel('時間');
ylabel('功率');
title('光伏發電出力')figure(3);
plot( pg3,'b-')
xlim([1 24])
grid figure(4);
plot( pg4,'k-')
xlim([1 24])
grid %plot (P_load,'-')
grid on
8 代碼獲取
