consfmin.m

MATPOWER 一款基于MATLAB的电力系统潮流计算及优化的程序

function [g, geq, dg, dgeq] = consfmin(x, baseMVA, bus, gen, gencost, branch, areas, Ybus, Yf, Yt, mpopt, parms, ccost, N, fparm, H, Cw)
%CONSFMIN  Evaluates nonlinear constraints and their Jacobian for OPF.
%   [g, geq, dg, dgeq] = consfmin(x, baseMVA, bus, gen, gencost, ...
%                           branch, areas, Ybus, Yf, Yt, mpopt, parms, ccost, N, fparm, H, Cw)

%   MATPOWER
%   $Id: consfmin.m,v 1.11 2006/03/14 21:32:06 ray Exp $
%   by Carlos E. Murillo-Sanchez, PSERC Cornell & Universidad Autonoma de Manizales
%   and Ray Zimmerman, PSERC Cornell
%   Copyright (c) 1996-2006 by Power System Engineering Research Center (PSERC)
%   See http://www.pserc.cornell.edu/matpower/ for more info.

%%----- initialize -----
%% define named indices into data matrices
[GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
    MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
    QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
[F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
    TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
    ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
[PW_LINEAR, POLYNOMIAL, MODEL, STARTUP, SHUTDOWN, NCOST, COST] = idx_cost;

%% constant
j = sqrt(-1);

%% unpack needed parameters
nb = parms(1);
ng = parms(2);
nl = parms(3);
ny = parms(4);
% nx = parms(5);
% nvl = parms(6);
nz = parms(7);
% nxyz = parms(8);
thbas = parms(9);
thend = parms(10);
vbas = parms(11);
vend = parms(12);
pgbas = parms(13);
pgend = parms(14);
qgbas = parms(15);
qgend = parms(16);
% ybas = parms(17);
% yend = parms(18);
% zbas = parms(19);
% zend = parms(20);
% pmsmbas = parms(21);
% pmsmend = parms(22);
% qmsmbas = parms(23);
% qmsmend = parms(24);
% sfbas = parms(25);
% sfend = parms(26);
% stbas = parms(27);
% stend = parms(28);

%% grab Pg & Qg
Pg = x(pgbas:pgend);            %% active generation in p.u.
Qg = x(qgbas:qgend);            %% reactive generation in p.u.

%% put Pg & Qg back in gen
gen(:, PG) = Pg * baseMVA;      %% active generation in MW
gen(:, QG) = Qg * baseMVA;      %% reactive generation in MVAr
 
%% rebuild Sbus
Sbus = makeSbus(baseMVA, bus, gen); %% net injected power in p.u.

%% ----- evaluate constraints -----
%% reconstruct V
Va = zeros(nb, 1);
Va = x(thbas:thend);
Vm = x(vbas:vend);
V = Vm .* exp(j * Va);

%% evaluate power flow equations
mis = V .* conj(Ybus * V) - Sbus;

%% compute branch power flows
Sf = V(branch(:, F_BUS)) .* conj(Yf * V);   %% complex power injected at "from" bus (p.u.)
St = V(branch(:, T_BUS)) .* conj(Yt * V);   %% complex power injected at "to" bus (p.u.)

%%----- evaluate constraint function values -----
%% first, the equality constraints (power flow)
geq = [ real(mis);              %% active power mismatch for all buses
        imag(mis) ];            %% reactive power mismatch for all buses

%% then, the inequality constraints (branch limits)
if mpopt(24) == 1   %% P limit (Pan Wei)
  g = [ real(Sf) - branch(:, RATE_A)/baseMVA;   %% branch real power limits (from bus)
        real(St) - branch(:, RATE_A)/baseMVA ]; %% branch real power limits (to bus)
else                %% |S| limit
  g = [ abs(Sf) - branch(:, RATE_A)/baseMVA;    %% branch apparent power limits (from bus)
        abs(St) - branch(:, RATE_A)/baseMVA ];  %% branch apparent power limits (to bus)
end

%%----- evaluate partials of constraints -----
if nargout > 2
  %% compute partials of injected bus powers
  [dSbus_dVm, dSbus_dVa] = dSbus_dV(Ybus, V);               %% w.r.t. V
  dSbus_dPg = sparse(gen(:, GEN_BUS), 1:ng, -1, nb, ng);    %% w.r.t. Pg
  dSbus_dQg = sparse(gen(:, GEN_BUS), 1:ng, -j, nb, ng);    %% w.r.t. Qg
  
  %% compute partials of line flows w.r.t. V
  [dSf_dVa, dSf_dVm, dSt_dVa, dSt_dVm, Sf, St] = dSbr_dV(branch, Yf, Yt, V);

  %% line limits are w.r.t apparent power, so compute partials of apparent power
  [dAf_dVa, dAf_dVm, dAt_dVa, dAt_dVm] = ...
            dAbr_dV(dSf_dVa, dSf_dVm, dSt_dVa, dSt_dVm, Sf, St);

  %% construct Jacobian of equality constraints (power flow) and transpose it
  dgeq = [
    %% equality constraints
    real(dSbus_dVa), real(dSbus_dVm), ...
          real(dSbus_dPg), real(dSbus_dQg), sparse(nb, ny+nz);  %% P mismatch
    imag(dSbus_dVa), imag(dSbus_dVm), ...
          imag(dSbus_dPg), imag(dSbus_dQg), sparse(nb, ny+nz);  %% Q mismatch
   ]'; 
  %% construct Jacobian of inequality constraints (branch limits)
  %% and transpose it so fmincon likes it
  if mpopt(24) == 1     %% P limit (Pan Wei)
    dg = [
      real(dSf_dVa), real(dSf_dVm), sparse(nl,2*ng+ny+nz);  %% Pf limit
      real(dSt_dVa), real(dSt_dVm), sparse(nl,2*ng+ny+nz);  %% Pt limit
    ]';
  else                  %% |S| limit
    dg = [
      dAf_dVa, dAf_dVm, sparse(nl,2*ng+ny+nz);  %% |Sf| limit
      dAt_dVa, dAt_dVm, sparse(nl,2*ng+ny+nz);  %% |St| limit
    ]';
  end  

  %% the following lines can be removed if/when fmincon
  %% supports sparse matrices for non-linearly constrained problems
  dgeq = full(dgeq);
  dg = full(dg);
end

return;