目錄

?

１配電網故障定位原理

1.1?故障信息的數學表示

1.2 二進制PSO算法的配電網故障定位原理?

1.3?構造配電網故障定位評價函數?

1.4 求期望函數

2 二進制粒子群算法的簡單舉例------背包問題

2.1?算例

2.2?Python?實現

2.3 Matlab?實現

3??二進制粒子群算法的配電網故障定位(Matlab實現)

---

１配電網故障定位原理

1.1?故障信息的數學表示

在下圖中K表示斷路器，每一個斷路器上均有一個FTU裝置，可以反饋斷路器開關是否過流，用表示上傳的故障信息，反映的是各分段開關處是否流過故障電流有故障電流為1,否則為0）。即：

　　　　　　　?

因為FTU上傳的信息可分為有故障信息及無故障信息兩類,對于分段區間來講也只能是有故障及無故障兩種情況,所以我們可以用二進制編碼規則對配電網故障定位問題進行數學建模。以上圖所示輻射狀配電網為例,系統擁有12個分段開關,我們可以用一串12位的二進制代碼表示FTU的上傳信息,作為程序的輸入,1代表對應的開關有過流信息,0代表對應的開關無過流信息。同時用另一串12位的二進制代碼作為程序的輸出,代表對應饋線區間發生故障,代表無故障。

1.2 二進制PSO算法的配電網故障定位原理?

應用二進制PSO算法求解配電網故障定位問題,粒子的位置代表配電網中饋線區段的狀態,粒子的維數代表配電網的饋線區段總數。每一饋線區段存在0和1兩種狀態,0表示正常狀態,1表示故障狀態,饋線區段的狀態為待求量。因此,N段饋線區段的狀態求解就轉化成N維粒子群優化求解,每個粒子的N維位置都表示為配電網N段饋線區段的潛在狀態。每次迭代過程中,通過評價函數評價各粒子位置優劣,更新粒子的當前最優位置和全體粒子的最優位置,進而更新粒子的速度和位置,直到滿足程序終止條件為止。最終得出的粒子群的全局最優位置就是所求的各饋線區段的實際狀態。

離散粒子群
? ? ? 基本的粒子群算法是在連續域中搜索函數極值的有力工具，在此之后，提出了離散二進制粒子群，將離散問題空間映射到連續粒子運動空間，并適當修改粒子群算法來求解，在計算上仍保留經典粒子群算法速度-位置更新運算法則，粒子在狀態空間的取值和變化只限于0和1兩個值，而速度的每一維Vij代表位置每一位xij的取值為1的可能性，，因此在連續粒子群中的vij更新公式仍保持不變，但是Pbest和Gbest只在【0,1】的取值，其位置更新等式為：

(式中：r是0-1之間的隨機數。)

? ? ? ? ? ? ? ? ? ? ? ???

1.3?構造配電網故障定位評價函數?

基于待求各饋線區段實際狀態下所對應的信息應與實際上傳故障信息偏差最小的原則，構造如下評價函數：

? ? ? ? ? ? ? ? ? ? ? ? ??

表達式的值為每個潛在解對應的適應度值,值越小表示解越優良,因此評價函數應取極小值。

式中：為第j個開關FTU上傳的故障信息,取值為1認為該開關流過了故障電流,為0則未流過；為各開關節點的期望狀態,若該開關流過了故障電流,則期望狀態為1,相反,期望狀態為0,表達式為各區段狀態的函數。N為配電網中饋線區段的總數，為配電網中各設備狀態,為1表明設備故障,取0則設備正常。表示一個權系數乘以故障設備數,。是根據故障診斷理論中“最小集”概念設置的權重系數,取值范圍介于與1之間,表明故障區間數越少解越優,避免出現誤診斷,本文中權系數取0.5。

1.4 求期望函數

針對不同的網絡，期望函數的表達式不同，無法用統一的格式表示出來。僅以上述算例為例來說明如何求取期望函數。

在輻射型配電網中,任何區間發生故障,必將導致其上游開關流過故障過電流。這是求取期望函數的原則。假設圖一的饋線區段6處發生故障,那么分段開關K1,K2,K3 ,K6,都將有過電流若饋線區段11處發生故障,則K1,K2,K3 ,K6,K10都將有過電流。由此類推我們可得到各個分段開關的期望函數:?

? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ??

2 二進制粒子群算法的簡單舉例------背包問題

2.1?算例

背包問題。有N件物品和一個容量為V的背包。第i件物品的體積是c(i)，價值是w(i)。求解將哪些物品放入背包可使物品的體積總和不超過背包的容量，且價值總和最大。假設物品數量為10，背包的容量為300。每件物品的體積為[95，75，23， 73，50， 22，6，57，89， 98]， 價值為[89， 59，19， 43,100，72，44，16，7, 64]。

2.2?Python?實現

import numpy as np import matplotlib.pyplot as plt from tqdm import tqdm#進度條設置 import matplotlib.pyplot as plt import matplotlib; matplotlib.use('TkAgg') from pylab import * mpl.rcParams['font.sans-serif'] = ['SimHei'] # 指定默認字體 mpl.rcParams['axes.unicode_minus'] = False # 解決保存圖像是負號'-'顯示為方塊的問題class BackPack(object):def __init__(self):self.Weight = [95,75,23,73,50,22,6,57,89,98] # 物品體積self.Value = [89, 59, 19, 43, 100, 72, 44, 16, 7, 64] # 物品價值self.N = 100 # 群體粒子個數self.D = len(self.Weight) # 粒子維數self.T = 100 # 最大迭代次數self.c1 = 1.5 # 學習因子1self.c2 = 1.5 # 學習因子2self.w = 1 # 慣性因子，一般取1self.V_max = 10 # 速度最大值self.V_min = -10 # 速度最小值self.Weight_max = 300 # 背包容量self.afa = 10 # 懲罰系數self.G=100 #迭代次數#=====初始化種群=========def init_x(self):""":return: 隨機生成的種群(二維list)"""X = np.random.choice([0, 1], size=(self.N, self.D))return X#======初始化速度==========def init_v(self):""":return: 速度"""V = np.random.random(size=(self.N, self.D)) ##10維return V#=====適應度值========def fitness_func(self,X):""":param X: 種群:return:""""""計算粒子的的適應度值，也就是目標函數值，X 的維度是 size * 10 """size = np.shape(X)[0] #種群個數result = np.zeros((size, 1))for i in range(size): # 遍歷每一個粒子fit = sum(X[i] * self.Value) # 物品價值TotalSize = sum(X[i] * self.Weight) # 物品重量if TotalSize <= self.Weight_max: # 如果小于限制重量fit = fitresult[i] = fitelse:#fit = fit - self.afa * (TotalSize - self.Weight_max) # 對fit進行懲罰fit=0result[i] = fitreturn result # 我們要求result越大#=====速度更新公式========def velocity_update(self,V, X, pbest, gbest):"""根據速度更新公式更新每個粒子的速度種群size=100:param V: 粒子當前的速度矩陣，size*10 的矩陣:param X: 粒子當前的位置矩陣，size*10 的矩陣:param pbest: 每個粒子歷史最優位置，size*10 的矩陣:param gbest: 種群歷史最優位置，1*10 的矩陣"""r1 = np.random.random((self.N, 1)) #(粒子個數,1)r2 = np.random.random((self.N, 1))V = self.w * V + self.c1 * r1 * (pbest - X) + self.c2 * r2 * (gbest - X) # 直接對照公式寫就好了#===防止越界處理===V[V < self.V_min] = np.random.random() * (self.V_max - self.V_min) + self.V_minV[V > self.V_max] = np.random.random() * (self.V_max - self.V_min) +self.V_minreturn V#==========位置更新公式==========def position_update(self,X, V):"""根據公式更新粒子的位置:param X: 粒子當前的位置矩陣，維度是 size*10:param V: 粒子當前的速度舉著，維度是 size*10return 更新后的種群"""for i in range(self.N): # 遍歷每一個粒子# 修改速度為sigmoid形式V[i, :] = 1. / (1 + np.exp(-np.array(V[i, :])))for j in range(self.D): # 遍歷粒子中的每一個元素rand = np.random.random() # 生成 0-1之間的隨機數if V[i, j] > rand:X[i, j] = 1else:X[i, j] = 0#對當前個體進行限制，放在超重while np.sum(X[i]*self.Weight)>self.Weight_max:#如果當前粒子超重X[i]=np.random.choice([0, 1], size=(1, self.D))return Xdef update_pbest(self,X, fitness, pbest, pbestfitness, m):"""更新個體最優:param X: 當前種群:param fitness: 當前每個粒子的適應度:param pbest: 更新前的個體最優解:param pbestfitness: 更新前的個體最優適應度:param m: 粒子數量:return: 更新后的個體最優解、個體最優適應度"""for i in range(m):if fitness[i] > pbestfitness[i]:pbest[i] = X[i]pbestfitness[i] = fitness[i]return pbest, pbestfitnessdef update_gbest(self,pbest, pbestfitness, gbest, gbestfitness, m):"""更新全局最優解:param pbest: 粒子群:param pbestfitness: 個體適應度(個體最優):param gbest: 全局最優解:param gbestfitness: 全局最優解:param m: 粒子數量:return: gbest全局最優，g對應的全局最優解"""for i in range(m):if pbestfitness[i] > gbestfitness:gbest = pbest[i]gbestfitness = pbestfitness[i]return gbest, gbestfitnessdef main(self):fitneess_value_list = [] # 記錄每次迭代過程中的種群適應度值變化x = self.init_x() # 初始化xv = self.init_v() # 初始化v# 計算種群各個粒子的初始適應度值p_fitness = self.fitness_func(x)# 計算種群的初始最優適應度值g_fitness = p_fitness.max()# 講添加到記錄中fitneess_value_list.append(g_fitness)# 初始的個體最優位置和種群最優位置pbest = xgbest = x[p_fitness.argmax()] ##===接下來就是不斷迭代了========for i in tqdm(range(self.G)):pbest=pbest.copy()#必須加，不然會被篡改值，造成結果錯p_fitness= p_fitness.copy()#必須加，不然會被篡改值，造成結果錯gbest=gbest.copy()#必須加，不然會被篡改值，造成結果錯g_fitness=g_fitness.copy()#必須加，不然會被篡改值，造成結果錯v = self.velocity_update(v, x, pbest, gbest) # 更新速度x = self.position_update(x, v) # 更新位置p_fitness2 = self.fitness_func(x) # 計算子代各個粒子的適應度# 更新每個粒子的歷史最優位置pbest, p_fitness=self.update_pbest(x, p_fitness2, pbest, p_fitness, self.N)#更新群體的最優位置gbest, g_fitness=self.update_gbest(pbest, p_fitness, gbest, g_fitness, self.N)# 記錄最優迭代結果fitneess_value_list.append(g_fitness)print("最優適應度是：%.5f" % fitneess_value_list[-1])print("最優解是", gbest)print('最優解對應的體積',np.sum(gbest*self.Weight))print('最優解對應的價值',np.sum(gbest*self.Value))plt.plot(fitneess_value_list,label='迭代曲線')plt.xlabel('迭代次數')plt.ylabel('適應度')plt.legend()plt.show()if __name__=='__main__':pso = BackPack()pso.main()

