//MAKESBUS   Builds the vector of complex bus power injections.
//   SBUS = MAKESBUS(BASEMVA, BUS, GEN)
//   SBUS = MAKESBUS(BASEMVA, BUS, GEN, MPOPT, VM)
//   SBUS = MAKESBUS(BASEMVA, BUS, GEN, MPOPT, VM, SG)
//   returns the vector of complex bus power injections, that is, generation
//   minus load. Power is expressed in per unit. If the MPOPT and VM arguments
//   are present it evaluates any ZIP loads based on the provided voltage
//   magnitude vector. If VM is empty, it assumes nominal voltage. If SG is
//   provided, it is a complex ng x 1 vector of generator power injections in
//   p.u., and overrides the PG and QG columns in GEN, using GEN only for
//   connectivity information.
//
//   [SBUS, DSBUS_DVM] = MAKESBUS(BASEMVA, BUS, GEN, MPOPT, VM)
//   With two output arguments, it computes the partial derivative of the
//   bus injections with respect to voltage magnitude, leaving the first
//   return value SBUS empty. If VM is empty, it assumes no voltage dependence
//   and returns a sparse zero matrix.
//
//   See also MAKEYBUS.

//   MATPOWER
//   Copyright (c) 1996-2016, Power Systems Engineering Research Center (PSERC)
//   by Ray Zimmerman, PSERC Cornell
//
//   This file is part of MATPOWER.
//   Covered by the 3-clause BSD License (see LICENSE file for details).
//   See https://matpower.org for more info.

// define named indices into bus, gen matrices

fn make_sbus(baseMVA, bus, gen) {
    nb = size(bus, 0);
    Sd = make_sdzip(baseMVA, bus, []);
    on = find(slice(gen, [0], [GEN_STATUS-1, GEN_STATUS]) > 0);  // which generators are on?
    gbus = get_multi(slice(gen, [0], [GEN_BUS-1, GEN_BUS]), on) - 1; // what buses are they at?
    ngon = size(on, 0);
    Cg = full(sparse(gbus, range(0, ngon), 1, nb, ngon));
    Gp = slice(gen, [0], [PG-1, PG]);
    Gq = slice(gen, [0], [QG-1, QG]);
    Sbusg = Cg * (get_multi(Gp, on)  + c(0,1) * get_multi(Gq, on)) / baseMVA;
    z = slice(Sd, [0], [0,1]);
    i = slice(Sd, [0], [1,2]);
    p = slice(Sd, [0], [2,3]);
    Vm = ones(nb, 1);
    Sbusd = p + i .* Vm + z .* Vm.^2;
    // complex power injection at each bus
    return Sbusg' - Sbusd;
}