diff --git a/run.asv b/run.asv
index 8ff07d4..da0a6aa 100644
--- a/run.asv
+++ b/run.asv
@@ -1,13 +1,81 @@
 addpath('.\Powerflow')
-[~, ~, ~, ~,Volt,Vangle,Y,Yangle,newwordParameter]=pf('ieee4.dat', '0');
+[~, ~, ~, ~,Volt,Vangle,Y,Yangle,r,c,newwordParameter,PG,QG,PD,QD,Balance]=pf('ieee30.dat', '0');
 %% 量测量
-% 电压 节点电流 支路 节点功率 支路功率
+% 电压 节点电流 支路电流 节点功率 支路功率
 %%
 %% 状态量
 % 电压 相角
 %%
 %% 开始生成量测量
 sigma=0.05;% 标准差
-%电压
-mVolt=Volt.*(normrnd(0,sigma,length(Volt),1)+1);
-
+%% 电压
+%电压幅值
+rVolt=Volt; %幅值
+mVolt=rVolt.*(normrnd(0,sigma,length(Volt),1)+1);%电压量测量
+%% 电流
+%注入电流
+cmpY=Y.*exp(1j*sparse(r,c,Yangle,length(Y),length(Y)));%复数导纳矩阵
+cmpV=Volt.*exp(1j*Vangle); %复数电压
+cmpI=cmpY*cmpV;% 注入电流
+rI=abs(cmpI); %注入电流量测量要的是电流幅值
+mI=rI.*(normrnd(0,sigma,length(rI),1)+1);%电流量测量
+%% 支路电流
+% 支路电流
+lineI=newwordParameter.line.lineI;
+lineJ=newwordParameter.line.lineJ;
+lineR=newwordParameter.line.lineR;
+lineX=newwordParameter.line.lineX;
+lineB2=newwordParameter.line;
+cmpBranchI=(cmpV(lineI)-cmpV(lineJ))./(lineR+1j*lineX);%复数支路电流
+rBranchI=abs(cmpBranchI);% 支路电流幅值
+mBranchI=rBranchI.*(normrnd(0,sigma,length(rBranchI),1)+1);%支路电流量测量
+%% 支路功率
+rBranchP=real((cmpV(lineI)-cmpV(lineJ)).*conj(cmpBranchI));
+mBranchP=rBranchP.*(normrnd(0,sigma,length(rBranchP),1)+1);%支路功率量测量
+rBranchQ=imag((cmpV(lineI)-cmpV(lineJ)).*conj(cmpBranchI));
+mBranchQ=rBranchQ.*(normrnd(0,sigma,length(rBranchQ),1)+1);%支路功率量测量
+%% 注入功率
+rPD=PD(PD~=0);
+PDi=find(PD~=0);
+rQD=QD(QD~=0);
+QDi=find(QD~=0);
+rPG=PG(PG~=0);
+PGi=find(PG~=0);
+rQG=QG(QG~=0);
+QGi=find(QG~=0);
+mPD=rPD.*(normrnd(0,sigma,length(rPD),1)+1);
+mQD=rQD.*(normrnd(0,sigma,length(rQD),1)+1);
+mPG=rPG.*(normrnd(0,sigma,length(rPG),1)+1);
+mQG=rQG.*(normrnd(0,sigma,length(rQG),1)+1);
+%% 冗余度计算
+stateVarCount=2*length(Volt);
+measurements=length(mVolt)+length(mI)+length(mBranchI)+length(mBranchP)+length(mBranchQ)+length(mPG)+length(mQG);
+fprintf('冗余度 %f\n',measurements/stateVarCount);
+%% 进入状态估计计算
+SEVolt=sdpvar(length(Volt),1);
+SEVAngel=sdpvar(length(Vangle),1);
+Objective=(SEVolt-mVolt)'*(1./sigma^2*eye(length(mVolt)))*(SEVolt-mVolt);%%电压
+%% 支路电流
+cmpSEV=SEVolt.*exp(1j*SEVAngel); %复数电压
+cmpSEBranchI=(cmpV(lineI)-cmpV(lineJ))./(lineR+1j*lineX);%复数支路电流
+SEBranchI=abs(cmpSEBranchI);% 支路电流幅值
+Objective=Objective+(SEBranchI-mBranchI)'*(1./sigma^2*eye(length(mBranchI)))*(SEBranchI-mBranchI);%%电流
+%% 支路功率
+SEBranchP=real((cmpSEV(lineI)-cmpSEV(lineJ)).*conj(SEBranchI));
+SEBranchQ=imag((cmpSEV(lineI)-cmpSEV(lineJ)).*conj(SEBranchI));
+Objective=Objective+(SEBranchP-mBranchP)'*(1./sigma^2*eye(length(mBranchP)))*(SEBranchP-mBranchP);
+Objective=Objective+(SEBranchQ-mBranchQ)'*(1./sigma^2*eye(length(mBranchQ)))*(SEBranchQ-mBranchQ);
+%% 发电机注入功率
+%% 0注入节点
+zerosInjectionIndex=1:length(Volt);
+zerosInjectionIndex=zerosInjectionIndex( ~(PD~=0|QD~=0|PG~=0|QG~=0) );
+%% YALMIP求解
+PQ=diag(SEVolt)*cmpY*SEVolt;
+zeros
+Constraints=[SEVAngel(Balance)==0,];
+% Constraints=[[zeros(length(c)) A' -eye(length(lbounds))]*x==-c;[A zeros(length(b)) zeros(length(b),length(lbounds))]*x<=b;10>=x>=0];
+options = sdpsettings('verbose',2,'solver','ipopt');
+solvesdp(Constraints,Objective,options)
+double(Objective)
+fprintf('相对误差\n');
+(abs(rVolt-double(SEVolt)))./(rVolt)
diff --git a/run.m b/run.m
index 5cfa25b..3ee6150 100644
--- a/run.m
+++ b/run.m
@@ -1,5 +1,5 @@
 addpath('.\Powerflow')
-[~, ~, ~, ~,Volt,Vangle,Y,Yangle,r,c,newwordParameter]=pf('ieee4.dat', '0');
+[~, ~, ~, ~,Volt,Vangle,Y,Yangle,r,c,newwordParameter,PG,QG,PD,QD,Balance]=pf('ieee30.dat', '0');
 %% 量测量
 % 电压 节点电流 支路电流 节点功率 支路功率
 %%
@@ -28,4 +28,55 @@ lineX=newwordParameter.line.lineX;
 lineB2=newwordParameter.line;
 cmpBranchI=(cmpV(lineI)-cmpV(lineJ))./(lineR+1j*lineX);%复数支路电流
 rBranchI=abs(cmpBranchI);% 支路电流幅值
-mBranchI=rBranchI.*(normrnd(0,sigma,length(rBranchI),1)+1);%支路电流量测量
\ No newline at end of file
+mBranchI=rBranchI.*(normrnd(0,sigma,length(rBranchI),1)+1);%支路电流量测量
+%% 支路功率
+rBranchP=real((cmpV(lineI)-cmpV(lineJ)).*conj(cmpBranchI));
+mBranchP=rBranchP.*(normrnd(0,sigma,length(rBranchP),1)+1);%支路功率量测量
+rBranchQ=imag((cmpV(lineI)-cmpV(lineJ)).*conj(cmpBranchI));
+mBranchQ=rBranchQ.*(normrnd(0,sigma,length(rBranchQ),1)+1);%支路功率量测量
+%% 注入功率
+rPD=PD(PD~=0);
+PDi=find(PD~=0);
+rQD=QD(QD~=0);
+QDi=find(QD~=0);
+rPG=PG(PG~=0);
+PGi=find(PG~=0);
+rQG=QG(QG~=0);
+QGi=find(QG~=0);
+mPD=rPD.*(normrnd(0,sigma,length(rPD),1)+1);
+mQD=rQD.*(normrnd(0,sigma,length(rQD),1)+1);
+mPG=rPG.*(normrnd(0,sigma,length(rPG),1)+1);
+mQG=rQG.*(normrnd(0,sigma,length(rQG),1)+1);
+%% 冗余度计算
+stateVarCount=2*length(Volt);
+measurements=length(mVolt)+length(mI)+length(mBranchI)+length(mBranchP)+length(mBranchQ)+length(mPG)+length(mQG);
+fprintf('冗余度 %f\n',measurements/stateVarCount);
+%% 进入状态估计计算
+SEVolt=sdpvar(length(Volt),1);
+SEVAngel=sdpvar(length(Vangle),1);
+Objective=(SEVolt-mVolt)'*(1./sigma^2*eye(length(mVolt)))*(SEVolt-mVolt);%%电压
+%% 支路电流
+cmpSEV=SEVolt.*exp(1j*SEVAngel); %复数电压
+cmpSEBranchI=(cmpV(lineI)-cmpV(lineJ))./(lineR+1j*lineX);%复数支路电流
+SEBranchI=abs(cmpSEBranchI);% 支路电流幅值
+Objective=Objective+(SEBranchI-mBranchI)'*(1./sigma^2*eye(length(mBranchI)))*(SEBranchI-mBranchI);%%电流
+%% 支路功率
+SEBranchP=real((cmpSEV(lineI)-cmpSEV(lineJ)).*conj(SEBranchI));
+SEBranchQ=imag((cmpSEV(lineI)-cmpSEV(lineJ)).*conj(SEBranchI));
+Objective=Objective+(SEBranchP-mBranchP)'*(1./sigma^2*eye(length(mBranchP)))*(SEBranchP-mBranchP);
+Objective=Objective+(SEBranchQ-mBranchQ)'*(1./sigma^2*eye(length(mBranchQ)))*(SEBranchQ-mBranchQ);
+%% 发电机注入功率
+%% 0注入节点
+zerosInjectionIndex=1:length(Volt);
+zerosInjectionIndex=zerosInjectionIndex( ~(PD~=0|QD~=0|PG~=0|QG~=0) );
+%% YALMIP求解
+PQ=diag(SEVolt)*cmpY*SEVolt;
+zeroInjP=real(PQ(zerosInjectionIndex));%% 0注入节点
+zeroInjQ=imag(PQ(zerosInjectionIndex));%% 0注入节点
+Constraints=[SEVAngel(Balance)==0,zeroInjP==0,zeroInjP==0];
+% Constraints=[[zeros(length(c)) A' -eye(length(lbounds))]*x==-c;[A zeros(length(b)) zeros(length(b),length(lbounds))]*x<=b;10>=x>=0];
+options = sdpsettings('verbose',2,'solver','ipopt');
+solvesdp(Constraints,Objective,options)
+double(Objective)
+fprintf('相对误差\n');
+(abs(rVolt-double(SEVolt)))./(rVolt)