makePTDF

MAKEPTDF   Builds the DC PTDF matrix for a given choice of slack.
   H = MAKEPTDF(BASEMVA, BUS, BRANCH, SLACK) returns the DC PTDF
   matrix for a given choice of slack. The matrix is nbr x nb, where
   nbr is the number of branches and nb is the number of buses. The SLACK
   can be a scalar (single slack bus) or an nb x 1 column vector of
   weights specifying the proportion of the slack taken up at each bus.
   If the SLACK is not specified the reference bus is used by default.

   Examples:
       H = makePTDF(baseMVA, bus, branch);
       H = makePTDF(baseMVA, bus, branch, 1);
       slack = rand(size(bus, 1), 1);
       H = makePTDF(baseMVA, bus, branch, slack);

   See also MAKELODF.

0001 function H = makePTDF(baseMVA, bus, branch, slack)
0002 %MAKEPTDF   Builds the DC PTDF matrix for a given choice of slack.
0003 %   H = MAKEPTDF(BASEMVA, BUS, BRANCH, SLACK) returns the DC PTDF
0004 %   matrix for a given choice of slack. The matrix is nbr x nb, where
0005 %   nbr is the number of branches and nb is the number of buses. The SLACK
0006 %   can be a scalar (single slack bus) or an nb x 1 column vector of
0007 %   weights specifying the proportion of the slack taken up at each bus.
0008 %   If the SLACK is not specified the reference bus is used by default.
0009 %
0010 %   Examples:
0011 %       H = makePTDF(baseMVA, bus, branch);
0012 %       H = makePTDF(baseMVA, bus, branch, 1);
0013 %       slack = rand(size(bus, 1), 1);
0014 %       H = makePTDF(baseMVA, bus, branch, slack);
0015 %
0016 %   See also MAKELODF.
0017 
0018 %   For convenience, SLACK can also be an nb x nb matrix, where each
0019 %   column specifies how the slack should be handled for injections
0020 %   at that bus.
0021 
0022 %   MATPOWER
0023 %   Copyright (c) 2006-2015 by Power System Engineering Research Center (PSERC)
0024 %   by Ray Zimmerman, PSERC Cornell
0025 %
0026 %   $Id: makePTDF.m 2644 2015-03-11 19:34:22Z ray $
0027 %
0028 %   This file is part of MATPOWER.
0029 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0030 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0031 
0032 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0033     VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0034 
0035 %% use reference bus for slack by default
0036 if nargin < 4
0037     slack = find(bus(:, BUS_TYPE) == REF);
0038     slack = slack(1);
0039 end
0040 
0041 %% set the slack bus to be used to compute initial PTDF
0042 if length(slack) == 1
0043     slack_bus = slack;
0044 else
0045     slack_bus = 1;      %% use bus 1 for temp slack bus
0046 end
0047 
0048 nb = size(bus, 1);
0049 nbr = size(branch, 1);
0050 noref   = (2:nb)';      %% use bus 1 for voltage angle reference
0051 noslack = find((1:nb)' ~= slack_bus);
0052 
0053 %% check that bus numbers are equal to indices to bus (one set of bus numbers)
0054 if any(bus(:, BUS_I) ~= (1:nb)')
0055     error('makePTDF: buses must be numbered consecutively in bus matrix')
0056 end
0057 
0058 %% compute PTDF for single slack_bus
0059 [Bbus, Bf, Pbusinj, Pfinj] = makeBdc(baseMVA, bus, branch);
0060 H = zeros(nbr, nb);
0061 H(:, noslack) = full(Bf(:, noref) / Bbus(noslack, noref));
0062         %%    = full(Bf(:, noref) * inv(Bbus(noslack, noref)));
0063 
0064 %% distribute slack, if requested
0065 if length(slack) ~= 1
0066     if size(slack, 2) == 1  %% slack is a vector of weights
0067         slack = slack/sum(slack);   %% normalize weights
0068         
0069         %% conceptually, we want to do ...
0070         %%    H = H * (eye(nb,nb) - slack * ones(1, nb));
0071         %% ... we just do it more efficiently
0072         v = H * slack;
0073         for k = 1:nb
0074             H(:, k) = H(:, k) - v;
0075         end
0076     else
0077         H = H * slack;
0078     end
0079 end