genform.m

该程序是计算最优潮流的matlab工具箱。可以很好的求解目标函数不同的最优潮流问题。

%GENFORM: Help file describing the generalized OPF formulation used by
%         the fmincopf and MINOPF solvers.
%
%*-*-*-*-*
%CONTENTS
%*-*-*-*-*
%
%1: General OPF problem formulation
%2: Problem data transformations & general linear restrictions
%3: Piecewise linear convex cost formulation using constrained cost variables
%4: Dispatchable loads
%5: Additional information on the generalized formulation structure
%
%
%*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%1. GENERAL OPF PROBLEM FORMULATION
%*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
%The generalized OPF formulation used by the fmincon and MINOPF solvers
%can be written as follows:
%
%   min  sum( f1i(Pgi) + f2i(Qgi) ) + c * [x; y; z]
%  x,y,z
%
%subject to
%  g(x) <=> 0          (nonlinear constraints: bus power balance equations
%                                              & branch flow limits)
%  l <= A*[x;y;z] <= u (general linear constraints)
%  xmin <= x <= xmax   (variable bounds: voltage limits, generation limits)
%
%where x = [ Theta; V; Pg; Qg ], and (Pg, Qg) are the ng active and ng
%reactive injections into the network and (V,Theta) are the voltage
%magnitudes and voltage phase angles at each of the nb buses in the
%network. The objective is a generator-wise additive function f(Pg,Qg)
%plus a linear term, including (y, z) which are described below. The
%nonlinear constraints g() are grouped as
%
%  g(x) = [ gp(x);    % nb bus real power mismatch constraints
%           gq(x);    % nb bus reactive power mismatch constraints
%           gsf(x);   % nl branch MVA injection upper limits at "from" end
%           gst(x) ]; % nl branch MVA injection upper limits at "to" end
%
%The additional variables y and z can optionally augment the problem to
%allow casting more general constraints and is also the mechanism though
%which piecewise linear convex costs are modeled. The linear cost vector
%c allows the addition of a general linear cost on the whole set [x;y;z]
%of variables. With this setup, mixed polynomial/piecewise linear costs
%can be dealt with, as well as dispatchable or price-sensitive loads
%requiring a constant power factor.
%
%
%*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%2. PROBLEM DATA TRANSFORMATIONS AND GENERAL LINEAR RESTRICTIONS
%*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
%To add general linear constraints of ones own, it is necessary to
%understand the standard transformations performed on the input data
%(bus, gen, branch, areas and gencost tables) before the problem is
%solved in order to know where the optimization variables end up in the x
%vector. All of these transformations are reversed after solving the
%problem so that output data is in the right place in the tables.
%
%The first step filters out inactive generators and branches; original
%tables are saved for data output.
%
%  comgen   = find(gen(:,GEN_STATUS) > 0);     % find online generators
%  onbranch = find(branch(:,BR_STATUS) ~= 0);  % find online branches
%  gen      = gen(comgen, :);
%  branch   = branch(onbranch, :);
%
%The second step is a renumbering of the bus numbers in the bus table so
%that the resulting table contains consecutively-numbered buses starting
%from 1:
%
%  [i2e, bus, gen, branch, areas] = ext2int(bus, gen, branch, areas);
%
%where i2e is saved for inverse reordering at the end. Finally,
%generators are further reordered by bus number:
%
%  ng = size(gen,1);                 % number of generators or injections
%  [tmp, igen] = sort(gen(:, GEN_BUS));
%  [tmp, inv_gen_ord] = sort(igen);  % save for inverse reordering at the end
%  gen  = gen(igen, :);
%  if ng == size(gencost,1)          % This is because gencost might have
%    gencost = gencost(igen, :);     % twice as many rows as gen if there
%  else                              % are reactive injection costs.
%    gencost = gencost( [igen; igen+ng], :);
%  end
%
%Having done this, the variables inside the x vector now have the same
%ordering as in the bus, gen tables:
%
%  x = [ Theta ;    % nb bus voltage angles
%          V   ;    % nb bus voltage magnitudes
%          Pg  ;    % ng active power injections (p.u.) (ascending bus order)
%          Qg ];    % ng reactive power injections (p.u.)(ascending bus order)
%
%and the nonlinear constraints have the same order as in the bus, branch
%tables
%
%  g = [ gp;        % nb real power flow mismatches (p.u.)
%        gq;        % nb reactive power flow mismatches (p.u.)
%        gsf;       % nl "from" end apparent power injection limits (p.u.)
%        gst ];     % nl "to" end apparent power injection limits (p.u.)
%
%With this setup, box bounds on the variables are applied as follows: the
%reference angle is bounded above and below with the value specified for
%it in the original bus table. The V section of x is bounded above and
%below with the corresponding values for VMAX and VMIN in the bus table.
%The Pg and Qg sections of x are bounded above and below with the
%corresponding values for PMAX, PMIN, QMAX and QMIN in the gen table. The
%nonlinear constraints are similarly setup so that gp and gq are equality
%constraints (zero RHS) and the limits for gsf, gst are taken from the
%RATE_A column in the branch table.
%
%The following example illustrates how an additional general linear
%constraint can be added to the problem formulation. In the standard
%solution to case9.m, the voltage angle for bus 7 lags the voltage angle
%in bus 2 by 6.09 degrees. Suppose we want to limit that lag to 5 degrees
%at the most. A linear restriction of the form
%
%  Theta(2) - Theta(7) <= 5 degrees
%
%would do the trick. We have nb = 9 buses, ng = 3 generators and nl = 9
%branches. Therefore the first 9 elements of x are bus voltage angles,
%elements 10:18 are bus voltage magnitudes, elements 19:21 are active
%injections corresponding to the generators in buses 1, 2 and 3 (in that
%order) and elements 22:24 are the corresponding reactive injections.
%Note that in this case the generators in the original data already
%appear in ascending bus order, so no permutation with respect to the
%original data is necessary. Going back to the angle restriction, we see
%that it can be cast as
%
%  [ 0 1 0 0 0 0 -1 0 0 zeros(1,nb+ng+ng) ] * x  <= 5 degrees
%
%We can set up the problem as follows:
%
%  A = sparse([1;1], [2;7], [1;-1], 1, 24);
%  l = -Inf;
%  u = 5 * pi/180;
%  mpopt = mpoption('OPF_ALG', 520); % use fmincon w/generalized formulation
%  opf('case9', A, l, u, mpopt)
%
%which indeed restricts the angular separation to 5 degrees. NOTE: In
%this example, the total number of variables is 24, but if there are any
%piecewise linear cost functions, there may be additional "helper"
%variables used by the solver and in that case the number of columns in A
%may need to be larger. The next section describes how this is done. If
%all costs are polynomial, however, no extra variables are needed.
%
%
%*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%3. PIECEWISE LINEAR CONVEX COST FORMULATION USING CONSTRAINED COST VARIABLES
%*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%
%The generalized formulation allows for an easy way to model any
%piecewise linear costs. Such a cost curve looks like
%
% cost axis
%     ^
%c2   |                           *
%     |                        *
%     |                     *
%c1   |                  *
%     |          *
%c0   |   *
%     |------------------------------->   Pg axis
%         x0            x1         x2
%        Pmin                     Pmax
%
%This non-differentiable cost can be modeled using one helper cost
%variable for each such cost curve and additional restrictions on this
%variable and Pg, one restriction for each segment of the curve. The
%restrictions build a convex "basin" and they are equivalent to saying
%that the cost variable must lie in the epigraph of the cost curve. When
%the cost is minimized, the cost variable will be pushed against this
%basin. If the helper cost variable is y, then the contribution of the
%generators' cost to the total cost is exactly y, and in the above case
%the two restrictions needed are
%
%  1)  y >= m1*(Pg - x0) + c0   (y must lie above the first segment)
%
%  2)  y >= m2*(Pg - x1) + c1   (y must lie above the second segment)
%
%where, m1 and m2 are the slopes of the two segments. Also needed, of
%course, are the box restrictions on Pg:  Pmin <= Pg <= Pmax. The