psat.m

基于PSAT 软件的多目标最优潮流计算用于中小型电力系统的分析和管理

%      _____      _____     _____      ______
%     /\  _  \   /\  ___\  /\  _  \   /\__  __\
%     \ \ \_\ \  \ \ \___  \ \ \_\ \  \/_/\ \_/
%      \ \  ___\  \ \___ \  \ \  __ \    \ \ \
%       \ \ \__/   \/___\ \  \ \ \ \ \    \ \ \
%        \ \_\      /\_____\  \ \_\ \_\    \ \_\
%         \/_/      \/_____/   \/_/\/_/     \/_/
%
% Copyright (C) 2002-2006 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.

echo off
if exist('command_line_psat') ~= 1
  clear all
  clear global
  command_line_psat = 0;
  if exist('OCTAVE_VERSION')
    command_line_psat = 1;
  end
end

%PSAT Power System Analysis Toolbox
%
%Author:    Federico Milano
%Date:      11-Nov-2002
%Update:    18-Feb-2003
%Update:    22-Feb-2003
%Update:    09-Jul-2003
%Update:    30-Aug-2003
%Update:    22-Feb-2004
%Update:    09-Mar-2004
%Version:   1.2.0
%
%E-mail:    [e-mail]
%Web-site:  http://thunderbox.uwaterloo.ca/~fmilano

% check host version
a = version;
a = str2num(a(1:3));
if exist('OCTAVE_VERSION')
  disp('PSAT 2 cannot run on Octave')
  return
elseif ~exist('OCTAVE_VERSION') & a < 7.0
  disp('PSAT needs a Matlab version >= 7.0')
  return
end

psatver = '2.0.0-b1';
psatdate = 'March 24,  2006';

failed = 0;
disp(' ')
disp('                        < P S A T >')
disp('          Copyright (C) 2002-2006 Federico Milano')
disp(['                     Version  ',psatver])
disp(['                      ',psatdate])
disp('  ')
disp('PSAT comes with ABSOLUTELY NO WARRANTY; type ''gnuwarranty''')
disp('for details. This is free software, and you are welcome to')
disp(['redistribute it under certain conditions; type' ...
      ' ''gnulicense'''])
disp('for details.')
disp('  ')
if exist('OCTAVE_VERSION')
  disp(['Host:         Octave ',version])
else
  disp(['Host:         Matlab ',version])
end
disp(['Session:      ', datestr(now,0)])
if command_line_psat
  disp('Usage:        Command Line')
else
  disp('Usage:        Graphical User Interface')
end
disp(['Path:         ',pwd])
disp(' ')

% searching the "PSAT" path
if exist('OCTAVE_VERSION')
  if isempty(findstr(LOADPATH,'psat'))
    path([':/',pwd,'//']);
  end
else
  psatdir = what('psat');
  if isempty(psatdir)
    psatdir = struct('path',pwd);
  end
  if isempty(findstr(path,[psatdir.path,pathsep]))
    addpath(psatdir.path)
    disp('PSAT folders have been added to the Matlab path.')
  end
  addpath([psatdir.path,filesep,'filters'])
end

% check for "lipsol" package
%lipdir = which('lipsol');
%if isempty(lipdir)
%  addpath([psatdir.path,filesep,'lipsol'])
%  disp('LIPSOL folder has been added to the Matlab path.')
%end
% initialize lipsol
%lipconfig

% check for "lp_solve" package
%lipdir = which('lp_solve');
%if isempty(lipdir)
%  addpath([psatdir.path,filesep,'lp_solve'])
%  disp('LP_SOLVE folder has been added to the Matlab path.')
%end

fm_var

Path = struct('data','','local','','psat','','pert','', ...
              'build','','images','','themes','', ...
              'filters','','temp','');

if command_line_psat
  clpsat.init = 1;
  clpsat.mesg = 1;
  clpsat.refresh = 1;
  clpsat.refreshsim = 0;
  clpsat.readfile = 1;
  clpsat.showopf = 0;
  clpsat.pq2z = 1;
  clpsat.viewrep = 0;
