diff --git a/run.m b/run.m index 578125f..c555055 100644 --- a/run.m +++ b/run.m @@ -33,6 +33,8 @@ lineX=newwordParameter.line.lineX; lineB2=newwordParameter.line.lineB2; lineG=real(1./(lineR+1j*lineX)); lineB=imag(1./(lineR+1j*lineX)); +% 处理线路电阻或电抗为0的情况，即消除NaN +zerosRXInd=union(find(abs(lineR)<1e-5),find(abs(lineX)<1e-5)); cmpBranchI=BranchI( cmpV,lineI,lineJ,lineR,lineX );%复数支路电流 rBranchI=abs(cmpBranchI);% 支路电流幅值 mBranchI=rBranchI.*(normrnd(0,sigma,length(rBranchI),1)+1);%支路电流量测量 @@ -47,6 +49,9 @@ mBranchP=rBranchP.*(normrnd(0,sigma,length(rBranchP),1)+1);%支路 rBranchQ=BranchQ( cmpV,cmpBranchI,lineI,lineB2 ); % rBranchQ=rBranchQ(abs(rBranchQ)>1e-5); mBranchQ=rBranchQ.*(normrnd(0,sigma,length(rBranchQ),1)+1);%支路功率量测量 +% 处理线路电阻或电抗为0的情况，即消除NaN +mBranchP=mBranchP(setdiff( 1:length(mBranchP),zerosRXInd)); +mBranchQ=mBranchQ(setdiff( 1:length(mBranchQ),zerosRXInd)); %% 注入功率 rPD=PD; PDi=find(PD~=0); @@ -106,7 +111,7 @@ Iteration=0; while maxD>1e-3; dV_dV=sparse(1:length(mVolt),1:length(mVolt),1,length(mVolt),length(mVolt));%电压量测量的微分 dV_dTyta=sparse(length(mVolt),length(mVolt)); -dLPij_dVi=sparse(1:length(lineI),lineJ, ... +dLPij_dVi=sparse(1:length(lineI),lineI, ... -SEVolt(lineJ).*( ... lineG.*cos(SEVAngel(lineI)-SEVAngel(lineJ)) +lineB.*sin(SEVAngel(lineI)-SEVAngel(lineJ))... )... @@ -117,7 +122,7 @@ dLPij_dVj=sparse(1:length(lineI),lineJ, ... lineG.*cos(SEVAngel(lineI)-SEVAngel(lineJ)) +lineB.*sin(SEVAngel(lineI)-SEVAngel(lineJ))... ) ... ,length(lineI),length(mVolt));%线路的 -dLPij_dThetai=sparse(1:length(lineI),lineJ, ... +dLPij_dThetai=sparse(1:length(lineI),lineI, ... SEVolt(lineI).*SEVolt(lineJ).*( ... lineG.*sin(SEVAngel(lineI)-SEVAngel(lineJ)) -lineB.*cos(SEVAngel(lineI)-SEVAngel(lineJ))... )... @@ -128,9 +133,9 @@ dLPij_dThetaj=sparse(1:length(lineI),lineJ, ... )... ,length(lineI),length(mVolt));%线路的 -dLQij_dVi=sparse(1:length(lineI),lineJ, ... +dLQij_dVi=sparse(1:length(lineI),lineI, ... -SEVolt(lineJ).*( ... - lineG.*cos(SEVAngel(lineI)-SEVAngel(lineJ)) -lineB.*sin(SEVAngel(lineI)-SEVAngel(lineJ))... + lineG.*sin(SEVAngel(lineI)-SEVAngel(lineJ)) -lineB.*cos(SEVAngel(lineI)-SEVAngel(lineJ))... )... -2*(lineB+lineB2).*SEVolt(lineI) ... ,length(lineI),length(mVolt));%线路的 @@ -139,7 +144,7 @@ dLQij_dVj=sparse(1:length(lineI),lineJ, ... lineG.*sin(SEVAngel(lineI)-SEVAngel(lineJ)) -lineB.*cos(SEVAngel(lineI)-SEVAngel(lineJ))... ) ... ,length(lineI),length(mVolt));%线路的 -dLQij_dThetai=sparse(1:length(lineI),lineJ, ... +dLQij_dThetai=sparse(1:length(lineI),lineI, ... -SEVolt(lineI).*SEVolt(lineJ).*( ... lineG.*cos(SEVAngel(lineI)-SEVAngel(lineJ)) +lineB.*sin(SEVAngel(lineI)-SEVAngel(lineJ))... ) ... @@ -150,7 +155,6 @@ dLQij_dThetaj=sparse(1:length(lineI),lineJ, ... ) ... ,length(lineI),length(mVolt));%线路的 % 处理线路电阻或电抗为0的情况，即消除NaN -zerosRXInd=union(find(abs(lineR)<1e-5),find(abs(lineX)<1e-5)); dLPij_dVi=dLPij_dVi(setdiff( 1:size(dLPij_dVi,1),zerosRXInd),:); dLPij_dVj=dLPij_dVj(setdiff( 1:size(dLPij_dVj,1),zerosRXInd),:); dLPij_dThetai=dLPij_dThetai(setdiff( 1:size(dLPij_dThetai,1),zerosRXInd),:); @@ -160,16 +164,26 @@ dLQij_dVj=dLQij_dVj(setdiff( 1:size(dLQij_dVj,1),zerosRXInd),:); dLQij_dThetai=dLQij_dThetai(setdiff( 1:size(dLQij_dThetai,1),zerosRXInd),:); dLQij_dThetaj=dLQij_dThetaj(setdiff( 1:size(dLQij_dThetaj,1),zerosRXInd),:); % 对量测值做同样处理 -mBranchP=mBranchP(setdiff( 1:length(mBranchP),zerosRXInd)); -mBranchQ=mBranchQ(setdiff( 1:length(mBranchQ),zerosRXInd)); + %% 进入迭代 -H=[dV_dV ,dV_dTyta; - dLPij_dVi+dLPij_dVj, dLPij_dThetai+dLPij_dThetaj; - dLQij_dVi+dLQij_dVj, dLQij_dThetai+dLQij_dThetaj];%jacobi -% H=[dV_dV]; -h=[SEVolt;mBranchP;mBranchQ]; -z=[mVolt;mBranchP;mBranchQ]; -W=sparse(1:length(mVolt),1:length(mVolt),1/sigma.^2,length(mVolt),length(mVolt)); +% H=[dV_dV ,dV_dTyta; +% dLPij_dVi+dLPij_dVj, dLPij_dThetai+dLPij_dThetaj; +% dLQij_dVi+dLQij_dVj, dLQij_dThetai+dLQij_dThetaj];%jacobi + +H=[dV_dV; + dLPij_dVi+dLPij_dVj ; + dLQij_dVi+dLQij_dVj ];%jacobi + +SEBranchI=BranchI( SEVolt.*exp(1j*SEVAngel),lineI,lineJ,lineR,lineX );%复数支路电流 +SEBranchP=BranchP( SEVolt.*exp(1j*SEVAngel),SEBranchI,lineI,lineB2 ); +SEBranchQ=BranchQ( SEVolt.*exp(1j*SEVAngel),SEBranchI,lineI,lineB2 ); +SEBranchP=SEBranchP(setdiff( 1:length(SEBranchP),zerosRXInd)); +SEBranchQ=SEBranchQ(setdiff( 1:length(SEBranchQ),zerosRXInd)); +% h=[SEVolt;SEBranchP;SEBranchQ]; +% z=[mVolt;mBranchP;mBranchQ]; +h=[SEVolt;SEBranchP;SEBranchQ]; +z=[mVolt;mBranchP;mBranchP]; +W=sparse(1:length(h),1:length(h),1/sigma.^2,length(h),length(h)); G=H'*W*H; g=-H'*W*(z-h); % 平衡节点相角恒定; @@ -182,11 +196,10 @@ maxD=max(abs(dX)) % 更新变量 SEVolt=SEVolt+dX(1:length(mVolt)); Iteration=Iteration+1; -SEVAngel=SEVAngel+dX(length(mVolt)+1:end); +% SEVAngel=SEVAngel+dX(length(mVolt)+1:end); end %% 输出结果 fprintf('迭代%d次\n',Iteration); -h=SEVolt; fval=full((z-h)'*W*(z-h)); fprintf('目标函数为:%f\n',fval); fprintf('相对误差\n'); diff --git a/鍏紡/鍏紡.tex b/鍏紡/鍏紡.tex index 7a861f7..781bc6b 100644 --- a/鍏紡/鍏紡.tex +++ b/鍏紡/鍏紡.tex @@ -25,7 +25,7 @@ \end{equation} \begin{equation} \begin{aligned} -\dot{S}^*_{12}&=V_1e^{j \theta_1} +\dot{S}_{12}&=V_1e^{j \theta_1} \dot{I}^*_{12} \\ &=V_1e^{j \theta_1} (V_1e^{-j \theta_1} - V_2e^{-j \theta_2}) @@ -50,15 +50,74 @@ 鏈夊姛浼犺緭鍔熺巼 \begin{equation} \begin{aligned} -\dot{P}_{ij}&=[V_1^2-V_1V_2cos(\theta_1 - \theta_2)]G_{ij}-V_1V_2sin (\theta_1 - \theta_2)B_{ij} \\ +P_{ij}&=[V_1^2-V_1V_2cos(\theta_1 - \theta_2)]G_{ij}-V_1V_2sin (\theta_1 - \theta_2)B_{ij} \\ &=V_1^2G_{ij}-V_1V_2[cos(\theta_1 - \theta_2)G_{ij}+sin (\theta_1 - \theta_2)B_{ij}] \end{aligned} \end{equation} 鏃犲姛浼犺緭鍔熺巼 \begin{equation} \begin{aligned} -\dot{Q}_{ij}&=-[V_1^2-V_1V_2cos(\theta_1 - \theta_2)]B_{ij}-V_1V_2sin (\theta_1 - \theta_2)G_{ij} \\ +Q_{ij}&=-[V_1^2-V_1V_2cos(\theta_1 - \theta_2)]B_{ij}-V_1V_2sin (\theta_1 - \theta_2)G_{ij} \\ &=-V_1^2B_{ij}-V_1V_2[sin(\theta_1 - \theta_2)G_{ij}-cos (\theta_1 - \theta_2)B_{ij}] \end{aligned} \end{equation} +绾胯矾鏈夊姛鍔熺巼Jacobi +\begin{equation} +\begin{aligned} +\frac{\partial P_{ij}}{\partial V_1}= +2V_1G_{ij}-V_2[cos(\theta_1 - \theta_2)G_{ij}+sin (\theta_1 - \theta_2)B_{ij}] +\end{aligned} +\end{equation} + +\begin{equation} +\begin{aligned} +\frac{\partial P_{12}}{\partial V_2}= +-V_1[cos(\theta_1 - \theta_2)G_{ij}+sin (\theta_1 - \theta_2)B_{ij}] +\end{aligned} +\end{equation} + +\begin{equation} +\begin{aligned} +\frac{\partial P_{12}}{\partial \theta_1}&= +-V_1V_2[-sin(\theta_1 - \theta_2)G_{ij}+cos (\theta_1 - \theta_2)B_{ij}] \\ +&=V_1V_2[sin(\theta_1 - \theta_2)G_{ij}-cos (\theta_1 - \theta_2)B_{ij}] +\end{aligned} +\end{equation} + +\begin{equation} +\begin{aligned} +\frac{\partial P_{12}}{\partial \theta_2}&= +-V_1V_2[sin(\theta_1 - \theta_2)G_{ij}-cos (\theta_1 - \theta_2)B_{ij}] \\ +\end{aligned} +\end{equation} + +\begin{equation} +\begin{aligned} +\frac{\partial Q_{12}}{\partial V_1}&= +-2V_1B_{12}-V_2[sin(\theta_1 - \theta_2)G_{ij}-cos (\theta_1 - \theta_2)B_{ij}] +\end{aligned} +\end{equation} + +\begin{equation} +\begin{aligned} +\frac{\partial Q_{12}}{\partial V_2}&= +-V_1[sin(\theta_1 - \theta_2)G_{ij}-cos (\theta_1 - \theta_2)B_{ij}] +\end{aligned} +\end{equation} + +\begin{equation} +\begin{aligned} +\frac{\partial Q_{12}}{\partial \theta_1}&= +-V_1V_2[cos(\theta_1 - \theta_2)G_{ij}+sin (\theta_1 - \theta_2)B_{ij}] +\end{aligned} +\end{equation} + +\begin{equation} +\begin{aligned} +\frac{\partial Q_{12}}{\partial \theta_2}&= +-V_1V_2[-cos(\theta_1 - \theta_2)G_{ij}-sin (\theta_1 - \theta_2)B_{ij}] \\ +&=V_1V_2[cos(\theta_1 - \theta_2)G_{ij}+sin (\theta_1 - \theta_2)B_{ij}] +\end{aligned} +\end{equation} + \end{document}