Диссертация (1144079), страница 14
Текст из файла (страница 14)
in[Bmin Bmax] !!!continueelsedelta=abs((tempo(1)-Data(k,1))/tempo(1))+abs((tempo(2)Data(k,2))/tempo(2));if delta <= gamma(n)h = h+1;Data(h,:) = tempo;% lay x1, x2, x3,... moi o phan tren de tinhtoan tiep theo% define variable power flows% define branch's voltagefor nn = 1:size(A,2)Data(h,nvar+nn)=(Ku(1,nn)+Ku(2,nn)*Data(h,1)+Ku(3,nn)*Data(h,2)+Ku(4,nn)*Data(h,1)*Data(h,2))/(1+Ki(3,Lo(ncase,1))*Data(h,1)+Ki(4,Lo(ncase,1))*Data(h,2)+Ki(5,Lo(ncase,1))*Data(h,1)*Data(h,2));endclear nn;% define branch's currentsfor nb=1:size(A,2)if nb == Lo(ncase,1)Data(h,nvar+size(A,2)+nb)=(Ki(1,nb)+Ki(2,nb)*Data(h,2))/(1+Ki(3,Lo(ncase,1))*Data(h,1)+Ki(4,Lo(ncase,1))*Data(h,2)+Ki(5,Lo(ncase,1))*Data(h,1)*Data(h,2));elseif nb == Lo(ncase,2)Data(h,nvar+size(A,2)+nb)=(Ki(1,nb)+Ki(2,nb)*Data(h,1))/(1+Ki(3,Lo(ncase,1))*Data(h,1)+Ki(4,Lo(ncase,1))*Data(h,2)+Ki(5,Lo(ncase,1))*Data(h,1)*Data(h,2));elseData(h,nvar+size(A,2)+nb)=(Ki(1,nb)+Ki(2,nb)*Data(h,1)+Ki(3,nb)*Data(h,2)+Ki(4,nb)*Data(h,1)*Data(h,2))/(1+Ki(3,Lo(ncase,1))*Data(h,1)+Ki(4,Lo(ncase,1))*Data(h,2)+Ki(5,Lo(ncase,1))*Data(h,1)*Data(h,2));end108end; clear nbData(h,nvar+2*size(A,2)+1:nvar+3*size(A,2))=Data(h,nvar+1:nvar+size(A,2)).*Data(h,nvar+size(A,2)+1:nvar+2*size(A,2));%PvarData(h,nvar+3*size(A,2)+1)=sum(Data(h,nvar+2*size(A,2)+1:nvar+3*size(A,2)));%sum PvarData(h,nvar+3*size(A,2)+2:nvar+3*size(A,2)+4)(S\[Data(h,16); Data(h,20); Data(h,24)]).';% I secquenceData(h,nvar+3*size(A,2)+5:nvar+3*size(A,2)+7)(S\[Data(h,17); Data(h,21); Data(h,25)]).';Data(h,nvar+3*size(A,2)+8:nvar+3*size(A,2)+10)(S\[Data(h,18); Data(h,22); Data(h,26)]).';Data(h,nvar+3*size(A,2)+11:nvar+3*size(A,2)+13)(S\[Data(h,19); Data(h,23); Data(h,27)]).';====endendend%Result(:,:,n,ncase,p) = Data;endData(:,nvar+3*size(A,2)+14)=(WeiCoe(WeC,1)*abs(Data(:,nvar+3*size(A,2)+8))).^2+(WeiCoe(WeC,1)*abs(Data(:,nvar+3*size(A,2)+9))).^2+(WeiCoe(WeC,1)*abs(Data(:,nvar+3*size(A,2)+10))).^2;Data = sortrows(Data,(nvar+3*size(A,2)+14));Result(ncase,:) = Data(1,:);% save all locationopti(p,:,ncase) = Data(1,:);clear Data;end% find the best locationIND=find(abs(Result(:,nvar+3*size(A,2)+14))<=min(abs(Result(:,nvar+3*size(A,2)+14))));optiLocal(p,:)= [Lo(IND,:) Result(IND,:)];clear Result IND;endfname = sprintf('ResultSec2%d%de%dWeC%dill.mat', Bus, Times,1000*eps, WeC);%clear fname;Прил.
9clear;clc;109format long e;tab = fopen('circuit.txt');TAB = fscanf(tab, '%d %d %d %d %d %d %d %d %d %d %d %d', [13 inf]);fclose(tab);A = TAB';NYG = 1; nBranch = (size(A,2)-NYG)/3; nNode = (size(A,1)-NYG)/3;S = [1 1 1; exp(-2*pi*i/3) exp(2*pi*i/3) 1; exp(2*pi*i/3) exp(-2*pi*i/3) 1];Fz = 18; kF =0.9;BUS = [1 2 4 5 6];KZ = [0.8 0.8 1 0.8 0.8 1; 0.7 0.7 1 0.7 0.7 1];Lim(1,1:6) = 210;%\Lim(2,1:6) = 225;%/Limit of V(t).Lim(3,1:6) = 0;%\Lim(4,1:6) = [1100 500 inf 1000 400 200]; %/Limit of I(t).Lim(5,1:6) = 0;%\Lim(6,1:6) = [2.5e5 1e5 inf 2e5 1e5 5e4]; %/Limit of P(t).Lim(7,1:6) = 0;%\Lim(8,1:6) = [800 380 inf 650 260 160]; %/Limit of I1(t).Lim(9,1:6) = 0;%\Lim(10,1:6) = [175 inf inf inf 55 60]; %/Limit of I2,I0(t).for nf = 1:length(BUS)Bus = BUS(nf)fname = sprintf('RealData2%d.txt',Bus);tab = fopen(fname);clear fname;colon = char(58);[FG,count1] = fscanf(tab,['%d' colon '%d' '%f' '%f' '%f' '%f' '%f' '%f' '%f''%f' '%f' '%f' '%f' '%f' '%f' '%f' '%f' '%f' '%f' '%f'],[20 inf]);fclose(tab);for p = 1:720%size(FG,2)ZA(1:nBranch)=[1e-9;0.001063+j*0.005905;KZ(1,Bus)*FG(18,p)*(FG(15,p)+j*sin(acos(FG(15,p)))); 1.31429-j*8.89675];ZB(1:nBranch)=[1e-9;0.001063+j*0.005905;KZ(2,Bus)*FG(19,p)*(FG(16,p)+j*sin(acos(FG(16,p)))); 1.31429-j*8.89675];ZC(1:nBranch)=[1e-9;0.001063+j*0.005905;FG(20,p)*(FG(17,p)+j*sin(acos(FG(17,p)))); 1.31429-j*8.89675];if NYG == 0ZN = [];elseZN(1:NYG) = 0.0017+j*0.0089;% if the circuit contain Zn!!!110endZ=[ZA ZB ZC ZN];Y=1./Z;clear ZA ZB ZC ZN;Y = diag(Y);E = zeros(size(A,2),1);E(1,1) = 220;%FG(3,p);E(5,1) = 220*exp(-2*pi*i/3);%FG(4,p)*exp(-2*pi*i/3);E(9,1) = 220*exp(2*pi*i/3);%FG(5,p)*exp(2*pi*i/3);%E(1,1) = FG(3,p)*(1+abs(Z(2))/abs(Z(3)));%E(5,1) = FG(4,p)*(1+abs(Z(6))/abs(Z(7)))*exp(-2*pi*i/3);%E(9,1) = FG(5,p)*(1+abs(Z(10))/abs(Z(11)))*exp(2*pi*i/3);Z1 = Z;Y1 = diag(1./Z1);% Solve circuitU0 = -(A*Y1*A')\(A*Y1*E);Us = A'*U0;U = Us+E;I = Y1*U;Data(1:size(A,2),p) = U;Data(1+size(A,2):2*size(A,2),p) = I;Data(1+2*size(A,2):3*size(A,2),p) = U.*I;Data(3*size(A,2)+1:3*size(A,2)+3,p) = (S\[Data(14,p);Data(22,p)])';% I secquenceData(3*size(A,2)+4:3*size(A,2)+6,p) = (S\[Data(15,p);Data(23,p)])';Data(3*size(A,2)+7:3*size(A,2)+9,p) = (S\[Data(16,p);Data(24,p)])';Data(3*size(A,2)+10:3*size(A,2)+12,p) = (S\[Data(17,p);Data(25,p)])';endfigure;subplot(3,1,1);plot(abs(Data(3,:)));hold on;grid on;ylim([Lim(1,Bus) Lim(2,Bus)]);plot(abs(Data(7,:)));plot(abs(Data(11,:)));legend('Phase A','Phase B','Phase C');title('Voltages [V]');set(gca,'TitleFontWeight','normal');set(gca,'XTickLabel',[''] );xlabel('24 hours');set(gca,'FontSize', Fz );Data(18,p);Data(19,p);Data(20,p);Data(21,p);111set(gca,'TitleFontSizeMultiplier',kF);set(gca,'LabelFontSizeMultiplier',kF);subplot(3,1,2)plot(abs(Data(16,:)));hold on;grid on;ylim([Lim(3,Bus) Lim(4,Bus)]);plot(abs(Data(20,:)));plot(abs(Data(24,:)));%legend('Phase A','Phase B','Phase C');title('Currents [A]');set(gca,'TitleFontWeight','normal');set(gca,'XTickLabel',[''] );xlabel('24 hours');set(gca,'FontSize', Fz );set(gca,'TitleFontSizeMultiplier',kF);set(gca,'LabelFontSizeMultiplier',kF);subplot(3,1,3)P(:,1) = abs(Data(3,:)).*abs(Data(16,:));%.*FG(15,:); Apparent powerP(:,2) = abs(Data(7,:)).*abs(Data(20,:));%.*FG(16,:);P(:,3) = abs(Data(11,:)).*abs(Data(24,:));%.*FG(17,:);plot(P(:,1));hold on;grid on;ylim([Lim(5,Bus) Lim(6,Bus)]);plot(P(:,2));plot(P(:,3));%legend('Phase A','Phase B','Phase C');title('Apparent powers [V.A]');set(gca,'TitleFontWeight','normal');set(gca,'XTickLabel',[''] );xlabel('24 hours');set(gca,'FontSize', Fz );set(gca,'TitleFontSizeMultiplier',kF);set(gca,'LabelFontSizeMultiplier',kF);fname = sprintf('T%d without compensation device',BUS(nf));print(fname, '-dpng', '-r300'); clear fname;figureyyaxis leftplot(abs(Data(46,:)));hold on;grid on;ylim([Lim(7,Bus) Lim(8,Bus)]);yyaxis rightplot(abs(Data(47,:)));plot(abs(Data(48,:)));ylim([Lim(9,Bus) Lim(10,Bus)]);112legend('I1','I2','I0');title('Sequence components of Currents [A]');set(gca,'TitleFontWeight','normal');set(gca,'XTickLabel',[''] );xlabel('24 hours');set(gca,'FontSize', Fz );set(gca,'TitleFontSizeMultiplier',kF);set(gca,'LabelFontSizeMultiplier',kF);fname = sprintf('T%d I120 without compensation device',Bus);print(fname, '-dpng', '-r300'); clear fname;meanI0I1 = mean(abs(Data(48,:))./abs(Data(46,:)))endПрил.
10clear;clc;Fz = 18; kF =0.9;EPS = [0 0.005 0.01 0.025];BUS = [1 2 4 5 6];TIMES = [4 3 3 4 3 3;];%Good resultLX = [3 1 1 3 1 1;];%Optimal LocalWEICOE = [20];% weight coefficient = (1 1 sqrt(2))for nWeC = 1:length(WEICOE)WeC = WEICOE(nWeC);for nBus = 3%1:length(BUS)Bus = BUS(nBus); % 1 2 4 5 6;Times = TIMES(nWeC,Bus);lx = LX(nWeC, Bus);Lim(1,1:6) = 210;%\Lim(2,1:6) = 225;%/Limit of V(t).Lim(3,1:6) = 0;%\Lim(4,1:6) = [1100 500 inf 1000 400 200]; %/Limit of I(t).Lim(5,1:6) = 0;%\Lim(6,1:6) = [2.5e5 1e5 inf 2e5 1e5 5e4]; %/Limit of P(t).Lim(7,1:6) = 0;%\Lim(8,1:6) = [800 380 inf 650 260 160]; %/Limit of I1(t).Lim(9,1:6) = 0;%\Lim(10,1:6) = [175 inf inf inf 55 60]; %/Limit of I2,I0(t).fname = sprintf('ResultSec2%d%de05k1010%di.mat',Bus,Times,WeC)load(fname);clear fname;mopti = abs(opti(:,:,lx));%figure113%subplot(3,1,1)%plot(mopti(:,5));%hold on;grid on;%plot(mopti(:,9));%plot(mopti(:,13));%subplot(3,1,2)%plot(mopti(:,18));%hold on;grid on;%plot(mopti(:,22));%plot(mopti(:,26));%subplot(3,1,3)%hold on;grid on;%plot(mopti(:,42));G=mopti;period = 5;for m=1%:15for n = 1:(size(G,1)-period)G(n,1:size(G,2)) = mean(G(n:n+period,1:size(G,2)));end;for m = 1:periodG(n+m,1:size(G,2)) = G(n,1:size(G,2));endendfigure;subplot(3,1,1)plot(abs(G(1:720,5)));hold on;grid on;ylim([Lim(1,Bus) Lim(2,Bus)]);plot(abs(G(1:720,9)));plot(abs(G(1:720,13)));legend('Phase A','Phase B','Phase C');title('Voltages [V]');set(gca,'TitleFontWeight','normal');set(gca,'XTickLabel',[''] );xlabel('24 hours');set(gca,'FontSize', Fz );set(gca,'TitleFontSizeMultiplier',kF);set(gca,'LabelFontSizeMultiplier',kF);subplot(3,1,2)plot(abs(G(1:720,18)));hold on;grid on;ylim([Lim(3,Bus) Lim(4,Bus)]);plot(abs(G(1:720,22)));plot(abs(G(1:720,26)));114%legend('Phase A','Phase B','Phase C');title('Currents [A]');set(gca,'TitleFontWeight','normal');set(gca,'XTickLabel',[''] );xlabel('24 hours');set(gca,'FontSize', Fz );set(gca,'TitleFontSizeMultiplier',kF);set(gca,'LabelFontSizeMultiplier',kF);subplot(3,1,3)plot(G(1:720,31));hold on;grid on;ylim([Lim(5,Bus) Lim(6,Bus)]);plot(G(1:720,35));plot(G(1:720,39));%legend('Phase A','Phase B','Phase C');title('Apparent powers [V.A]');set(gca,'TitleFontWeight','normal');set(gca,'XTickLabel',[''] );xlabel('24 hours');set(gca,'FontSize', Fz );set(gca,'TitleFontSizeMultiplier',kF);set(gca,'LabelFontSizeMultiplier',kF);fname = sprintf('T%d 1010%d',Bus,WeC);print(fname, '-dpng', '-r300'); clear fname;%print('T4 101020', '-dpng', '-r300')%plot(G(:,44));figure;yyaxis leftplot(G(1:720,49));hold on;grid on;ylim([Lim(7,Bus) Lim(8,Bus)]);yyaxis rightylim([Lim(9,Bus) Lim(10,Bus)]);plot(G(1:720,50));plot(G(1:720,51));legend('I1','I2','I0');title('Sequence components of Currents [A]');set(gca,'TitleFontWeight','normal');set(gca,'XTickLabel',[''] );xlabel('24 hours');set(gca,'FontSize', Fz );set(gca,'TitleFontSizeMultiplier',kF);set(gca,'LabelFontSizeMultiplier',kF);fname = sprintf('T%d I120 1010%d',Bus,WeC);print(fname, '-dpng', '-r300'); clear fname;115figure;plot(G(1:720,1));ylim([0 inf]);hold on;grid on;plot(G(1:720,2));legend('X1','X2');title('Reactance values [Ohm]');set(gca,'TitleFontWeight','normal');set(gca,'XTickLabel',[''] );xlabel('24 hours');set(gca,'FontSize', Fz );set(gca,'TitleFontSizeMultiplier',kF);set(gca,'LabelFontSizeMultiplier',kF);fname = sprintf('T%d value of X 1010%d',Bus,WeC)print(fname, '-dpng', '-r300'); clear fname;meanI0I1 = mean(G(1:720,51)./G(1:720,49))meanI0I12 = mean(mopti(1:720,51)./mopti(1:720,49))endendПрил.
11clear; clc;format long e;EPS = [0 0.005];BUS = [1 2 3 4 5];XMAX1 = [20 20 20 20 20];XMAX2 = [20 20 20 20 20];XMIN1 = [0 0 0 0 0];XMIN2 = [0 0 0 0 0];TIMES = 1;KZ = [0.8 0.8 0.8 0.8 1; 0.7 0.7 0.7 0.7 1];for nBus = 1:length(BUS)Bus = BUS(nBus ); % 1 2 4 5 6;for nTimes = 1:TIMESTimes = nTimes;for ie = 1:length(EPS)SymPar;end;end;end;Прил. 12format long e;profile on;tab = fopen('circuit.txt');TAB = fscanf(tab, '%d %d %d %d %d %d %d %d %d %d %d %d', [13 inf]);116fclose(tab);A = TAB';NYG = 1; nBranch = (size(A,2)-NYG)/3; nNode = (size(A,1)-NYG)/3;tab = fopen('local.txt');Lo = fscanf(tab, '%d %d', [2 inf]);fclose(tab);Lo=Lo';Xmax=[XMAX1(Bus) XMAX2(Bus)];% limit of reactive conductivityXmin=[XMIN1(Bus) XMIN2(Bus)];eps = EPS(ie);% loadfname = sprintf('RealData%d.txt', Bus);tab = fopen(fname);clear fname;colon = char(58);[FG,count1] = fscanf(tab,['%d' colon '%d' '%f' '%f' '%f' '%f' '%f' '%f' '%f' '%f''%f' '%f' '%f' '%f' '%f' '%f' '%f' '%f' '%f' '%f'],[20 inf]);fclose(tab);for p = 1:size(FG,2)%ZA(1:nBranch)=[1e-9;0.001063+j*0.005905;KZ(1,Bus)*FG(18,p)*(FG(15,p)-j*sin(acos(FG(15,p)))); 1.31429-j*8.89675];%ZB(1:nBranch)=[1e-9;0.001063+j*0.005905;KZ(2,Bus)*FG(19,p)*(FG(16,p)-j*sin(acos(FG(16,p)))); 1.31429-j*8.89675];%ZC(1:nBranch) = [1e-9; 0.001063+j*0.005905; FG(20,p)*(FG(17,p)j*sin(acos(FG(17,p)))); 1.31429-j*8.89675];ZA(1:nBranch)=[1e-9;0.001063+j*0.005905;KZ(1,Bus)*FG(18,p)*FG(15,p)-KZ(1,Bus)*FG(18,p)*j*sin(acos(FG(15,p)));1.31429-j*8.89675];ZB(1:nBranch)=[1e-9;0.001063+j*0.005905;KZ(1,Bus)*FG(18,p)*FG(15,p)-KZ(2,Bus)*FG(19,p)*j*sin(acos(FG(16,p)));1.31429-j*8.89675];ZC(1:nBranch)=[1e-9;0.001063+j*0.005905;KZ(1,Bus)*FG(18,p)*FG(15,p)-FG(20,p)*j*sin(acos(FG(17,p)));1.31429j*8.89675];if NYG == 0ZN = [];elseZN(1:NYG) = 0.0017+j*0.0089;% if the circuit contain Zn!!!endZ=[ZA ZB ZC ZN];Y=1./Z;clear ZA ZB ZC ZN;117Y = diag(Y);E = zeros(size(A,2),1);E(1,1) = 220;E(5,1) = 220*exp(-2*pi*i/3);E(9,1) = 220*exp(2*pi*i/3);ncase = 1;nvar = 2;for ncase = 1:size(Lo,1)X1 = [1 1 2 5 7 9 11 15];X2 = [10 7 4 1 1 15 13 9];for nex = 1:size(X1,2)Z1 = Z;Z1(Lo(ncase,1)) = (Z1(Lo(ncase,1))*j*X1(nex))/(Z1(Lo(ncase,1)) +j*X1(nex));Z1(Lo(ncase,2)) = (Z1(Lo(ncase,2))*j*X2(nex))/(Z1(Lo(ncase,2)) +j*X2(nex));Y1 = diag(1./Z1);% Solve circuitU0 = -(A*Y1*A')\(A*Y1*E);Us = A'*U0;U = Us+E;I = Y1*U;Uex(:,nex) = U.*(1 + eps*2*(rand(size(A,2),1) - 0.5));Iex(:,nex) = I.*(1 + eps*2*(rand(size(A,2),1) - 0.5));endfor nf = 1:size(X1,2)Vx(nf,:)=[1X2(nf)-Iex(Lo(ncase,1),nf)*X1(nf)Iex(Lo(ncase,1),nf)*X2(nf) -Iex(Lo(ncase,1),nf)*X1(nf)*X2(nf)];Tx(nf,1) = [Iex(Lo(ncase,1),nf)];end;Ki(:,Lo(ncase,1)) = (Vx'*Vx)\(Vx'*Tx);clear Vx Tx nf;for nf = 1:3Vy(nf,:) = [1 X1(nf)];Ty(nf,1)=[(1+Ki(3,Lo(ncase,1))*X1(nf)+Ki(4,Lo(ncase,1))*X2(nf)+Ki(5,Lo(ncase,1))*X1(nf)*X2(nf))*Iex(Lo(ncase,2),nf)];end;Ki(1:2,Lo(ncase,2)) = (Vy'*Vy)\(Vy'*Ty);clear Vy Ty nf;% tim he so cua bieu thuc dong dien cac nhanh con laifor ni = 1:size(A,2)if ni == Lo(ncase,1) || ni == Lo(ncase,2)118continueelsefor nf = 1:6V(nf,:) = [1 X1(nf) X2(nf) X1(nf)*X2(nf)];T(nf,1)=[(1+Ki(3,Lo(ncase,1))*X1(nf)+Ki(4,Lo(ncase,1))*X2(nf)+Ki(5,Lo(ncase,1))*X1(nf)*X2(nf))*Iex(ni,nf)];end;Ki(1:4,ni) = (V'*V)\(V'*T);clear nf V T;endend; clear ni% tim he so cua bieu thuc dien ap tren cac phan tu tieu thufor nu = 1:size(A,2)for nf = 1:6V(nf,:) = [1 X1(nf) X2(nf) X1(nf)*X2(nf)];T(nf,1)=[(1+Ki(3,Lo(ncase,1))*X1(nf)+Ki(4,Lo(ncase,1))*X2(nf)+Ki(5,Lo(ncase,1))*X1(nf)*X2(nf))*Uex(nu,nf)];end;Ku(1:4,nu) = (V'*V)\(V'*T);clear nf V T;end; clear nu;% Mathematical Optimization% First step, we choose 'n_row' random value of com.dev.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