else
  clpsat.init = 0;
  clpsat.mesg = 1;
  clpsat.refresh = 0;
  clpsat.refreshsim = 0;
  clpsat.readfile = 1;
  clpsat.showopf = 1;
  clpsat.pq2z = 1;
  clpsat.viewrep = 1;
end

Path.local = [pwd,filesep];
if exist('OCTAVE_VERSION')
  Path.psat = strrep(which('psat'),'psat.m','');
else
  if strcmp(computer,'GLNX86')
    Path.local = strrep(Path.local,getenv('HOME'),'~');
  end
  Path.psat = [psatdir(1).path,filesep];
  if strcmp(computer,'GLNX86')
    Path.psat = strrep(Path.psat,getenv('HOME'),'~');
  end
end
Path.build   = [Path.psat,'build',  filesep];
Path.images  = [Path.psat,'images', filesep];
Path.themes  = [Path.psat,'themes', filesep];
Path.filters = [Path.psat,'filters',filesep];

if exist(Path.build) ~= 7
  disp('No "build" folder found in the current PSAT path.')
  failed = ~mkdir(Path.psat,'build');
  if failed,
    disp('The folder "build" could not be created')
  end
end
if exist(Path.images) ~= 7
  disp('No "images" folder found in the current PSAT path.')
  failed = 1;
end
if exist(Path.themes) ~= 7
  disp('No "themes" folder found in the current PSAT path.')
  failed = ~mkdir(Path.psat,'themes');
  if failed,
    disp('The folder "themes" could not be created')
  end
end
if exist(Path.filters) ~= 7
  disp('No "filters" folder found in the current PSAT path.')
  failed = 1;
end
if failed
  disp('The PSAT folder seems to be uncomplete and needs to be installed again.')
  disp('PSAT initialization failed.')
  return
end

if ~command_line_psat
  hdl = fm_enter(psatdir(1).path,psatver,psatdate);
  pause(0.1)
  hold on
  % the size of the splash window are 346 x 410 pixels
  hdlp = patch(8+[0 0 33 33 0],410-[5 15 15 5 5], ...
               [0.50 0.45 0.95],'EdgeColor',[1 1 1]);
  set(hdlp,'EraseMode','none');
  idx = 33;
  pause(0.1)
else
  idx = 0;
end

File = struct('data','','pert','','temp','','modify','');
Source = struct('data','','pert','','description','');

jay = sqrt(-1);

History = struct('text','', ...
                 'index',1, ...
                 'string','', ...
                 'workspace',0, ...
                 'Max',500, ...
                 'FontName', ...
                 'Courier', ...
                 'FontSize',12, ...
                 'FontAngle','normal', ...
                 'FontWeight','normal', ...
                 'BackgroundColor','w', ...
                 'ForegroundColor','k');

History.text = {['PSAT version ',psatver, ...
		 ', Copyright (C) 2002-2006 Federico Milano']; ...
		' '; ['Session  ', datestr(now)]};

if ~command_line_psat
  set(hdlp,'XData',8+idx+[0 0 33 33 0])
  idx = idx + 33;
  pause(0.1)
end

try
  settings
