function [z,l,u,w,x,y,PG,QR,Pg,Qg,e,f,T,G,B,capK,meritF]=snd_solve(HSB,const,mu,r,nodeNum,balNode,pgNum,pvNum,pgvNum,xNum,... transNum,pgNode,pvNode,m,Lz,Lw,LZ,UW,dgxT,L1,U1,transI,transJ,transG,transB,transG1,transB1,Tn,capKn,capNum,... capI,Tn1,capKn1,x,y,z,l,u,w,G,B,PG,QR,Pg,Qg,opfGoal,Gama,Arlfa,meritF,k2,T,capK,lineNum,lineI,lineJ,Pd,Qd,gxx,gsx,mut,mub) %程序功能:光滑处理过程,求解修正方程 %编写时间:2010年11月 %% 求解方程,求出各变量的修正量 dX=-HSB\const; dx=dX(1:xNum,1); dy=dX(xNum+1:xNum+m,1); dl=dgxT*dx+Lz; du=-dgxT*dx-Lw; Ll=(l.*z-mu*ones(r,1)); Lu=(u.*w+mu*ones(r,1)); diagZ=sparse(1:r,1:r,z,r,r); diagW=sparse(1:r,1:r,w,r,r); dz=-LZ*dgxT*dx-L1*(diagZ*Lz+Ll); dw=UW*dgxT*dx+U1*(diagW*Lw-Lu); L3=-L1*(diagZ*Lz+Ll); L4=U1*(diagW*Lw-Lu); L5=Lz; L6=Lw; f_NCP1=[const;L3;L4;L5;L6]; f_NCPsum1=(f_NCP1'*f_NCP1)/2; meritF=[meritF;f_NCPsum1]; minusdz=find(dz<0); dzmin=min(-z(minusdz)./dz(minusdz)); minusdl=find(dl<0); dlmin=min(-l(minusdl)./dl(minusdl)); minusdw=find(dw>0); dwmin=min(-w(minusdw)./dw(minusdw)); minusdu=find(du<0); dumin=min(-u(minusdu)./du(minusdu)); ap=min(dlmin,dumin); ad=min(dzmin,dwmin); stepp=0.9995*min(ap,1); stepd=0.9995*min(ad,1); %% 一维线性搜索,求步长 [mk,meritF,stepp,stepd]=linesearch(r,nodeNum,balNode,pgNum,pvNum,pgvNum,... xNum,transNum,pgNode,pvNode,m,x,y,z,l,u,w,G,B,transI,transJ,transG,transB,transG,transB,... Tn,capKn,capNum,capI,Tn1,capKn1,PG,QR,Pg,Qg,dx,dl,du,dy,dz,dw,... opfGoal,Gama,Arlfa,meritF,k2,mu,T,capK,lineNum,lineI,lineJ,Pd,Qd,gxx,gsx,mut,mub,stepp,stepd); %% 更新变量 x=x+stepp*dx; l=l+stepp*dl; u=u+stepp*du; y=y+stepd*dy; z=z+stepd*dz; w=w+stepd*dw; %保留上一次计算的变比,用于修正节点导纳矩阵 Tlast=T; capKlast=capK; %给控制变量、状态变量重新赋值 T=x(1:transNum,1); capK=x(transNum+pgvNum+1:transNum+pgvNum+capNum); PG(pgNode,1)=x(transNum+1:transNum+pgNum,1); QR(pvNode,1)=x(transNum+pgNum+1:transNum+pgNum+pvNum,1); e=x(transNum+pgNum+pvNum+capNum+1:transNum+pgNum+pvNum+capNum+nodeNum,1); f=x(transNum+pgNum+pvNum+capNum+nodeNum+1:xNum,1); Pg(pgNode)=PG(pgNode); Qg(pvNode)=QR(pvNode); %% 变比变化后需修正节点导纳矩阵 transG=transG1; transB=transB1; G=G-sparse(transI,transI,(T.*T-Tlast.*Tlast).*transG,nodeNum,nodeNum)... +sparse(transI,transJ,transG.*(T-Tlast),nodeNum,nodeNum)... +sparse(transJ,transI,transG.*(T-Tlast),nodeNum,nodeNum); B=B-sparse(transI,transI,(T.*T-Tlast.*Tlast).*transB,nodeNum,nodeNum)... +sparse(transI,transJ,transB.*(T-Tlast),nodeNum,nodeNum)... +sparse(transJ,transI,transB.*(T-Tlast),nodeNum,nodeNum)... +sparse(capI,capI,capK-capKlast,nodeNum,nodeNum);