fm_pl.m

电力系统的psat

function fm_pl(flag)
% FM_PL define Polinomial Loads
%
% FM_PL(FLAG)
%       FLAG -> 0 initialization
%       FLAG -> 1 algebraic equations
%       FLAG -> 2 algebraic Jacobians
%
%Author:    Federico Milano
%Date:      11-Nov-2002
%Version:   1.0.0
%
%E-mail:    [e-mail]
%Web-site:  http://thunderbox.uwaterloo.ca/~fmilano
%
% Copyright (C) 2002-2005 Federico Milano
%
% This toolbox is free software; you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation; either version 2.0 of the License, or
% (at your option) any later version.
%
% This toolbox is distributed in the hope that it will be useful, but
% WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANDABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
% General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this toolbox; if not, write to the Free Software
% Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307,
% USA.

global Pl PQ DAE Bus Settings

V = DAE.V(Pl.bus);

switch flag
 case 0
  i = find(Pl.con(:,11));
  if i
    for j = 1:length(i)
      k = i(j);
      idx = find(PQ.bus == Pl.bus(k));
      if isempty(idx)
	fm_disp(['No PQ load found for initializing ZIP load ', ...
		 'at bus #',int2str(Pl.bus(k))])
      else
	PQ.con(idx,4) = PQ.con(idx,4) - ...
	    PQ.P0(idx)*sum(Pl.con(k,[5:7]))/100;
	PQ.con(idx,5) = PQ.con(idx,5) - ...
	    PQ.Q0(idx)*sum(Pl.con(k,[8:10]))/100;
	Pl.con(k,5)  = Pl.con(k,5)*PQ.P0(idx)/V(k)/V(k)/100;
	Pl.con(k,6)  = Pl.con(k,6)*PQ.P0(idx)/V(k)/100;
	Pl.con(k,7)  = Pl.con(k,7)*PQ.P0(idx)/100;
	Pl.con(k,8)  = Pl.con(k,8)*PQ.Q0(idx)/V(k)/V(k)/100;
	Pl.con(k,9)  = Pl.con(k,9)*PQ.Q0(idx)/V(k)/100;
	Pl.con(k,10) = Pl.con(k,10)*PQ.Q0(idx)/100;
	if PQ.con(idx,4) == 0 & PQ.con(idx,5) == 0,
	  PQ.con(idx,:) = [];
	  PQ.bus(idx) = [];
	  PQ.n = PQ.n - 1;
	  PQ.P0(idx) = [];
	  PQ.Q0(idx) = [];
	end
      end
    end
  end
  fm_disp('Initialization of ZIP loads completed.')
 case 1
  if Settings.init
    DAE.gp = DAE.gp + sparse(Pl.bus,1,(Pl.con(:,5).*V + ...
             Pl.con(:,6)).*V + Pl.con(:,7),Bus.n,1);
    DAE.gq = DAE.gq + sparse(Pl.bus,1,(Pl.con(:,8).*V + ...
             Pl.con(:,9)).*V + Pl.con(:,10),Bus.n,1);
  elseif ~isempty(Pl.init)
    i = Pl.init;
    DAE.gp = DAE.gp + sparse(Pl.bus(i),1,(Pl.con(i,5).*V(i) + ...
             Pl.con(i,6)).*V(i) + Pl.con(i,7),Bus.n,1);
    DAE.gq = DAE.gq + sparse(Pl.bus(i),1,(Pl.con(i,8).*V(i) + ...
             Pl.con(i,9)).*V(i) + Pl.con(i,10),Bus.n,1);
    end
 case 2
  if Settings.init
    DAE.J12 = DAE.J12 + sparse(Pl.bus,Pl.bus,2*Pl.con(:,5).*V + ...
              Pl.con(:,6),Bus.n,Bus.n);
    DAE.J22 = DAE.J22 + sparse(Pl.bus,Pl.bus,2*Pl.con(:,8).*V + ...
              Pl.con(:,9),Bus.n,Bus.n);
  elseif ~isempty(Pl.init)
    i = Pl.init;
    DAE.J12 = DAE.J12 + sparse(Pl.bus(i),Pl.bus(i),2*Pl.con(i,5).*V(i) + ...
			       Pl.con(i,6),Bus.n,Bus.n);
    DAE.J22 = DAE.J22 + sparse(Pl.bus(i),Pl.bus(i),2*Pl.con(i,8).*V(i) + ...
			       Pl.con(i,9),Bus.n,Bus.n);
  end
end