catch
  Theme = struct('color01',[0.800 0.800 0.800], ...
                 'color02',[0.753 0.753 0.753], ...
                 'color03',[0.502 0.502 0.502], ...
                 'color04',[1.000 1.000 1.000], ...
                 'color05',[0.000 0.000 0.502], ...
                 'color06',[1.000 1.000 0.502], ...
                 'color07',[0.502 0.000 0.000], ...
                 'color08',[0.750 0.390 0.503], ...
                 'color09',[1.000 1.000 1.000], ...
                 'color10',[1.000 1.000 1.000], ...
                 'color11',[0.920 0.920 0.920], ...
                 'font01','Courier','hdl',zeros(18,1));

  Settings = struct('dlf',0, ...
                    'dac',0, ...
                    'pq2z',0, ...
                    'pv2pq',0, ...
                    'showlf',0, ...
                    'init', 0, ...
                    'status', 1, ...
                    'conv', 1, ...
                    'red',0, ...
                    'plot',0, ...
                    'plottype',1, ...
                    'method',2, ...
                    'show',1, ...
                    'vs',0, ...
                    'ok',0, ...
                    'pfsolver', 1, ...
                    'deltadelta',180, ...
                    'checkdelta',0, ...
                    'deltat',1e-5, ...
                    'deltatmax',1e-3, ...
                    'deltatmin',1e-5, ...
                    'chunk',100, ...
                    'maxvar',1500, ...
                    'mv',0, ...
                    'iter',0, ...
                    'static', 0, ...
                    'zoom','', ...
                    'freq', 50, ...
                    'beep',0, ...
                    'dyntol',1e-5, ...
                    'dynmit',20, ...
                    'lftol',1e-5, ...
                    'lfmit',20, ...
                    'lftime', 0, ...
                    't0',0, ...
                    'tf',30, ...
                    'mva',100, ...
                    'rad',2*3.14159265358979*50, ...
                    'distrsw', 0, ...
                    'refbus',0, ...
                    'color',[], ...
                    'fixt',0, ...
                    'tstep',0.001, ...
                    'xlabel','time (s)', ...
                    'locksnap',0, ...
                    'tviewer',[], ...
                    'absvalues', 'off', ...
                    'shuntvalues', 'off', ...
                    'violations', 'off', ...
		    'export','txt', ...
                    'noarrows',1, ...
                    'octave',0, ...
                    'local',1, ...
                    'simtd',0, ...
                    'robust',0);

end

Settings.color = [0      0      1.0000;
                  0      0.5000 0;
                  1.0000 0      0;
                  0      0.7500 0.7500;
                  0.7500 0      0.7500;
                  0.7500 0.7500 0;
                  0.2500 0.2500 0.2500];

if ~command_line_psat
  set(0,'DefaultUicontrolBackgroundColor',Theme.color01)
end
a = version;
Settings.hostver = str2num(a(1:3));
Settings.init = 0;
Settings.version = psatver;
Settings.date = psatdate;
if exist('OCTAVE_VERSION')
  Settings.octave = 1;
else
  Settings.octave = 0;
end

if ~isfield(Settings,'local')
  Settings.local = 1;
end

if isunix & strcmp(Settings.export,'xls')
  Settings.export = 'txt';
end

switch computer
 case 'PCWIN'
  Settings.tviewer = '!notepad ';
 otherwise
  Settings.tviewer = '!cat ';
end

Hdl = struct('hist',0, ...
             'text',0, ...
             'status',0, ...
             'frame',0, ...
             'bar',0, ...
             'pert','');

if Settings.hostver >= 6
  Hdl.pert = str2func('pert');
else
  Hdl.pert = 'pert';
end

Fig = struct('main',0, ...
             'plot',0, ...
             'update',0, ...
             'hist',0, ...
             'pset',0, ...
             'make',0, ...
             'lib',0, ...
             'comp',0, ...
             'cset',0, ...
             'sset',0, ...
             'xset',0, ...
             'eigen',0, ...
             'matrx',0, ...
             'theme',0, ...
             'simset',0, ...
             'setting',0, ...
             'snap',0, ...
             'stat',0, ...
             'pmu',0, ...
             'dir',0, ...
             'cpf',0, ...
             'snb',0, ...
             'about',0, ...
             'author',0, ...
             'opf',0, ...
             'line',0, ...
             'clock',0, ...
             'license',0, ...
	     'warranty',0, ...
             'tviewer',0, ...
             'gams',0, ...
             'uwpflow',0, ...
             'laprint',0, ...
             'plotsel',0);

% Component Structures
% clear Comp
Comp.funct = cell(5,1);
Comp.number = cell(5,1);
Comp.prop = [];
Comp.n = 0;
Comp.init = 0;
Comp.descr = '';
Comp.name = '';
Comp.shunt = 1;

if ~command_line_psat
  set(hdlp,'XData',8+idx+[0 0 33 33 0])
  idx = idx + 33;
  pause(0.1)
