Matlab电动汽车选址

clear; clc; close all; %% 基础数据 % 1 2 3 4 5 % b=[需求坐标，需求负荷，成本系数，区域] % 充电需求点坐标 b(:,1) = 1.0e+003 .* ([0.1092, 0.1197, 0.2578, 0.4259, 0.1257, 0.2803, 0.4439, 0.5505, 0.5610, 0.5700, 0.1332, 0.3013, 0.4559, 0.5850, 0.1452, 0.3163, 0.4739, 0.5880, 0.3193, 0.4784, 0.6015, 0.6736, 0.8657, 1.0308, 0.6781, 0.8327, 1.0188, 0.6811, 0.8192, 1.0083, 0.6691, 0.6976, 0.8191, 1.0098]); b(:,2) = 1.0e+003 .* [0.1348; 0.2399; 0.1724; 0.1739; 0.3420; 0.3375; 0.3360; 0.1108; 0.2024; 0.3075; 0.4591; 0.4455; 0.4380; 0.4260; 0.6092; 0.5341; 0.5341; 0.5341; 0.6407; 0.6452; 0.6467; 0.1574; 0.1649; 0.1634; 0.2924; 0.2909; 0.2939; 0.4425; 04576; 0.4546; 0.5431; 0.6392; 0.6377; 0.6347]; % 充电需求点常规电力负荷点负荷b(:,3)（kW） b(:,3) = [2480; 2480; 8680; 11400; 890; 2340; 4160; 560; 1670; 5010; 2670; 8280; 7400; 1430; 7500; 4840; 3400; 4290; 3840; 3680; 2560; 7000; 14800; 8960; 3160; 7000; 5000; 2280; 10360; 10000; 760; 6000; 7040; 5600]; % 集中充电站坐标,两个，各区域一个 bcs = [937.7296 379.5010; 310.3141 238.4076]; % 选址数量 Tn = 7; %% 成本参数 na = 4500; alp = 0.1; b(:,4) = round(alp .* b(:,3) ./ sum(b(:,3)) .* na); b(23,4) = 37; ns = 4; mui = 0.6; Nchz = round(mui .* sum(b(:,4)) ./ ns); bm = 1.0e+003 .* [0.0086, 0.0088; 1.1734, 0.0088; 1.1734, 0.7412; 0.0086, 0.7412; 0.0086, 0.0088]; BL = sqrt(8.2 * 1.0e6 ./ ((max(bm(:,1)) - min(bm(:,1))) * (max(bm(:,2)) - min(bm(:,2))))); %% 区域划分 Area1 = 1.0e+003 .* [0.9377 -1.0860; 1.1734 0.0088; 1.1734 0.7412; 0.3103 1.7040; 0.9377 -1.0860]; Area1 = [Area1, 1 .* ones(size(Area1,1),1)]; Area2 = 1.0e+003 .* [0.0086 0.0088; 0.9377 -1.0860; 0.3103 1.7040; 0.0086 0.7412; 0.0086 0.0088]; Area2 = [Area2, 2 .* ones(size(Area2,1),1)]; % 区域分界 vv = [Area1; Area2]; for k = 1:size(bcs,1) Ai = find(vv(:,3) == k); xx = vv(Ai,1); yy = vv(Ai,2); kk = convhull(xx,yy); % 判断点是否在多边形内 in = inpolygon(b(:,1), b(:,2), xx(kk), yy(kk)); b(in,5) = k; end % 区域归属 Ep = []; for i = 1:size(bcs,1) gb = b(b(:,5) == i, :); Ep = [Ep; [sum(gb(:,4)), round(mui .* sum(gb(:,4)) ./ ns), i]]; end %% 粒子群算法参数 PopSize = 30; % 种群数量（优化参数） MaxIter = 500; % 最大迭代次数（优化参数） % 学习因子 c1s = 2.5; c2s = 0.5; c1e = 0.5; c2e = 2.5; w_start = 0.9; w_end = 0.4; Iter = 1; % 上下限约束 xmax = max(Area1(:,1)); xmin = min(Area1(:,1)); ymax = max(Area1(:,2)); ymin = min(Area1(:,2)); x = xmin + (xmax - xmin) .* rand(Tn, PopSize); y = ymin + (ymax - ymin) .* rand(Tn, PopSize); %% 初始化 X = [x; y]; V = rand(Tn * 2, PopSize); Vmax = 0.1 * max((xmax - xmin), (ymax - ymin)); inAr1 = find(b(:,5) == 1); bb = [b(inAr1,1:2), b(inAr1,4)]; % 初始寻优 for pk = 1:PopSize [FX(pk),~,~,~,~,~,~,~,~,~] = VorCostCDEV(X(1:Tn,pk), X(Tn+1:end,pk), bb, bcs(1,:), BL); end % 输出初始最优 PBest = X; FPBest = FX; [Fgbest, r] = min(FX); Fm(Iter) = Fgbest; CF = Fgbest; Best = X(:,r); FBest(Iter) = Fgbest; FgNum = 0; %% 粒子群寻优 while (Iter <= MaxIter) Iter = Iter + 1; w_now = ((w_start - w_end) * (MaxIter - Iter) / MaxIter) + w_end; A = repmat(X(:,r), 1, PopSize); R1 = rand(Tn * 2, PopSize); R2 = rand(Tn * 2, PopSize); c1 = c1e + (c1s - c1e) * (1 - acos(-2 * Iter / (MaxIter + 1) + 1) / pi); c2 = c2e + (c2s - c2e) * (1 - acos(-2 * Iter / (MaxIter + 1) + 1) / pi); % 更新速度和位置 V = w_now * V + c1 * R1 .* (PBest - X) + c2 * R2 .* (A - X); changeRows = V > Vmax; V(changeRows) = Vmax; changeRows = V < -Vmax; V(changeRows) = -Vmax; X = X + 1.0 * V; % 计算适应度 for pk = 1:PopSize [FX(pk),~,~,~,~,~,~,~,~,~] = VorCostCDEV(X(1:Tn,pk), X(Tn+1:end,pk), bb, bcs(1,:), BL); end % 更新局部最优解 P = FX < FPBest; FPBest(P) = FX(P); PBest(:,P) = X(:,P); % 更新全局最优解 [Fgbest, r] = min(FPBest); Fm(Iter) = Fgbest; if Fgbest < CF [FBest, gr] = min(FPBest); Best = PBest(:,gr); CF = Fgbest; FgNum = 0; else FgNum = FgNum + 1; end % 全局最优扰动策略 if mod(Iter, 100) == 0 % 每100次迭代进行一次扰动 Best = Best + 0.01 * randn(size(Best)); end % 粒子重置策略 if FgNum > 20 % 连续20次未更新则重置 SubX = r; while SubX == r || SubX == 0 SubX = round(rand * (PopSize - 1)) + 1; end r = SubX; FgNum = 0; end end %% 输出最优值 % 坐标 x = Best(1:Tn); y = Best(Tn+1:end); [f1, f2, f3, f4, Num_chI, num_Epc, f5, f6, f7, f8] = VorCostCDEV(x, y, bb, bcs(1,:), BL); FBest Best %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% 作图 %% 迭代图 figure(1); Iter_t = 1:1:MaxIter + 1; plot(Iter_t, Fm, '.-'); legend('每次迭代值'); title('迭代图'); %% 选址图 figure(2); hold on; [vxT, vyT] = VoronoiT(bcs(:,1), bcs(:,2), 0); % 集中式充电站 plot(bcs(:,1), bcs(:,2), 'gs', 'linewidth', 12); % 区域分界 plot(vxT, vyT, 'r -', 'linewidth', 3); % 需求点 plot(b(:,1), b(:,2), 'k*', 'linewidth', 5); % 周边 plot(bm(:,1), bm(:,2), 'k-', 'linewidth', 3); % 需求点排序 for k = 1:length(b(:,4)) str = num2str(b(k,4)); text(b(k,1) - 15, b(k,2) + 25, str, 'FontSize', 23, 'color', 'black'); end axis equal; [vx, vy] = voronoi(x, y); plot(x, y, 'b^', 'linewidth', 3); % 最优充电站排序 for k = 1:length(x) str = num2str(k); text(x(k), y(k), str, 'FontSize', 20, 'color', 'red'); end axis([min(bm(:,1)) - 3 max(bm(:,1)) + 3 min(bm(:,2)) - 3 max(bm(:,2)) + 3]); legend('集中充电站', '区域分界', '充电需求点', '周边界限', '充电站最优选址区'); title('规划选址图'); %% 选址图 figure(3); [vxT, vyT] = VoronoiT(bcs(:,1), bcs(:,2), 0); % 集中式充电站 plot(bcs(:,1), bcs(:,2), 'gs', 'linewidth', 12); % 区域分界 plot(vxT, vyT, 'r -', 'linewidth', 3); % 需求点 plot(b(:,1), b(:,2), 'k*', 'linewidth', 5); % 周边 plot(bm(:,1), bm(:,2), 'k-', 'linewidth', 3); % 需求点排序 for k = 1:length(b(:,4)) str = num2str(b(k,4)); text(b(k,1) - 15, b(k,2) + 25, str, 'FontSize', 23, 'color', 'black'); end