2.3 Matlab?實現

%% 離散粒子群算法解決0-1背包問題 %% 初始化 clear all; %清除所有變量 close all; %清圖 clc; %清屏 N=100; %群體粒子個數 D=10; %粒子維數 c1=1.5; %學習因子1 c2=1.5; %學習因子2 Wmax=0.8; %慣性權重最大值 Wmin=0.4; %慣性權重最小值 Vmax=10; %速度最大值 Vmin=-10; %速度最小值 V = 300; %背包容量 C = [95,75,23,73,50,22,6,57,89,98]; %物品體積 W = [89,59,19,43,100,72,44,16,7,64]; %物品價值 afa = 2; %懲罰函數系數 %% 初始化種群個體（限定位置和速度） x=randi(N,D); %隨機獲得二進制編碼的初始種群，randint在后期版本中會被randi代替 v=rand(N,D) * (Vmax-Vmin)+Vmin; %% 初始化個體最優位置和最優值 p=x; pbest=ones(D,1); for i=1:Dpbest(i)= func(x(i,:),C,W,V,afa); end %% 初始化全局最優位置和最優值 g=ones(1,D); gbest=-10000; for i=1:Dif(pbest(i)>gbest)g=p(i,:);gbest=pbest(i);end end gb=ones(1,D); %% 按照公式依次迭代直到滿足精度或者迭代次數 while 1>0 %無限迭代for j=1:D%% 更新個體最優位置和最優值if (func(x(j,:),C,W,V,afa)>pbest(j))p(j,:)=x(j,:);pbest(j)=func(x(j,:),C,W,V,afa);end%% 更新全局最優位置和最優值if(pbest(j)>gbest)g=p(j,:);gbest=pbest(j);end%% 計算動態慣性權重值w=Wmax-(Wmax-Wmin)*i/D;%% 跟新位置和速度值v(j,:)=w*v(j,:)+c1*rand*(p(j,:)-x(j,:))...+c2*rand*(g-x(j,:));%% 邊界條件處理for ii=1:Dif (v(j,ii)>Vmax) | (v(j,ii)< Vmin)v(j,ii)=rand * (Vmax-Vmin)+Vmin;endendvx(j,:)=1./(1+exp(-v(j,:)));for jj=1:Dif vx(j,jj)>randx(j,jj)=1;elsex(j,jj)=0;endendend%% 記錄歷代全局最優值gb(i)=gbest;if g.*C<=300 %如果滿足條件 退出breakendg; %最優個體 figure plot(gb) xlabel('迭代次數'); ylabel('適應度值'); title('適應度進化曲線') end %% 適應度函數 function result = func(f,C,W,V,afa) fit = sum(f.*W); TotalSize = sum(f.*C); if TotalSize <= Vfit = fit; elsefit = fit - afa * (TotalSize - V); end result = fit;

