function Va = dcpf(B, Pbus, Va0, ref, pv, pq)%DCPF  Solves a DC power flow.%   [Va, success] = dcpf(B, Pbus, Va0, ref, pv, pq) solves for the bus%   voltage angles at all but the reference bus, given the full system%   B matrix and the vector of bus real power injections, the initial%   vector of bus voltage angles (in radians), and column vectors with%   the lists of bus indices for the swing bus, PV buses, and PQ buses,%   respectively. Returns a vector of bus voltage angles in radians.%   MATPOWER%   $Id: dcpf.m,v 1.4 2004/08/23 20:56:12 ray Exp $%   by Carlos E. Murillo-Sanchez, PSERC Cornell & Universidad Autonoma de Manizales%   and Ray Zimmerman, PSERC Cornell%   Copyright (c) 1996-2004 by Power System Engineering Research Center (PSERC)%   See http://www.pserc.cornell.edu/matpower/ for more info.%% initialize result vectorVa = Va0;%% Since pulling off rows of a large sparse matrix like B%% can be slow in Matlab 5, we do this ...temp = B(:, [pv; pq])';B_pvpq_rows = temp(:, [pv; pq])';temp = B(:, ref)';B_ref_row   = temp(:, [pv; pq])';Va([pv; pq]) = B_pvpq_rows \ (Pbus([pv; pq]) - B_ref_row * Va0(ref));%% ... instead of this ...% Va([pv; pq]) = B([pv; pq], [pv; pq]) \ (Pbus([pv; pq]) - B([pv; pq], ref) * Va0(ref));return;