diff --git a/OPF.m b/OPF.m
index a29d989..88de316 100644
--- a/OPF.m
+++ b/OPF.m
@@ -74,10 +74,10 @@ while(k<=kmax+10 && maxD> EPS)
     SA=VoltpA.*conj(VoltpA.*sparse(cap.capNode,1,1j*cap.capB(:,1),busNum,1));
     SB=VoltpB.*conj(VoltpB.*sparse(cap.capNode,1,1j*cap.capB(:,2),busNum,1));
     SC=VoltpC.*conj(VoltpC.*sparse(cap.capNode,1,1j*cap.capB(:,3),busNum,1));
-    %先不要电容
-    %     SA=0;
-    %     SB=0;
-    %     SC=0;
+    %先不要电容，已经加到序的导纳矩阵中了。20150405 By 杜孟远
+        SA=0;
+        SB=0;
+        SC=0;
     iterPD=PD+real(SA+SB+SC)/3;
     iterQD=QD+imag(SA+SB+SC)/3;
     iterPhaseASpotLoadP=phaseASpotLoadP+real(SA);
diff --git a/dataRead.m b/dataRead.m
index 0b37306..4135a1e 100644
--- a/dataRead.m
+++ b/dataRead.m
@@ -27,6 +27,18 @@ fsY1=fsY1+fs31;
 fsY2=fsY2+fs32;
 phaseABCY=phaseABCY+retphaseABCY;    
 end
+%% 把电容加到矩阵中
+cap=data(zeroEntries(2)+1:zeroEntries(3)-1,:);
+if ~isempty(cap)
+capNode=nodeNum(ismember(setIJ,cap(:,1)));
+end
+%电容引起的电流
+capG=1j*cap(:,2)/1000;
+% fs3Y=3/3*diag([Zl+2*Zm,Zl-Zm,Zl-Zm]);%这里不是1/3，IEEE上的公式错了。
+fsY0=fsY0+sparse(capNode,capNode,capG,length(fsY0),length(fsY0));
+fsY1=fsY1+sparse(capNode,capNode,capG,length(fsY1),length(fsY1));
+fsY2=fsY2+sparse(capNode,capNode,capG,length(fsY2),length(fsY2));
+
 %% spot load
 busNum=length(nodeNum);
 spotloads=data(zeroEntries(3)+1:zeroEntries(4)-1,:);
diff --git a/feeder13/data1.txt b/feeder13/data1.txt
index 0bb5f24..da177f9 100644
--- a/feeder13/data1.txt
+++ b/feeder13/data1.txt
@@ -12,11 +12,10 @@
 602	684	611	300
 601	692	675	500
 0
-652 20 20 0
-671 20 220 120
+611 100 100 100
+652 200 200 200
 684 200 200 200
-675 120 220 220
-634 50 130 120
+692 200 200 200
 0
 634	120	90	120	90	160	110
 645	40	30	50	40	60	50
@@ -24,7 +23,7 @@
 652	20	10	60	30	40	40
 671	385	220	385	220	385	220
 675	68	60	290	212	285	190
-692	60	50	40	40	70	40
+692	260	200	40	40	70	40
 611	80	20	60	40	30	20
 0
 600	632	645	500