3??二進制粒子群算法的配電網故障定位(Matlab實現)

clc;clear all;%% 初始化參數global Ij; %全局變量 global Vmax;%允許誤差MaxNum=120;%粒子最大迭代次數D=12;%粒子的維數代表配電網的饋線區段總數 Particlesize=30;%粒子群規模c1=2.1;%每個粒子的個體學習因子，也稱為加速常數c2=2.1;%每個粒子的社會學習因子，也稱為加速常數w=1;%慣性因子Vmax=4;%粒子的最大飛翔速度S_b=randi([0,1],Particlesize,D);%粒子所在的位置,randint在后期版本中會被randi代替V=Vmax*rand(Particlesize,D);%粒子的飛翔速度Ij=[1 1 1 1 0 1 1 1 0 0 1 1];%FTU 反饋的值這里就是輸入，請修改這里，看結果fitness=F(S_b');%這里注意下，要轉置%% inline定義的適應度函數會使程序運行速度大大降低for i=1:Particlesizef(i)=F(S_b(i,:));%S_b配電網中各設備狀態,為1表明設備故障,取0則設備正常。end%% 個體最優與總群最優personalbest_x=S_b;%這里是種群里的每個粒子都認為自己是最好的personalbest_faval=f;[globalbest_favali]=min(personalbest_faval);%全局最優globalbest_x=personalbest_x(i,:);k=1;Start=0;%% 這里就是把整個種群中的所有粒子都檢查了一遍，保證當前最優while k<=MaxNumRemember=globalbest_favali;for i=1:Particlesizef(i)=F(S_b(i,:));%計算適應度函數值if f(i)<personalbest_faval(i) %判斷當前位置是否是歷史上最佳位置personalbest_faval(i)=f(i);personalbest_x(i,:)=S_b(i,:);end end%% 更新種群最優[globalbest_favali]=min(personalbest_faval);globalbest_x=personalbest_x(i,:);%% 更新粒子群里每個個體的最新位置for i=1:Particlesize V(i,:)=w*V(i,:)+c1*rand*(personalbest_x(i,:)-S_b(i,:))...+c2*rand*(globalbest_x-S_b(i,:));%更新了速度for j=1:length(Ij)if rand(1)<Sigmoid(V(i,j))%這里體現了速度決定位置S_b(i,j)=1;elseS_b(i,j)=0;endendend%% 粒子的位置代表配電網中饋線區段的狀態globalbest_xif globalbest_favali-Remember==0Start=Start+1;if Start==20disp(globalbest_x);break;endendk=k+1; end

完整代碼：

?

總結

以上是生活随笔為你收集整理的改进粒子群算法的配电网故障定位(PythonMatlab代码实现)的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。