end

% Output Variables and Names Structures
Varout = struct('t',[],'vars',[],'idx',[]);
Varname = struct('bus','','comp','','uname','','fname','','nvars',0, ...
                 'uvars','','fvars','','fcomp','','format','', ...
                 'idx',[],'custom',0,'fixed',1,'x',1,'V',1,'PQ',0,'Pij',0);

% Snapshots
Snapshot = struct('name','', 'time',Settings.t0, 'V',[], ...
                  'ang',[], 'x', [], 'Y', [], 'Pg', [], ...
                  'Qg', [], 'Pl', [], 'Ql', [], 'vfd', [], ...
                  'pmech', [], 'Jlf', [], 'Jlfv', [], ...
                  'Fx', [], 'Fy', [], 'Gx', [], 'Ploss', [], ...
                  'Qloss', [], 'it', 0);

filemode = 'rt';

% Existing component variable names
fid = fopen([Path.psat,'namevar.ini'],filemode);
if fid == -1,
  disp('#Error: File "namevar.ini" cannot be open.')
  failed = 1;
else
  nname = 0;
  while 1
    sline = fgetl(fid);
    if ~ischar(sline), break; end
    try
      Varname.uname{nname+1,1} = deblank(sline(1:20));
      Varname.fname{nname+1,1} = deblank(sline(21:40));
      Varname.comp{nname+1,1}  = deblank(sline(41:end));
      if Settings.octave & strcmp(Varname.comp{nname+1,1}(end),'\r')
        Varname.comp{nname+1,1} = Varname.comp{nname+1,1}(1:end-1);
      end
      nname = nname + 1;
    catch
      % nothing to do ...
    end
  end
  count = fclose(fid);
end
fid = fopen([Path.psat,'history.ini'],filemode);
if fid == -1,
  disp('#Error: File "history.ini" cannot be open.')
  failed = 1;
else
  while 1
    sline = fgetl(fid);
    if ~ischar(sline), break, end
    try
      eval(['History.',sline(1:19),' = ',sline(20:end),';']);
    catch
      % nothing to do ...
    end
  end
  count = fclose(fid);
end

if ~command_line_psat
  set(hdlp,'XData',8+idx+[0 0 33 33 0])
  idx = idx + 33;
  pause(0.1)
end

% Existing component function names
fid = fopen([Path.psat,'fcomp.ini'],filemode);
if fid == -1,
  disp('#Error: File "fcomp.ini" cannot be open.')
  failed = 1;
else
  tipi = 0;
  while 1
    sline = fgetl(fid);
    if ~ischar(sline),
      break
    end
    try
      tipi = tipi + 1;
      Varname.fcomp{tipi,1} = deblank(sline(1:15));
      Varname.format{tipi,1} = sline(16:end);
      if Settings.octave & strcmp(Varname.format{tipi,1}(end),'\r')
        Varname.format{tipi,1} = Varname.format{tipi,1}(1:end-1);
      end
    catch
      % nothing to do ...
    end
  end
  fclose(fid);
end

if ~command_line_psat
  set(hdlp,'XData',8+idx+[0 0 33 33 0])
  idx = idx + 33;
  pause(0.1)
end

% Existing components functions
fid = fopen([Path.psat,'comp.ini'],filemode);
if fid == -1,
  disp('#Error: File "comp.ini" cannot be open.')
  failed = 1;
else
  ncomp=0;
  while 1
    sline = fgetl(fid);
    if ~ischar(sline), break; end
    try
      Comp.funct{ncomp+1,1} = deblank(sline(1:23));
      Comp.number{ncomp+1,1} = [deblank(sline(25:44)),'.n'];
      Comp.prop(ncomp+1,:) = str2num(sline(46:60));
      ncomp=ncomp+1;
    catch
      % nothing to do ...
    end
  end
  count = fclose(fid);
  Comp.funct{ncomp+1}  = 'fm_pv';
  Comp.number{ncomp+1} = 'PV.n';
  Comp.prop(ncomp+1,:) = [1 1 0 0 0 1 0 0];