function [Bbus, Bf, Pbusinj, Pfinj] = makeBdc(baseMVA, bus, branch)%MAKEBDC   Builds the B matrices and phase shift injections for DC power flow.%   [Bbus, Bf, Pbusinj, Pfinj] = makeBdc(baseMVA, bus, branch) returns the%   B matrices and phase shift injection vectors needed for a DC power flow.%   The bus real power injections are related to bus voltage angles by%       P = Bbus * Va + Pbusinj%   The real power flows at the from end the lines are related to the bus%   voltage angles by%       Pf = Bf * Va + Pfinj%   Does appropriate conversions to p.u.%   MATPOWER%   $Id: makeBdc.m,v 1.10 2007/06/26 16:16:17 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.%% constantsj = sqrt(-1);nb = size(bus, 1);          %% number of busesnl = size(branch, 1);       %% number of lines%% define named indices into bus, branch matrices[PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...    VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;[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;%% check that bus numbers are equal to indices to bus (one set of bus numbers)if any(bus(:, BUS_I) ~= [1:nb]')    error('makeBdc: buses must be numbered consecutively in bus matrix')end%% for each branch, compute the elements of the branch B matrix and the phase%% shift "quiescent" injections, where%%%%      | Pf |   | Bff  Bft |   | Vaf |   | Pfinj |%%      |    | = |          | * |     | + |       |%%      | Pt |   | Btf  Btt |   | Vat |   | Ptinj |%%stat = branch(:, BR_STATUS);                    %% ones at in-service branchesb = stat ./ branch(:, BR_X);                    %% series susceptancetap = ones(nl, 1);                              %% default tap ratio = 1i = find(branch(:, TAP));                       %% indices of non-zero tap ratiostap(i) = branch(i, TAP);                        %% assign non-zero tap ratiosb = b ./ tap;%% build Bbusf = branch(:, F_BUS);                           %% list of "from" busest = branch(:, T_BUS);                           %% list of "to" busesCf = sparse(f, 1:nl, ones(nl, 1), nb, nl);      %% connection matrix for line & from busesCt = sparse(t, 1:nl, ones(nl, 1), nb, nl);      %% connection matrix for line & to busesBbus =  Cf * spdiags(b, 0, nl, nl) * Cf' + ...  %% Bff term of branch admittance        Cf * spdiags(-b, 0, nl, nl) * Ct' + ... %% Bft term of branch admittance        Ct * spdiags(-b, 0, nl, nl) * Cf' + ... %% Btf term of branch admittance        Ct * spdiags(b, 0, nl, nl) * Ct';       %% Btt term of branch admittance%% build phase shift injection vectorsPfinj = b .* (-branch(:, SHIFT) * pi/180);      %% injected at the from bus ...    % Ptinj = -Pfinj;                           %% ... and extracted at the to busPbusinj = (Cf - Ct) * Pfinj;                    %% Pbusinj = Cf * Pfinj + Ct * Ptinj;%% Build Bf such that Bf * Va is the vector of real branch powers injected%% at each branch's "from" busif nargout > 1    i = [[1:nl]'; [1:nl]'];     %% double set of row indices        Bf = sparse(i, [f; t], [b; -b]);endreturn;