Hi Jarek,
I don’t think any code changes should be necessary for solving the problem you
want to solve. That is, if you put the constant impedance loads in the GS and
BS columns of the bus matrix, the calculations should be correct, including all
of the resulting voltages and flows. The only issue is that the power delivered
to these elements is not reported anywhere as “load”, so you’ll have to do that
manually after the fact, but that is a trivial calculation once you have the
bus voltages. For example, if the nominal (i.e. at 1 p.u. voltage) active and
reactive constant impedance loads are given in the nb x 1 vectors Pdz0 and
Qdz0, then you can calculate the constant impedance loads at the solution as
follows ...
define_constants;
mpc = loadcase(...);
mpc.bus(:, GS) = mpc.bus(:, GS) + Pdz0;
mpc.bus(:, BS) = mpc.bus(:, BS) - Qdz0;
r = runpf(mpc);
Pdz = Pdz0 ./ r.bus(:, VM) .^ 2;
Qdz = Qdz0 ./ r.bus(:, VM) .^ 2;
The constant power loads can be found using total_load() …
[Pds, Qds] = total_load(r.bus, r.gen, 'bus');
So, no need for any of your steps beyond a) under #2, except I suppose f) could
be done by simply doing …
r.bus(:, PD) = r.bus(:, PD) + Pdz;
r.bus(:, QD) = r.bus(:, QD) + Qdz;
Ray
> On Apr 5, 2016, at 1:43 AM, Jaroslaw Krata <j.kr...@uq.edu.au> wrote:
>
> Dear Matpower authors,
>
> I wish to use the Matpower software in my PhD research.
> I need to calculate power flow for the distribution system that consists of
> many constant impedance and constant power loads.
>
> After I read in the manual that " (...) the constant impedance portions can
> be modeled as a shunt element (...)",
> I thought naively that it's gonna be easy ;-).
> I went through other messages from this mailing list and I understand that I
> need to changes the code of the Power Flow calculation.
>
> In particular, I was trying use this advise
> https://www.mail-archive.com/matpower-l%40cornell.edu/msg03254.html
> <https://www.mail-archive.com/matpower-l%40cornell.edu/msg03254.html>
> but it seems that since 2014 the code of Matpower has changed and this advice
> is not easy applicable to the newest version of the Matpower.
>
> Thus, I would like to ask you for help of the required code changes in the
> Matpower 5.1.
> I wish to achieve two targets:
> 1. correct power flow calculation assuming const P and const Z loads
> 2. output branch and power matrices showing updated (real) values of power
> from const Z loads.
>
> I started as follows:
> a) I assume that whole constant power loads are assigned to the bus matrix -
> Pd and Qd columns
> while all constant impedance loads are assigned to Gs and Bs columns
>
> b) I modified makeSbus() function by adding lines related to const P and Z
> line1: [Sbus, Sbus_gen, Sbus_constP, Sbus_constZ] = makeSbus(baseMVA, bus,
> gen) %adding Sbus_gen, Sbus_consP and Sbus_constZ vectors
> line36 to 38:
> Sbus = ( Cg * (gen(on, PG) + 1j * gen(on, QG)) ... %% power injected by
> generators
> - (bus(:, PD) + 1j * bus(:, QD)) ... %% plus power
> injected by const P loads
> ) / baseMVA; %% converted to p.u.
>
> replaced by:
> Sbus = ( Cg * (gen(on, PG) + 1j * gen(on, QG)) ... %% power injected by
> generators
> - (bus(:, PD) + 1j * bus(:, QD)) ... %% plus power injected
> by const P loads
> - (bus(:, GS) + 1j * bus(:, BS)) )/ ... %% plus power injected
> by const Z loads
> baseMVA; %% converted to p.u.
>
> Sbus_gen = (Cg * (gen(on, PG) + 1j * gen(on, QG))) / baseMVA;
> Sbus_constP = (bus(:, PD) + 1j * bus(:, QD)) / baseMVA;
> Sbus_constZ = (bus(:, GS) + 1j * bus(:, BS)) / baseMVA;
>
> c) I pass all Sbus vectors to newtonpf() function
> [V, converged, i] = newtonpf(Ybus, Sbus, Sbus_gen, Sbus_constP, Sbus_constZ,
> V0, ref, pv, pq, mpopt)
>
> d) In newtonpf(), I calculate 'mis' variable as
> line 54 and line 106: mis = V .* conj(Ybus * V) - (Sbus_gen - Sbus_constP -
> Sbus_constZ.*(Vm.^2))
> (replacing the default mis = V .* conj(Ybus * V) - Sbus;)
>
> e) How should I change the dSbus_dV() function??
> f) How to update bus matrix to show updated value of const Z power?
> g) Should I change other functions as well (e.g. makeYbus() function)?
>
>
> I will be grateful for your help and suggestions,
> BR,
> Jarek
>
>
