form_jac.m

这个程序是采用模块化的方式对潮流进行求解

function J=form_jac(bus,Y)  % 该子程序是用来求取jacci矩阵
global n m V G p q B delta 
    
    for i=1: n-1                    
         for j=1: n-1
            if i==j
               J(i,j)=0;
               for t=1:n  
                   if t~=i
                       J(i,j)= J(i,j)+V(t)*(G(i,t)*sin(delta(i)-delta(t))-B(i,t)*cos(delta(i)-delta(t)));
                   end
               end
                J(i,j)= J(i,j)*V(i);
            else 
                J(i,j)=-V(i)*V(j)*(G(i,j)*sin(delta(i)-delta(j))-B(i,j)*cos(delta(i)-delta(j)));
            end
         end
     end        %计算雅可比矩阵 H；   
   
     for i=1: n-1
       for j=1: m
            if i==j
                J(i,j+n-1)=0;
               for t=1:n  
                   if t~=i
                       J(i,j+n-1)= J(i,j+n-1)+V(t)*(G(i,t)*cos(delta(i)-delta(t))+B(i,t)*sin(delta(i)-delta(t)));
                   end
               end
                J(i,j+n-1)= -J(i,j+n-1)*V(i)-2*V(i)*V(i)*G(i,i);
          else
            J(i,j+n-1)=-V(i)*V(j)*(G(i,j)*cos(delta(i)-delta(j))+B(i,j)*sin(delta(i)-delta(j)));
        end
       end
     end       %计算雅可比矩阵 N；

    for i=1: m
        for j=1: n-1
           if i==j
               J(i+n-1,j)=0;
               for t=1:n  
                   if t~=i
                       J(i+n-1,j)= J(i+n-1,j)+V(t)*(G(i,t)*cos(delta(i)-delta(t))+B(i,t)*sin(delta(i)-delta(t)));
                   end
               end
                J(i+n-1,j)= -J(i+n-1,j)*V(i);
           else
               J(i+n-1,j)=V(i)*V(j)*(G(i,j)*cos(delta(i)-delta(j))+B(i,j)*sin(delta(i)-delta(j)));
           end
        end
     end       %计算雅可比矩阵 K；

    for i=1: m
         for j=1: m
              if i==j
                  J(n-1+i,n-1+j)=0;
               for t=1:n  
                   if t~=i
                       J(n-1+i,n-1+j)= J(n-1+i,n-1+j)+V(t)*(G(i,t)*sin(delta(i)-delta(t))-B(i,t)*cos(delta(i)-delta(t)));
                   end
               end
                J(n-1+i,n-1+j)= -J(n-1+i,n-1+j)*V(i)+2*V(i)*V(i)*B(i,i);
             else
                  J(n-1+i,n-1+j)=-V(i)*V(j)*(G(i,j)*sin(delta(i)-delta(j))-B(i,j)*cos(delta(i)-delta(j)));
             end
         end
      end       %计算雅可比矩阵 L；