应用介绍

这是优化配电系统中太阳辐射年曲线的光伏渗透的matlab代码

clc;

clear all;

%bus[busno buscode voltage angle pd qd pg qg]

busdata = [ 1 1 1.06 0 0 0 0 0;2 0 1 0 0.5 0.2 0 0;3 0 1 0 0.4 0.3 0 0;4 0 1 0 0.3 0.1 0 0;]

%linedata = [startbus endbus r x shunt-y]

linedata = [1 2 2 8.0 0;1 3 1 4.0 0;2 3 0.666 2.664 0;2 4 1 4.0 0;3 4 2 8.0 0;]

       n1 =linedata(:,1);%linedata of startbus

       nr =linedata(:,2);%linedata of endbus

       nbr =length(n1);

       nbus =max(max(n1),max(nr));

       r =linedata(:,3);%line resistance

       x =linedata(:,4);%line reactance

       bc =linedata(:,5);

       y =complex(r,-x);

       epsilon =0.001;

       r =epsilon+1;

       epsilon;

       ybus =zeros(nbus,nbus);

       for k=1:nbr

           if n1(k)>0&nr(k)>0

               ybus(n1(k),nr(k))=-y(k);

               ybus(nr(k),n1(k))=-y(k);

           end

       end

       for n=1:nbus

           for k=1:nbr

               if n1(k)==n|nr(k)==n

                   ybus(n,n)=ybus(n,n)+y(k);

               else

               end

           end

       end

       y;

       ybus

       p=(busdata(:,7)-busdata(:,5));

       q=(busdata(:,8)-busdata(:,6));

       s=complex(p,q);

       v=busdata(:,3);

       void =v;

       iter=0;

       while(r>=epsilon)

       for k=2:nbus

           sum =0;

           for q=1:nbus

               if(q~=k)

               sum=sum+(ybus(k,q)*v(q));

               else

               end

           end

           sum;

           v(k,1)=((conj(s(k,1)/conj(v(k,1)))-sum)/ybus(k,k));

           iter+1;

           p;

           v;

       end

       v

       dv=abs(void-v);

       r=max(dv);

       void=v;

       iter=iter+1

       end

       for ab=1:nbus

           for ba=1:nbus

           if (ab~=ba)

               ibus(ab,ba)=v(ab,1)-v(ba,1)*(-ybus(ab,ba));

           else

               ibus(ab,ba)=0;

           end

           end

       end

       ibus

       for ab=1:nbus

           for ba=1:nbus

               if(ab~=ba)

                   lineflow(ab,ba)=(v(ab)*conj(ibus(ab,ba)));

               end

           end

       end

       lineflow

       islack=0;

       for ab=1:nbus

           for ba=1:nbus

               if(ab==1)

                   islack=islack+ibus(1,ba);

               else

               end

           end

       end

       islack

       slackpower =conj(v(1))*(islack)

       %powerloss

       for ab=1:nbus

           for ba=1:nbus

               if ba>ab

                   pv=ibus(ab,ba)*ibus(ab,ba)*(real(1/-ybus(ab,ba)));

               else

               end

           end

       end

       pv