radial_pf

RADIAL_PF  Solves the power flow using a backward-forward sweep method.

   [mpc, success, iterations] = radial_pf(mpc,mpopt)

   Inputs:
       mpc : MATPOWER case struct with internal bus numbering
       mpopt : MATPOWER options struct to override default options
           can be used to specify the solution algorithm, output options
           termination tolerances, and more.

   Outputs:
       mpc : results struct with all fields from the input MATPOWER case,
             with solved voltages, active and reactive power flows
             and generator active and reactive power output.
       success : success flag, 1 = succeeded, 0 = failed
       iterations : number of iterations

  See also CASEFORMAT, LOADCASE, MPOPTION.

0001 function [mpc, success, iterations] = radial_pf(mpc,mpopt)
0002 %RADIAL_PF  Solves the power flow using a backward-forward sweep method.
0003 %
0004 %   [mpc, success, iterations] = radial_pf(mpc,mpopt)
0005 %
0006 %   Inputs:
0007 %       mpc : MATPOWER case struct with internal bus numbering
0008 %       mpopt : MATPOWER options struct to override default options
0009 %           can be used to specify the solution algorithm, output options
0010 %           termination tolerances, and more.
0011 %
0012 %   Outputs:
0013 %       mpc : results struct with all fields from the input MATPOWER case,
0014 %             with solved voltages, active and reactive power flows
0015 %             and generator active and reactive power output.
0016 %       success : success flag, 1 = succeeded, 0 = failed
0017 %       iterations : number of iterations
0018 %
0019 %  See also CASEFORMAT, LOADCASE, MPOPTION.
0020 
0021 %% define named indices into bus, gen, branch matrices
0022 define_constants;
0023 %% branch ordering
0024 mpc = order_radial(mpc);
0025 if ~isempty(mpc.loop)
0026     fprintf('\nNumber of detected loops: %i\n', length(mpc.loop));
0027     if mpopt.verbose > 0
0028         fprintf('\nBranches forming loops\n')
0029         fprintf('LOOP# F_BUS T_BUS\n');
0030         for i = 1:length(mpc.loop)
0031             fprintf('%5i %5i %5i\n',i,mpc.order.bus.i2e(mpc.branch(mpc.loop(i),[F_BUS T_BUS])));
0032         end
0033     end
0034     error('radial_pf: power flow algorithm %s can only handle radial networks.', mpopt.pf.alg)
0035 end
0036 %% define vectors needed for backward-forward sweep method
0037 % branch and demand data
0038 [f, t, Zb, Yb, Sd, Ysh] = ...
0039     deal(mpc.branch(:,F_BUS),mpc.branch(:,T_BUS), ...
0040          mpc.branch(:,BR_R)+1j*mpc.branch(:,BR_X),1j*mpc.branch(:,BR_B), ...
0041          mpc.bus(:,PD)+1j*mpc.bus(:,QD),mpc.bus(:,GS)+1j*mpc.bus(:,BS));
0042 nl = size(mpc.branch,1);
0043 nb = size(mpc.bus,1);
0044 Sd  =  Sd/mpc.baseMVA;
0045 Ysh = Ysh/mpc.baseMVA;
0046 tap = ones(nl, 1);        % default tap ratio = 1
0047 i = find(mpc.branch(:, TAP)); % indices of non-zero tap ratios
0048 tap(i) = mpc.branch(i, TAP);  % assign non-zero tap ratios
0049 Ybf = Yb/2 + 1./tap .* (1./tap-1) ./ Zb;
0050 Ybt = Yb/2 + (1-1./tap) ./ Zb;
0051 [Ybf(mpc.br_reverse), Ybt(mpc.br_reverse)] = ...
0052     deal(Ybt(mpc.br_reverse), Ybf(mpc.br_reverse)); % reversed branches
0053 Zb = Zb .* tap;
0054 Yd = Ysh + (sparse(f, f, Ybf, nb, nb) + sparse(t, t, Ybt, nb, nb)) * ones(nb,1);
0055 % generator data (other than the slack bus)
0056 pv = mpc.gen(2:end,GEN_BUS);
0057 Pg = mpc.gen(2:end,PG)/mpc.baseMVA;
0058 Vg = mpc.gen(2:end,VG);
0059 %% calculate voltages and power flows
0060 Vslack = mpc.gen(1,VG);
0061 switch upper(mpopt.pf.alg);
0062     case 'PQSUM'
0063         [V, Qpv, Sf, St, Sslack, iterations, success] = calc_v_pq_sum(Vslack,nb,nl,f,Zb,Ybf,Ybt,Yd,Sd,pv,Pg,Vg,mpopt);
0064     case 'ISUM'
0065         [V, Qpv, Sf, St, Sslack, iterations, success] = calc_v_i_sum(Vslack,nb,nl,f,Zb,Ybf,Ybt,Yd,Sd,pv,Pg,Vg,mpopt);
0066     case 'YSUM'
0067         [V, Qpv, Sf, St, Sslack, iterations, success] = calc_v_y_sum(Vslack,nb,nl,f,Zb,Ybf,Ybt,Yd,Sd,pv,Pg,Vg,mpopt);
0068 end
0069 mpc.success = success;
0070 %% update data matrices with solution
0071 mpc.bus(:,VM) = abs(V);
0072 mpc.bus(:,VA) = angle(V)/pi*180;
0073 mpc.branch(:,PF) = real(Sf)*mpc.baseMVA;
0074 mpc.branch(:,QF) = imag(Sf)*mpc.baseMVA;
0075 mpc.branch(:,PT) = -real(St)*mpc.baseMVA;
0076 mpc.branch(:,QT) = -imag(St)*mpc.baseMVA;
0077 [mpc.branch(mpc.br_reverse,PF), mpc.branch(mpc.br_reverse,PT)] = ...
0078     deal(mpc.branch(mpc.br_reverse,PT), mpc.branch(mpc.br_reverse,PF));
0079 [mpc.branch(mpc.br_reverse,QF), mpc.branch(mpc.br_reverse,QT)] = ...
0080     deal(mpc.branch(mpc.br_reverse,QT), mpc.branch(mpc.br_reverse,QF));
0081 mpc.gen(1,PG) = real(Sslack)*mpc.baseMVA;
0082 mpc.gen(1,QG) = imag(Sslack)*mpc.baseMVA;
0083 if ~isempty(pv)
0084     mpc.gen(2:end,QG) = Qpv*mpc.baseMVA;
0085 end
0086 %% reverse bus and branch ordering
0087 mpc.bus = mpc.bus(mpc.bus_order_inv,:);
0088 mpc.bus(:,BUS_I) = mpc.bus_order(mpc.bus(:,BUS_I));
0089 [f, t] = deal(mpc.branch(:,F_BUS),mpc.branch(:,T_BUS));
0090 [f(mpc.br_reverse), t(mpc.br_reverse)] =  deal(t(mpc.br_reverse), f(mpc.br_reverse));
0091 mpc.branch(:,[F_BUS T_BUS]) = [mpc.bus_order(f) mpc.bus_order(t)];
0092 mpc.branch = mpc.branch(mpc.branch_order_inv,:);
0093 mpc.gen(:,GEN_BUS) = mpc.bus_order(mpc.gen(:,GEN_BUS));