% Manipulate, plot power flow example from class

V1 = 1.0; % slack bus
Zl = i*0.5;
Q2 = 0.0; % PQ (load) bus
P2 = 0.5;
pldim = 100;

S2 = P2 + i*Q2;
Yl = 1/Zl;
resid = zeros(pldim,pldim);
resid_inds = zeros(pldim,pldim);
Re_min = -0.5;
Re_max = 2.0;
Im_min = -1.0;
Im_max = 2.0;
Re_inc = (Re_max-Re_min)/pldim;
Im_inc = (Im_max-Im_min)/pldim;
Re_ind = [Re_min:Re_inc:Re_max-Re_inc];
Im_ind = [Im_min:Im_inc:Im_max-Im_inc];

for k=1:pldim
 for m = 1:pldim
	V2_re = Re_min+(k-1)*Re_inc;
	V2_im = Im_min+(m-1)*Im_inc;
	V2 = V2_re + i*V2_im;
	resid(k,m) = (conj(S2) + Yl*V1*conj(V2) - Yl*V2*conj(V2));
 end
end

figure(1)
mesh(Im_ind,Re_ind,abs(resid))
xlabel('Imaginary Part of V2');
ylabel('Real Part of V2');
title('Residual of PQ bus equation')
figure(2)
mesh(Im_ind,Re_ind,unwrap(angle(resid)))
xlabel('Imaginary Part of V2');
ylabel('Real Part of V2');
title('Phase of Residual of PQ bus equation')
%view(-37.5, -70);
