fex_nrloop.m

这是一个很适合研究和学习用的电力系统仿真软件

function fex_nrloop(tol,itermax)
% FEX_NRLOOP -> called by FEX_ABCD
% Perform Newton Raphson for non linear algebraic equation 
% solution
%
% Author:  Alberto Del Rosso
% Date:    June 2004
% Update:
% Version: preliminary version
% E-mail:  [e-mail]
% 
% I M P O R T A N T   N O T E
% ===========================
% This external function is not part of the original PSAT software
% package although it uses many of PSAT functions. It has not been
% approved yet by Dr. Federico Milano to be distributed along with
% PSAT files.

global Fig Settings Snapshot Hdl
global Bus File DAE Theme SW PV
global Varout Line Path clpsat

% initialize NR loop
iteration = 1;
incY(1) = 1;
if isempty(DAE.f), DAE.f = 0; end
fn = DAE.f;
% Newton-Raphson loop
while max(abs(incY)) > tol
  if (iteration > itermax),  break,  end
  drawnow
  if Fig.main
    if ~get(Fig.main,'UserData'), break, end
  end
  if isempty(Line.Y)
    DAE.gp = zeros(Bus.n,1);
    DAE.gq = zeros(Bus.n,1);
    DAE.J11 = sparse(Bus.n,Bus.n);
    DAE.J21 = sparse(Bus.n,Bus.n);
    DAE.J12 = sparse(Bus.n,Bus.n);
    DAE.J22 = sparse(Bus.n,Bus.n);
  end
  % DAE equations
  fm_call('i');
  DAE.Fy(:,SW.bus) = 0;
  DAE.Gx(SW.bus,:) = 0;
  DAE.Fy(:,Bus.n+SW.bus) = 0;
  DAE.Gx(Bus.n+SW.bus,:) = 0;
  DAE.Fy(:,Bus.n+PV.bus) = 0;
  DAE.Gx(Bus.n+PV.bus,:) = 0;
  
  % check for islanded buses
  if ~isempty(Bus.island)
    kkk = Bus.island;
    DAE.Jlfv(kkk,:) = 0;
    DAE.Jlfv(:,kkk) = 0;
    DAE.Jlfv(:,kkk+Bus.n) = 0;
    DAE.Jlfv(kkk+Bus.n,:) = 0;
    DAE.Jlfv(kkk,kkk) = speye(length(kkk));
    DAE.Jlfv(kkk+Bus.n,kkk+Bus.n) = speye(length(kkk));
    DAE.g(kkk) = 0;
    DAE.g(kkk+Bus.n) = 0;
    DAE.V(kkk) = 1e-6;
    DAE.a(kkk) = 0;
  end
  % solve for algebraic variables increment (linear system equation)
  % solution of linearized system within Newton-Raphson iterative
  % process
  incY=DAE.Jlfv\-DAE.g;   
  DAE.a = DAE.a + incY(1: Bus.n);
  DAE.V = DAE.V + incY(Bus.n+1: 2*Bus.n);
  iteration = iteration + 1;
end