Re: Convertion of CIM CGMES format into a MATPOWER case struct

[Jose Luis Marín]

I know you want to convert in the direction CIM CGMES ==> MATPOWER, but
just in case, the PowSyBl project has converters in the other direction:

*
 
https://nam12.safelinks.protection.outlook.com/?url=https%3A%2F%2Fgithub.com%2Fpowsybl%2Fpowsybl.github.io%2Fblob%2Fmain%2Fpages%2Fdocumentation%2Fdeveloper%2Frepositories%2Fpypowsybl.md%23features=05%7C01%7Cmatpower-l%40list.cornell.edu%7Caf35730a97bb41ea393b08dbe1b90a0a%7C5d7e43661b9b45cf8e79b14b27df46e1%7C0%7C0%7C638351957079434662%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C%7C%7C=8mGE74%2FPn5sF%2B3zCANWCsTORChJoj3tXWDfmXgIrYTY%3D=0
*
*   Grid modelling*

   - *We can load a network from a file. The supported formats are for the
   moment CGMES, MATPOWER, IEEE-CDF, PSS/E, UCTE and XIIDM.*
   - *We can save a network to a file. The supported formats are for the
   moment CGMES, UCTE, and XIIDM.*


-- 
Jose Luis Marin
Grupo AIA


El vie, 10 nov 2023 a las 7:10,  escribió:

> Dear Ray,
>
>
>
> thank you very much for your prompt feedback. Consequently, we will have
> to find another way.
>
>
>
> Best regards,
>
> Oliver
>
>
>
>
>
> *Von:* bounce-127879544-95614...@list.cornell.edu <
> bounce-127879544-95614...@list.cornell.edu> *Im Auftrag von *Ray Daniel
> Zimmerman
> *Gesendet:* Donnerstag, 9. November 2023 17:25
> *An:* MATPOWER-L 
> *Betreff:* Re: Convertion of CIM CGMES format into a MATPOWER case struct
>
>
>
> Unfortunately, I am not aware of any such converter and I am not in a
> position to create one. However, if someone else is, I’d be happy to
> include it in the MATPOWER distribution.
>
>
>
> Ray
>
>
>
> On Nov 9, 2023, at 1:43 AM, oliver.zehn...@elcom.admin.ch wrote:
>
>
>
> Dear all,
>
>
>
> ENTSO-E, the European Network of Transmission System Operators for
> Electricity, publishes datasets representing the power systems of the
> ENTSO-E members in the CIM CGMES format. Is there a converter available for
> MATPOWER, to convert data in the CIM CGMES format into a MATPOWER case
> struct (like e.g. psse2mpc)? And if not, are there any plans to build such
> a converter?
>
>
>
> Many thanks in advance for your feedback and best regards,
>
> Oliver Zehnder
>
>
>

Re: DC/DC converter model

[Jose Luis Marín]

I don't think DC grids can be accommodated in MATPOWER (only DC links can,
as far as I know).  But I'll try to answer your question because we
encountered this problem as well, a few years back.  The context was
steady-state power-flows in pure-DC microgrids.  What we found at the time
was that there seems to be no standardized modeling for the steady state of
DC/DC converters in the high-power, high-voltage domain.  But it's not
difficult to come up with one; we did so using this reference:

R. W. Erickson and D. Maksimovic, Fundamentals of Power Electronics (Second
Ed.), Kluwer Academic, 2004.
(see Chapters 2 & 3)

For our purposes, we did not need any more detailed modeling than that.
Your case may be different, I don't know.  Hope it helps.

-- 
Jose L. Marin
Grupo AIA



El vie., 15 mar. 2019 a las 10:35, Asma DABBABI (<
asma.dabb...@etu.univ-nantes.fr>) escribió:

> Hello everyone,
>
> I am looking for the power flow example of CIGRE B4 DC grid (see figure
> attached) to see the modeling of DC / DC converters. I need to do the
> power flow of an offshore wind farm connected all in DC on Mat AC/DC. In
> Mat AC/DC  i only find the modeling of AC / DC converter.
> Can you help me please to find how to do power flow in DC/DC converter
> case.
>
> Thank you
>
> Best regards

Re: negative line parameters from psse_convert_xfmr

[Jose Luis Marín]

I remember I did check extensively the conversion of two-winding
transformers a few years back; however, I did not check 3-winding ones.

You can check those formulas in this article by Siemens:

"Modeling of Three-Winding Voltage Regulating Transformers for Positive
Sequence Load Flow Analysis in PSS/E"
https://w3.usa.siemens.com/datapool/us/SmartGrid/docs/pti/2010July/PDFS/Modeling%20of%20Three%20Winding%20Voltage%20Regulating%20Transformers.pdf

-- 
Jose L. Marin
Grupo AIA



2018-08-31 0:12 GMT+02:00 Ray Zimmerman :

> This code was contributed by Yujia Zhu when he was at ASU working with Dan
> Tylavsky. I’m afraid I’m not familiar with the details of that code, but
> I’m pretty sure it was double-checked extensively against PSS/E’s internal
> models.
>
> Ray
>
>
> On Aug 6, 2018, at 4:50 AM, Aebi Stefanie 
> wrote:
>
> Dear all
>
>
> I am trying to convert a PSSE .raw network case to mpc format.
>
>
> I just discovered that I obtain negative branch parameters for some of the
> 3-winding transformers in my System, i.e. x < 0, r < 0. I tried to track
> the phenomenon and found the following section in psse_convert_xfmr.m:
>
>
> R1 = (R12+R31-R23) ./ 2;
> R2 = (R12+R23-R31) ./ 2;
> R3 = (R31+R23-R12) ./ 2;
> X1 = (X12+X31-X23) ./ 2;
> X2 = (X12+X23-X31) ./ 2;
> X3 = (X31+X23-X12) ./ 2;
>
>
> I've got a number of branches where the negative term (e.g. R23 in
> calculation of R1) is the biggest one, such that the resulting R or X is
> negative.
>
>
> Where does this definition come from?
> And is there a way I can avoid getting negative line parameters?
>
>
> Thanks for your help.
>
>
> Kind regards,
> Stefanie Aebi
>
>
>

Re: AC case network did not converge

[Jose Luis Marín]

Hello Asma,


I thank you very much for your help, I understood the solution that you
> proposed.  I decreased the cos phi to finally have a max_lambda> 1, after
> I saw the
> value of Q which corresponds to the value of max lambda > 1 and I put the
> Q for each turbine in negative like you did for P.
>

Please keep in mind that you probably cannot choose cos_phi freely for wind
turbines.  You will have to investigate what kind of technology
(inverters?) they're using, and then find out how much reactive power they
can really absorb/inject.  I'm not an expert in such models, but I bet
their operating range is not much wider than cos_phi 0.95 either way
(lagging/leading).

Additionally, if you trust your values of P to be correct, I think the only
way you can "torture your case until it solves" is by lowering the line
reactances.



> my questions are:
>
> 1/ why did you put the P and Q in negatives because I try to put them in
> positive and the calculation diverges,
>

Because the columns bus(:, [PD QD]) are power injection *demand* (see
http://www.pserc.cornell.edu//matpower/docs/ref/matpower6.0/idx_bus.html).
So they're supposed to be positive when they're loads, and negative when
they're generation.



> 2/ with max_lambda> 1 I got also mean voltage> Vmax = 1.1, is that
> logical, so I increased Vmax to 1.5 and it did not change the results
>

The columns bus(:, [VMAX VMIN]) play absolutely no role in *powerflow*;
they're only used for OPF.



> 3/ in my current case (see image and mfile) I have 2 transmission lines
> between 1002 and 3001 and between 1001 and 3001, my problem is when I
> increase the transmission distance (beyond 130 km) so I  increase the
> values of r x and b between these branches, the calculation diverges
> knowing that it worked for distances <130 km. and I also noticed that the
> Q  losses decrease by increasing the transmission distance which is not
> logical. what is the problem in your opinion?
>
> the solution that i did: I changed the value of system MVA base = 8MW *
> 100 (which was set to 100) so here my program converges without even
> injecting reactive power for the turbines and the value of Q increases by
> increasing the transmission distance.
> I thought to do that because I saw in other load flow calculation code
> that the basic value of the system = the production of each wind turbine
> (8MW) * nb of wind turbines (100) is my logic true and what is the effect
> of the base MVA value?
>

Increasing the values of X (while maintining everything else fixed) will
eventually lead to an infeasible case, that's no surprise. But the value of
net reactive losses (I^2 X) may vary in non-monotonic ways because it also
depends on other factors such as the line susceptances (Bshunt).

Now, did you change the MVAbase without making the corresponding change in
all of your R, X, B parameters?  If so, this is simply wrong. The MVA base
is just a choice of units for power; so if your case is infeasible in one
choice of units, it will still be infeasible in any other.  I'm guessing
you're obtaining convergence simply because the values of R, X, B are not
expressed correctly in the new MVAbase (they are 8 times lower), so you are
changing the case altogether!

If I were you I would first double check that your branch parameters R, X,
B are correct for the choice of MVAbase and the corresponding Voltage level.


and my last question:
> 4/how can I do the compensation of the reactive power  at nodes 1002 and
> 1001 (the 2 substations) and at pdl node 3001 and how can i  see the
> intermediate calculation of the reactive losses for each node because
> matpower gives us directly the results without seeing the details of the
> calculations.
>

I'm not sure I understand your question.  To find out the amount of
reactive compensation you need, here's one very simple thing you can do in
MATPOWER: introduce a new generator at the desired bus, with PG=0 (no real
power), and choose some sensible value for the desired voltage setpoint,
VG.  This will act as a "synchronous condenser", injecting or absorbing
just the right amount of reactive power Q so that the bus reaches the
voltage setpoint.  Solve the case and look at the solution. The Q injected
by this generator will give you exactly the amount of (local) reactive
support you need.  You can later convert this generator into a static
capacitor bank, if you prefer.  I mean, a bus capacitor whose BS value
would yield the same injection Q when the bus voltage is VG  (BS = Q /
VG^2).

-- 
Jose L. Marin
Grupo AIA

Re: CPF with Modal Analysis

[Jose Luis Marín]

Hi Mariusz,

I think the problem is that you have to provide a *solved* target case,
because the Jacobian should be evaluated *at the solution*.

Even for the base case, it wouldn't hurt to run a powerflow, since you
can't always trust that the input file is a solved case.

-- 
Jose L. Marin
Grupo AIA



2018-05-20 11:58 GMT+02:00 Mariusz Drabecki :

> Dear All,
>
>
>
> I am working on Modal Analysis. I would like to show to my students, that
> in the bifurcation point at least one of the reduced Jacobian eigenvalues
> is zero. This is why I calculate these values once for the base case and
> then for the nose point ( lambda max calculated by the CPF).
>
>
>
> Unfortunately there is no change in the results. Could anyone help?
>
>
>
> The results of the eigenvalues are as follows (for case14)
>
>
>
> Eigb=
>
> [
>
> 65,4207
>
> 39,9533
>
> 21,9939
>
> 18,9118
>
> 16,4492
>
> 11,3421
>
> 2,70658
>
> 5,56886
>
> 7,66373
>
> ]
>
>
>
> And *exactly the same ones* in the nose points.
>
>
>
> I would greatly appreciate any help from you as I got terribly stuck.
>
>
>
> Best,
>
> Mariusz Drabecki
>
>
>
>
>
> Here is my code:
>
>
>
> define_constants
>
> mpc=loadcase(case14);
>
> inc=2.5;
>
>
>
> %% CPF options
>
> mpopt = mpoption('out.all', 0, 'verbose', 2);
>
> mpopt = mpoption(mpopt, 'cpf.stop_at', 'NOSE', 'cpf.step', 0.1);
>
> mpopt = mpoption(mpopt, 'cpf.plot.level', 2);
>
> %% Base case (base solution operating points)
>
> %CPF
>
> mpcb=mpc;
>
> mpct = mpcb; % set up target case with
>
> mpct.gen(:, [PG QG]) = mpcb.gen(:, [PG QG]) * inc; % increased generation
>
> mpct.bus(:, [PD QD]) = mpcb.bus(:, [PD QD]) * inc; % and increased load
>
> cpf1 = runcpf(mpcb, mpct, mpopt);
>
> lambda_max=cpf1.cpf.max_lam;
>
>
>
> %Modal Analysis for Base Case
>
> [jb11,jb12,jb21,jb22,jacb,~,~,~] = makeJac(mpcb);
>
> jb11=full(jb11);
>
> jb12=full(jb12);
>
> jb21=full(jb21);
>
> jb22=full(jb22);
>
> Hb=jb11;
>
> Mb=jb12;
>
> Nb=jb21;
>
> Kb=jb22;
>
> % The reduced Jacobian for the base case (Q-V sensitivities)
>
> jacb_red=Kb-Nb*inv(Hb)*Mb;
>
>
>
> % The eigenvalues of the reduced Jacobian (Q-V sensitivities)
>
> [leftb eigb rightb]=(eig(jacb_red));
>
>
>
> %% Target case (nose point)
>
> mpct.gen(:, [PG QG]) = mpct.gen(:, [PG QG]) * lambda_max; % increased
> generation
>
> mpct.bus(:, [PD QD]) = mpct.bus(:, [PD QD])* lambda_max; % and increased
> load
>
> % The  Jacobian for the target case
>
> [jt11,jt12,jt21,jt22,jact,~,~,~] = makeJac(mpct);
>
> jt11=full(jb11);
>
> jt12=full(jb12);
>
> jt21=full(jb21);
>
> jt22=full(jb22);
>
> njact=size(jact,1);
>
> Ht=jt11;
>
> Mt=jt12;
>
> Nt=jt21;
>
> Kt=jt22;
>
> % The reduced Jacobian for the target case
>
> jact_red=Kt-Nt*inv(Ht)*Mt;
>
> % Eigenvalue decomposition of the reduced Jacobian for the target case
>
> [leftt eigt rightt]=(eig(jact_red));
>
>
>
>
>
>
>
>
>

Re: AC case network did not converge

[Jose Luis Marín]

Hello Asma,

I gave your case a quick look and I found several problems with the file
format (many extraneous lines in several places).  Anyway, I was courious
and I found that I could edit it by hand quickly because it's actually a
rather simple and repetitive topology structure.  I'm attaching my fixed
file, have a look at it (I also converted it to MATPOWER version 2).  The
second version contains an additional fix that I guessed (bus 1002
connecting to bus 3001, instead of 1001).  By the way, you can quickly
visualize the topology with this handy tool:
http://immersive.erc.monash.edu.au/stac/

So it's a *very* radial case, with no loads, only (wind) generators, and
the swing acting as the "sink" for all of the wind-farm power, as it is the
connection point to the grid.

>From the point of view of powerflow, your case has two kind of problems:

   1. The main one is that this network just cannot "evacuate" all that
   power. The case is infeasible, as you can see by running a CPF.
   2. A smaller problem is that all generators are running in PQ mode, with
   Q=0 (power factor 1.0), which means that you won't have enough "reactive
   support" to maintain voltage levels for such long radial sections.

Allowing the generators to inject at least some Mvars will let you transfer
a bit more power from the wind farm.  See the attached test script.  These
are the results I get (for caseborselle_fixed2):

   - Operating gens at cos_phi=1.0  ==>  max_lambda = 0.546875, mean
   voltage = 0.802666
   - Operating gens at cos_phi=0.95 ==>  max_lambda = 0.807189, mean
   voltage = 0.986915
   - Operating gens at cos_phi=0.90 ==> max_lambda = 0.968138, mean voltage
   = 1.104592

Hope it helps,

-- 
Jose L. Marin
Grupo AIA



2018-05-18 18:43 GMT+02:00 Asma DABBABI :

> Hello Jose,
> Thank you for rapid response, i'm quite in a rush to finish this study
> case for my thesis research.
>
> - i have a radial topology of wind turbines ( 100 wind turbines with 8 MW
> each) like i said i represent all my turbines nodes like a PQ node
> with P=8Mw et Q=0MVAR. and i have 3 AC offshores plateformes (ie 3
> transformers ) related to my wind farm (each transformer is related to one
> cluster in my wind farm) after that i have a long transmission line (
> distance about 60km) and finally i have a  terrestrial network which is
> composed of  one final transformer and the output node of the transformer
> is my slack bus ( my only generator).
>
> I put you my case file so you can more understand what i mean.
>
> from what I understood from your answer is that I have to add PV nodes and
> i have to add also Q for each turbine ?
>
> thank you again and i hope you can help me.
>
>
> Best regards
>
>
> > Hello Asma,
> >
> > Your description is not very clear, so let me ask you a couple of things
> > first:
> >
> >- You seem to be saying that your generators should have P=8 MW.
> > That's
> >fine, but then you total load should be balanced according to that (is
> >it?). Otherwise the imbalance, which has to be picked up by the slack
> > bus,
> >could render the case infeasible.
> >- You also seem to be saying that all your generators are PQ type.  A
> >case with no PV buses is quite prone to numerical ill-conditioning
> > (unless
> >the total number of buses is low).  Worse yet, if all gens have Q=0
> > then
> >all reactive power will have to be provided by the slack, and then
> > voltages
> >will quickly drop as you move away from the slack bus.  And if the
> > network
> >topology is dominantly radial, things are even worse yet.
> >
> > Maybe we can help you faster if you could share your case file here.
> >
> > --
> > Jose L. Marin
> > Grupo AIA
> >
> >
> >
> > 2018-05-18 17:31 GMT+02:00 Asma DABBABI  >:
> >
> >> the library Matpower works perfectly with small wind farms (30 wind
> >> turbines with
> >> 3 mw each, the turbines are represented as PQnode with P=3MW and Q=0MVAR
> >> and i have one slack bus in the terrestrial network)
> >> but when I tried to make an example of a larger farm such as Borsselle
> >> Park (100 wind turbines of 8 mw each) I had a divergence of the system.
> >>
> >> "Newton's method power flow did not converge in 10 iterations ",
> >> I looked on the internet and I found that in this case we have to do a
> >> continuation power flow that gradually increases the loading/generation
> >> like this:
> >>
> >> {define_constants;
> >> mpcbase = loadcase ('casefile');
> >> mpcbase.bus (:, PD) = 0;
> >> mpcbase.bus (:, QD) = 0;
> >> mpcbase.gen (:, PG) = 0;
> >> mpctarget = loadcase ('casefile');
> >> results = runcpf (mpcbase, mpctarget);
> >> results.cpf.max_lam}
> >>
> >> If the resulting value of results is greater than 1, it indicates that
> >> the
> >> load for the box is greater than that of the loading and unloading
> >> system
> >> at least by a factor of results.cpf.max_lam 

Re: AC case network did not converge

[Jose Luis Marín]

Hello Asma,

Your description is not very clear, so let me ask you a couple of things
first:

   - You seem to be saying that your generators should have P=8 MW.  That's
   fine, but then you total load should be balanced according to that (is
   it?). Otherwise the imbalance, which has to be picked up by the slack bus,
   could render the case infeasible.
   - You also seem to be saying that all your generators are PQ type.  A
   case with no PV buses is quite prone to numerical ill-conditioning (unless
   the total number of buses is low).  Worse yet, if all gens have Q=0 then
   all reactive power will have to be provided by the slack, and then voltages
   will quickly drop as you move away from the slack bus.  And if the network
   topology is dominantly radial, things are even worse yet.

Maybe we can help you faster if you could share your case file here.

-- 
Jose L. Marin
Grupo AIA



2018-05-18 17:31 GMT+02:00 Asma DABBABI :

> the library Matpower works perfectly with small wind farms (30 wind
> turbines with
> 3 mw each, the turbines are represented as PQnode with P=3MW and Q=0MVAR
> and i have one slack bus in the terrestrial network)
> but when I tried to make an example of a larger farm such as Borsselle
> Park (100 wind turbines of 8 mw each) I had a divergence of the system.
>
> "Newton's method power flow did not converge in 10 iterations ",
> I looked on the internet and I found that in this case we have to do a
> continuation power flow that gradually increases the loading/generation
> like this:
>
> {define_constants;
> mpcbase = loadcase ('casefile');
> mpcbase.bus (:, PD) = 0;
> mpcbase.bus (:, QD) = 0;
> mpcbase.gen (:, PG) = 0;
> mpctarget = loadcase ('casefile');
> results = runcpf (mpcbase, mpctarget);
> results.cpf.max_lam}
>
> If the resulting value of results is greater than 1, it indicates that the
> load for the box is greater than that of the loading and unloading system
> at least by a factor of results.cpf.max_lam to get a convergent power flow
> solution.
> I tried this and actually I found that it shows that results.cpf.max_lam
> is less than 1 so i must reduce the loads  to 1.35 mw for each turbine.
> in my case i don't need to reduce the loads, i must found the losses of
> the network with 8mw for each turbine.
> I've tried many solutions that i found in matpower page but until now I
> have not found the right solution to calculate with 8 mw.
>
> So can you help me please on that, I 'll be so grateful if you answer me.
>
> Best regards,
>
>
>
>

Re: Branch model - switch FROM and TO bus?

[Jose Luis Marín]

Hello Elis,

The ordering matters for the correct specification of the *tap ratio*.
>From the manual (p. 21--22):

"The transformer, whose tap ratio has magnitude $\tau$ and phase shift
angle $\theta_{shift}$, is *located at the from end of the branch*, as
shown in Figure 3-1".

If, for whatever reason, you wanted to swap FROM and TO labels, then the
values of the tap ratio would have to change accordingly (it's not hard to
figure out the formulas for this transformation).

Note also that the manual defines N as the *complex* value of the tap
ratio, so when you say N>1 or N<1, I assume you're thinking about its
modulus, tau.  This value can be greater or less than 1, either way.  And
MATPOWER doesn't care which side is the high-voltage or low-voltage one, as
long as your tap ratio is correctly modeled.

Hope it helps,

-- 
Jose L. Marin
Grupo AIA



2018-05-09 14:17 GMT+02:00 Elis Nycander :

> Hi,
>
> In the matpower branch model (Fig 3.1) the impedance is on the TO-bus
> side. However, since the impedance is in pu, for an ideal transformer we
> might as well move the pu impedance to the FROM-bus right? So in this case
> we can exchange the TO-bus and FROM-bus for any branch without effecting
> the results of a load flow.
>
> Q1. So when does the the ordering of the FROM and TO bus actually matter?
> Is there any case when we cannot swich the TO and FROM bus of a branch and
> still get the same result from a load flow or opf?
>
> Q2. From Figure 3.1 it looks as if the FROM bus is the high voltage bus
> (assuming N > 1). Is this always the case or can we have N < 1 just as well?
>
> Regards,
> Elis
>

Re: Power Flow in Matpower

[Jose Luis Marín]

Hello Mohammed,

The problem is in the *branch joining buses 16 and 17* (line 191 in the
file).  This branch had zero R and zero X.  Just edit the value of X to
something reasonable (I used 0.001), and the case solves just fine, in 4
iterations.

-- 
Jose L. Marin
Grupo AIA



2018-03-21 3:51 GMT+01:00 Mohammed Alhajri :

> i got
>
> 
>
> 19980   NaN
> 19981   NaN
> 19982   NaN
> 19983   NaN
> 19984   NaN
> 19985   NaN
> 19986   NaN
> 19987   NaN
> 19988   NaN
> 19989   NaN
> 19990   NaN
> 19991   NaN
> 19992   NaN
> 19993   NaN
> 19994   NaN
> 19995   NaN
> 19996   NaN
> 19997   NaN
> 19998   NaN
> 1   NaN
> 2   NaN
> Gauss-Seidel power flow did not converge in 2 iterations.
>
> >  Did NOT converge (47.99 seconds)  <
>
>
> On 20 March 2018 at 19:44, Ray Zimmerman  wrote:
>
>> Unfortunately, I do not have time to work on this myself. I was just
>> giving a suggestion for another direction to try if you want to understand
>> the issue that MATPOWER is having with your case. Could you post the output
>> (using verbose set to 2) of runpf() when using a MATPOIWER case file
>> that corresponds to the solved case from the Hadi Sadat code?
>>
>> And if you have any questions about the MATPOWER case format or MATPOWER
>> power flow options, feel free to ask.
>>
>> Ray
>>
>>
>>
>> On Mar 20, 2018, at 11:28 AM, Mohammed Alhajri  wrote:
>>
>> Ok, i have attached the case information in format of Hadi Saadat code,
>> can you please try it in MATPOWER?
>>
>> because we have spent more than three weeks checking the format, but
>> still dose not converge...
>>
>> Regards,,,
>>
>> بتاريخ ٢٠١٨/٠٣/٢٠ ٦:٢٦ م، كتب "Ray Zimmerman" :
>>
>>> It’s possible that the modeling is not identical or that there is some
>>> error in your conversion to MATPOWER format. You can check by talking the
>>> solved case from your other software, converting that solved case to
>>> MATPOWER and then trying the MATPOWER power flow. It should converge in a
>>> single iteration. If it does not, then you know that there is either a
>>> mistake somewhere or a difference in modeling.
>>>
>>>Ray
>>>
>>>
>>> On Mar 16, 2018, at 12:42 PM, Mohammed Alhajri 
>>> wrote:
>>>
>>> i tried that but unfortunately not work
>>>
>>> بتاريخ ٢٠١٨/٠٣/١٦ ٨:٣٠ م، كتب "Abhyankar, Shrirang G." >> >:
>>>
 See FAQ #5 



 Thanks,

 Shri

 Ph: (630) 252 0219 <(630)%20252-0219>

 www.mcs.anl.gov/~abhyshr







 *From: * on behalf of
 Mohammed Alhajri 
 *Reply-To: *MATPOWER discussion forum 
 *Date: Friday, March 16, 2018 at 11:26 AM*
 *To: MATPOWER discussion forum >*
 *Subject: Re: Power Flow in Matpower*



 any answer to this question?



 بتاريخ ٢٠١٨/٠٢/٢٥ ٧:١٨ م، كتب "Mohammed Alhajri" :

 Hello All,



 I did the power flow for a 89-bus network and it converges using Hadi
 Sadat code after 17080 iterations. The accuracy was 1e-8 and the method is
 Gauss-Seidel Method.



 But when I did the power flow using matpower it does not converge! I
 tried to increase the maximum iteration, I put it 10, and still did not
 cnvarge!



 I have attached the data according to Hadi Sadat Code, can any one try
 to do the power flow using matpower?




>>>
>>
>

Re: PF did not converge

[Jose Luis Marín]

Hi Mohammed,

If you can share the mpc file you have constructed with all this data,
maybe I could help.

-- 
Jose L. Marin
Grupo AIA



2017-12-14 18:55 GMT+01:00 Mohammed Alhajri :

> Hello
>
> i tried to do a PF for a system whose data is attached, but PF did not
> converge!
>
> Please help me, although i have tried my best to follow the manual, but
> the problem is still there.
>
> regards
>

Re: Dear

[Jose Luis Marín]

Actually, you have to change the setpoint VG of the generator (or
generators) attached to the slack bus:

define_constants;
mpc = loadcase('test'); % "test" is the name of your casefile
mpc.gen(k, VG) = 1.023456 ; % k is the generator index of the generator(s)
attached to your Slack bus

-- 
Jose L. Marin
Grupo AIA


2017-09-14 10:29 GMT+02:00 Jascha Korscheck :

> Thank you for the answer. This is the command I search for.
>
>
> best regards
>
> Jascha
>
>
>
> Am 13.09.2017 um 13:47 schrieb Ehsan Hejri:
>
> Dear Jascha Korscheck,
> I think you should change the slack bus voltage magnitude like this:
>
> define_constants;
> mpc = loadcase('test'); % "test" is the name of your casefile
> mpc.bus(1, VM) = 1.06 ;  % 1 is your Slack bus number %1.06 p.u. is just
> an example voltage magnitude for slack bus
>
> Best wishes,
> Ehsan
>
>
> On Wed, Sep 13, 2017 at 1:10 PM, Jascha Korscheck <
> jascha.korsch...@dutrain.de> wrote:
>
>> Dear MATPOWER-Team
>>
>> I've described a grid with the name "test". This grid has only three
>> lines with four buses.
>>
>> Bus 1 is connect over line 1 to bus 2, bus is connect over line 2 with
>> bus 3 and bus 3 is connect over line 3 with bus 4.
>>
>> Bus 1 is my slag.  Is the slag-voltage always 1 p.u. or can I change this
>> value?
>>
>> If I can change it, with  which command can I change it?
>>
>> I hope you can help me.
>>
>> --
>> with kindly regards
>>
>> Jascha Korscheck
>>
>> Tel.: +49 (0)2065 839555 – 7
>> Fax: +49 (0)2065 68999 – 8
>> jascha.korsch...@dutrain.de
>>
>> DUtrain GmbH
>> - Independent service & training centre
>>   for power system control
>> Dr.-Alfred-Herrhausen-Allee 16
>> 
>> D-47228 Duisburg
>> Germany
>>
>> www.dutrain.de
>>
>> Geschäftsführung: Dr.-Ing. Udo Spanel
>>   Dr.-Ing. Carsten Roggatz
>> Amtsgericht Duisburg - HRB 7009
>> VAT-ID: DE 1833 51 736
>>
>>
>>
>
>

Re: Regarding assigning initial values for load flow analysis in NR method in runpf() function

[Jose Luis Marín]

Hi Sayandev,

As you know, iterative methods need an initial guess.  But this is the
user's *choice*: you can use a flat start as you mention, or any other
information you may have (for instance, from previously solved cases that
are similar to the one at hand---as it happens when solving N-1 contingency
analysis cases).  MATPOWER cannot make this choice for you; so it will just
take the values of bus(:,[VM VA]) present in the input file as the initial
guess (except for VM in PV buses, where the value is replaced by the
regulation voltage VG of the corresponding generators).

So, if you want to test a "flat start", you have to explicitly set the
values of VM and VA before you call runpf().  By the way, setting the VM to
the swing's Vref is (arguably) slightly better than setting them to 1.  In
any case, note that "flat start" is by no means a guaranteed way to obtain
the correct solution (or a solution at all).

-- 
Jose L. Marin
Grupo AIA



2017-08-08 10:11 GMT+02:00 Sayandev Ghosh :

> Beg your pardon the missing snippet from the previous mail is as follows:
>
>
> %% initial state
>
> %V0= ones(size(bus, 1), 1);%% flat start
>
> V0  = bus(:, VM) .* exp(sqrt(-1) * pi/180 * bus(:, VA))
>
> vcb = ones(size(V0));   %% create mask of voltage-controlled
> buses
>
> vcb(pq) = 0;%% exclude PQ buses
>
> k = find(vcb(gbus));%% in-service gens at v-c buses
>
> V0(gbus(k)) = gen(on(k), VG) ./ abs(V0(gbus(k))).* V0(gbus(k));
>
>
> *SAYANDEV GHOSH*
> Assistant Professor
> Electrical Engineering Department
> Chaibasa Engineering College, Chaibasa, Jharkhand
> Mob. No. - 7292880077 / 9748509620
>
> On Tue, Aug 8, 2017 at 1:31 PM, Sayandev Ghosh 
> wrote:
>
>> Hello,
>> I am Sayandev Ghosh, presently teaching power system in an undergraduate
>> engineering college. I have done my M.Tech on Energy Science from Jadavpur
>> University. I have been using MATPOWER for load flow analysis of till 118
>> bus by Newton Raphson method. While doing that I am not able to understand
>> the following area:-
>>
>>
>>  In the runpf() function while using Newton Raphson method  the initial
>> values of voltage magnitude and voltage angle have been assigned by the
>> following lines (line 177 - 183):
>>
>>   [image: Inline image 2]
>> I am having a bit problem in understanding this section. Why are we
>> considering the final values of Voltage magnitude and Voltage angle as the
>> initial guess? In the case files for eg. 'case57' the Vm and Va are
>> probably the final voltage magnitudes and voltage angles.
>>
>> Shouldn't it be like assigning all the voltage magnitudes of PQ bus as 1.0
>> pu and voltage angles of all PV and PQ bus as 0 degree for a flat start.
>>
>> Depending on the answer of this question I would like to convey my other
>> doubts/questions.
>>
>> Thanks in advance.
>>
>>
>>
>> *SAYANDEV GHOSH*
>> Assistant Professor
>> Electrical Engineering Department
>> Chaibasa Engineering College, Chaibasa, Jharkhand
>>
>
>

Re: DC Power Flow

[Jose Luis Marín]

I don't mean to be disrespectful, but you should probably read this first:

   http://www.catb.org/esr/faqs/smart-questions.html

Your question is way too broad and vague.

-- 
Jose L. Marin
Grupo AIA



2017-07-17 15:21 GMT+02:00 Saeed Ahmed :

> I need to perform DC State Estimation and calculate the Real Power Flow in
> the lines. Need help please.
>
>
>
>
> Regards
> Saeed Ahmed
>
>

Re: Qestions on case file data in Matpower

[Jose Luis Marín]

Check if they have it in one of these repositories, which are both part of
ARPA-E's program "GRID DATA" (https://arpa-e.energy.gov/?q=
arpa-e-programs/grid-data):

   1. http://www.bettergrids.org/  Maintained by GridBright Inc and Utility
   Integration Solutions LLC (UISOL, a GE Company).
   2. https://egriddata.org/  Maintained by Pacific Northwest National Lab
   (PNNL) and the National Rural Electric Cooperative Association (NRECA).

Check also the NESTA archive (
https://research.csiro.au/data61/nicta-energy-system-test-case-archive-nesta/
).

-- 
Jose L. Marin
Grupo AIA



2017-07-11 6:58 GMT+02:00 Stephanie :

> Dear all,
>
>
> I am using Matpower with IEEE-57-bus system, and I’d like to test the case
> of 24h, does anyone knows where I can get load and generation data of 24h
> for IEEE-57-bus system, I greatly appreciate it if anyone can give some
> advice, thanks in advance.
>
> Best Regards
> Stephanie
>

Re: NREL-118 test case

[Jose Luis Marín]

Hi Zegal,

You can check if they have it in one of these repositories, which are both
part of ARPA-E's program "GRID DATA" (https://arpa-e.energy.gov/?q=
arpa-e-programs/grid-data):

   1. http://www.bettergrids.org/  Maintained by GridBright Inc and Utility
   Integration Solutions LLC (UISOL, a GE Company).
   2. https://egriddata.org/  Maintained by Pacific Northwest National Lab
   (PNNL) and the National Rural Electric Cooperative Association (NRECA).


-- 
Jose L. Marin
Grupo AIA



2017-06-29 20:50 GMT+02:00 Zegal, Walid :

> Hello everybody,
>
> I would like to know if there is a MatPower test case for NREL-118 test
> (An Extended IEEE 118-bus Test System with High Renewable Penetration).
>
> Thank you, best regards
>

Re: α and β of each transmission line

[Jose Luis Marín]

I think he may be referring to the real and imaginary parts of the line
propagation constant gamma (
https://en.wikipedia.org/wiki/Propagation_constant).  Note that this is a
distributed parameter (i.e., per unit length of the line), while the
parameters given in the *.m file are supposed to be lumped parameters.

Strictly speaking, you cannot recover the value of alpha and beta unless
you know how the lumped parameters were obtained from the distributed
ones.  If you assume that they were simply obtained by multiplying the
manufacturer's per-unit-length values by the length of the line (which is
only valid for "short" lines), *and* if you knew such length, then you
could obtain alpha and beta very easily.  It's  described in Kundur's book,
Chapter 6.1.4.

Schematically, in pseudo-code:

r = R / length
x = X / length
bshunt = Bshunt / length
alpha + j beta = sqrt ( (r + j x) * j bshunt)


-- 
Jose L. Marin
Grupo AIA


2017-05-30 22:22 GMT+02:00 Ray Zimmerman :

> The model used for the transmission line shown in Figure 3-1 in the User’s
> Manual. I’m not sure what you mean by the α and β of the transmission line.
> What do these parameters refer to?
>
>Ray
>
>
>
> > On May 28, 2017, at 10:04 PM, Ameeratul Balqis <
> ameeratulbal...@gmail.com> wrote:
> >
> > Hello and good day,
> >
> > I have been using command of 'runpf(casedata)' in MATPOWER in order to
> get the power flow of the system case. However, I realized that the result
> did not include the α and β of each transmission line.
> >
> > 1.How can I get the α and β of each transmission line? What command
> should I use?
> > 2.What is the .m file in order for me to able to learn more or edit the
> equation of the α and β in the lines.
> >
> > I'm new in using MATPOWER so I hope I would get help and guidance
> regarding the above problem I've stated.
> >
> > Thank you
>
>
>
>

Re: Power flow VM VA, CPF end point

[Jose Luis Marín]

1. Yes, VM VA are used both for input and output. Note one subtle point,
though: in runpf.m (Lines 177--183) the initial seed for the iteration is
first set to the values [VM VA] provided as input, but for
voltage-controlled buses (with active generators), the value of VM is
replaced by the setpoint VG of their respective bus generator(s) (the value
of VA is preserved). This makes sense because it guarantees that the seed
will be closer to the solution.

2. What you see here is to be expected. Basically, what happens is that the
basins of attraction of low-voltage volutions are usually smaller than
those of the operating solution. Also, note that there are *many* (for N
buses, somewhat of the order of 2^N) low voltage solutions, so chances are
that the iteration will converge to a solution that's different from the
one you followed by homotopy (CPF). Even homotopy methods can encounter
this same problem if their step-size is not small enough.


Here's some refs on the fractal nature of the problem:

@INPROCEEDINGS{KlumpOverbye00b,
author={Klump, R.P. and Overbye, T.J.},
booktitle={Power Engineering Society Summer Meeting},
title={A new method for finding low-voltage power flow solutions},
publisher={IEEE},
year={2000},
volume={1},
pages={593--597},
doi={10.1109/PESS.2000.867653},
ISSN={}
}

@INPROCEEDINGS{ThorpNaqavi89,
author={Thorp, J.S. and Naqavi, S.A.},
booktitle={Proceedings of the 28th IEEE Conference on Decision and Control},
title={Load flow fractals},
year={1989},
volume={2},
pages={1822--1827},
doi={10.1109/CDC.1989.70472}
}

@INPROCEEDINGS{ThorpNaqaviChiang90,
author={Thorp, J.S. and Naqavi, S.A. and Chiang, H.-D.},
booktitle={Decision and Control, 1990., Proceedings of the 29th IEEE
Conference on},
title={More load flow fractals},
year={1990},
month={dec},
volume={6},
pages={3028--3030},
doi={10.1109/CDC.1990.203339}
}

@ARTICLE{ThorpNaqavi97,
author={Thorp, J.S. and Naqavi, S.A.},
journal=IEEE_M_CAP,
title={Load-flow fractals draw clues to erratic behaviour},
year={1997},
month={jan},
volume={10},
number={1},
pages={59--62},
doi={10.1109/67.560872},
ISSN={0895-0156}
}

@INPROCEEDINGS{Mori00,
author={Mori, H.},
booktitle={IEEE International Symposium on Circuits and Systems (ISCAS)},
title={Chaotic behavior of the Newton-Raphson method with the optimal
multiplier for ill-conditioned power systems},
year={2000},
volume={4},
pages={237--240},
doi={10.1109/ISCAS.2000.858732}
}


-- 
Jose L. Marin
Grupo AIA




2017-05-19 11:30 GMT+02:00 Elis Nycander :

> Hi all matpower users!
>
> I have two questions:
>
> 1. In the bus matrix, the columns VM and VA are used both for the initial
> guess when solving the power flow, and to store the resulting voltages?
>
> 2. When solving a cpf, I get lam_max which corresponds to the nose point,
> i.e. maximum load/generation increase before "voltage collapse" happens. I
> can also find the power flow at the nose point just by running an ordinary
> power flow with flat start and conditions corresponding to the maximum
> load. However, I have tried to do the same thing for the lower part of the
> PV curve but failed to reproduce the solutions from the cpf. Basically I
> thought I could get the lower part of the PV curve by just solving a power
> flow using initial conditions which are close to the "unstable"/lower
> solutions from the cpf (instead of a flat start), but I get different
> solutions.
>
> Thanks,
> Elis
>

Re: Wrong solution and Losses on the branches

[Jose Luis Marín]

Hi Mahraz,

What is it exactly that makes you think the solution is wrong?  If you want
to perform an independent double-check on the solution, one way to do it
consists simply in verifying current mismatches at every bus (which is the
very expression of the power flow equations).

When calculating flows across branches, you've got to be aware of their
PI-model topology.  Line charging susceptances result in two shunt
capacitors---they typically inject reactive power, which has to be taken
into account when calculating flow balances accross the branch.

-- 
Jose L. Marin
Grupo AIA



2017-04-10 22:03 GMT+02:00 mahraz amini :

> Hi every one,
>
> I am trying to solve AC Power flow in matpower. I made my own model in
> matpower format and when I run it, it got converged in a few iterations and
> gave me the solutions. the problem is the given solution is not correct. I
> mean when I calculate the power flow based on given Vmag and Vang the
> mismatch is not only zero but also large numbers. I am wondering how it
> could happen.
>
>  Furthermore when I look at the flow at the branch, losses on the lines
> are not equal the difference between S_from and S_to, for example :
>
> Brnch  From  To   From Bus Injection   To Bus
> InjectionLoss (I^2 * Z)
>   #Bus   Bus   P (MW   Q (MVAr)  P (MW)
> Q (MVAr)  P (MW)   Q (MVAr)
>
>   1613   16   -168.57   26.10  168.93
>-47.47  0.357  3.82
>   1713   19122.42  -38.42 -122.07
>-4.820.349  3.64
>   1813   21107.36  -34.03 -107.17
> 3.360.189  1.97
>   1916   17-95.59   -2.66
> 95.74 -28.46  0.147  1.56
>
> I would really appreciate if you can help me understand why it is so.
>
> Thanks,
> Mahraz
>

Re: Load flow not convergent

[Jose Luis Marín]

By the looks of it, you may have bad data (or a badly formed) input case.
You seem to have one or more NaNs somewhere in bus columns PD, QD, or
BUS_TYPE.

Note that, unless you have wacky errors of this kind, the technique
mentioned in FAQ 5 (v) should always give you at least some solutions, at
some scaling.

-- 
Jose L. Marin
Grupo AIA



2017-03-31 16:49 GMT+02:00 Gamze Dogan :

> Hi,
>
> I am running a load flow on a very large grid but I cannot make it
> converge using Matpower.
> I have looked at all the discussions on the subject and test various
> things but neither the pf nor the opf function converge.
>
> When I look at the results (of the non-convergent solution) I see that
> Matpower could not compute the following values which are assigned to NaN:
>
> - Voltage magnitude and angle of buses
>
> - Reactive power output of generators
>
> - Reactive and active power injections at node from and to of each line
>
> When I tried to check the scaling factor (using the solution proposed on
> the website:5. Why does MATPOWER power flow not converge?
>  v.) I got the
> following error:
>
> ??? NaN's cannot be converted to logicals.
>
> Error in ==> printpf at 177
> nzld = find((bus(:, PD) | bus(:, QD)) & bus(:, BUS_TYPE) ~= NONE);
>
> Error in ==> runcpf at 392
> printpf(results, 1, mpopt);
>
>
> Given this, if anybody has an idea of what I should try, check, change to
> make the computation converge?
>
> Thank you,
>
> Gamze
>
>

Re: thanks & speeding up the Newton-Raphson power flow

[Jose Luis Marín]

Excellent news!  I don't have much time at the moment but, I definitely
want to try out the new code and experiment with other solvers, in
particular KLU from Tim Davis' SuiteSparse.

Just a few quick comments to summarize and clarify my previous post:

   - It's not obvious at all that MATLAB uses different sparse solvers
   dependng on the way you make the lu() call.   I re-read my post and
   realized that I probably didn't explain this in detail.  The problem is
   that MathWorks does not document this at all. Luckily, Tim Davis (the main
   guy behind SuiteSparse, and several of the sparse solvers in MATLAB)
   published how this works in his FACTORIZE paper
   <http://dl.acm.org/citation.cfm?id=2491498>.
   - The summary is: backslash and the four- and five-output syntax calls
   will invoke UMFPACK, while the two- and three-output syntax calls will
   invoke the GP solver.  Still, since all this is not officially documented,
   MathWorks may change it in future versions.
   - About the performance of sparse solvers, in general: to my knowledge,
   there are two very important factors:
  - The use of BLAS kernels. Methods with names such as "Supernodal"
  and "Multi-frontal" try to take advantage of dense BLAS calls as much as
  possible. In application areas such as finite-element methods, this is
  crucial for performance, because the matrices do exhibit many dense
  sub-blocks.  Not so in power systems matrices, which is the reason why
  "simplicial" methods are faster (using the BLAS would incurr in too much
  overhead -- too many BLAS calls for sub-blocks of very small size).
  - Fill-reducing reorderings. This affects all sparse methods,
  regardless of the use of the BLAS.  The problem is that the reordering
  algorithms have different performance for different types of matrices.
  Experimentally, people have found that AMD is the best for electrical
  circuits and power system matrices.


As you said, it is very likely that the choice of fill-reducing algorithm
has a greater influence than the choice of UMFPACK vs. GP.  I can imagine
that UMFPACK disables calls to the BLAS for very small blocks of, say,
3x3.  If that were the case, both methods would be similar for the kind of
matrices we typically have here.

On the other hand, Davis makes this comment in his paper:

*GP is faster than UMFPACK for small matrices and for some circuit
matrices, but in general the results in Figure 3 show that UMFPACK
([L,U,P,Q,R]=lu(A)) is the best choice for most sparse unsymmetric
matrices. UMFPACK does poorly for some circuit matrices because its
automatic ordering method selection makes the wrong choice (it selects
COLAMD, but AMD works much better for those matrices). *
So the problem seems to be that you can't tell MATLAB's UMFPACK what
reordering algorithm to use.

It's great to have the code you prepared, it'll be fun for experimenting.

-- 
Jose L. Marin
Grupo AIA




2017-02-22 21:12 GMT+01:00 Ray Zimmerman <r...@cornell.edu>:

> First, a big thanks to Richard Lincoln and Mirko Todorovski for their
> valuable contributions since we put MATPOWER on GitHub. Richard contributed
> the Travis CI integration <https://travis-ci.org/MATPOWER/matpower/> that
> now automatically runs the full test suite on every commit and pull request
> (supplying a little green check mark when everything passes). And Mirko
> contributed some distribution power flow algorithms and cases.
>
> Second, in response to the very informative post to MATPOWER-L back in Dec
> by Jose Luis Marin (see below), I made some modifications to the way the
> update step is computed in the Newton-Raphson power flow, resulting in
> significant speed improvements (up to 2x) for large models under Matlab
> (marginally better under Octave).
>
> You can see the details and comment in this pull request ...
>
> https://github.com/MATPOWER/matpower/pull/7
>
> … or by checking out the lu branch
> <https://github.com/MATPOWER/matpower/tree/lu> of the Git repository.
> Jose, I’m particularly interested in any further comments you might have.
> In particular, it appears that the vast majority of the speedup comes from
> the AMD reordering and not from the Gilbert-Peierls algorithm.
>
> Ray
>
>
>
>
>
> Begin forwarded message:
>
> *From: *Jose Luis Marín <mari...@aia.es>
> *Subject: **Re: RunPF vs GPU calculations*
> *Date: *December 1, 2016 at 12:52:05 PM EST
> *To: *MATPOWER discussion forum <matpowe...@list.cornell.edu>
> *Reply-To: *MATPOWER discussion forum <matpowe...@list.cornell.edu>
>
> Hello all,
>
> I've been meaning to report the following findings to Ray and the MATPOWER
> team for some time, but never got around doing it.  Well, here it is in
> summarized form:
>
>- There's a 

Re: RunPF vs GPU calculations

[Jose Luis Marín]

Ooops, sorry, they're *seconds* (for one runpf() run).

Here's the benchmarking script, it's pretty simple:

case_list = {'case300', 'case1354pegase', 'case1888rte', 'case1951rte', ...
'case2383wp', 'case2736sp', 'case2737sop', 'case2746wop', ...
'case2746wp', 'case2848rte', 'case2868rte', 'case2869pegase',
...
'case3012wp', 'case3120sp', 'case3375wp', ...
'case6468rte', 'case6470rte', 'case6495rte', 'case6515rte', ...
'case9241pegase', 'case13659pegase'};

opts = mpoption('pf.alg', 'NR', 'pf.tol', 1.0e-10, 'pf.nr.max_it', 100,
'out.all', 0, 'verbose', 0);
NREPEAT = 10;
timings = zeros(NREPEAT,1);


for k = 1:length(case_list)

cfile = char(case_list(k));
fprintf('*** BENCHMARKING %s:\t', cfile);
mpc = loadcase(cfile);
for i = 1:NREPEAT;
tic;
mpsol = runpf(mpc, opts);
timings(i) = toc;
end
fprintf('%11.3e\n', min(timings));
end


-- 
Jose L. Marin
Grupo AIA



2016-12-01 19:01 GMT+01:00 davor sutic <davor.su...@gmail.com>:

> Hi Jose,
> thanks for sharing. I have one quick question: What are the units of the
> your benchmark (seconds, ticks, etc)?
>
> Thanks a lot
>
> On Thu, Dec 1, 2016 at 6:52 PM, Jose Luis Marín <mari...@aia.es> wrote:
>
>> Hello all,
>>
>> I've been meaning to report the following findings to Ray and the
>> MATPOWER team for some time, but never got around doing it.  Well, here it
>> is in summarized form:
>>
>>- There's a very cheap way to speedup the J \ F step (actually, the
>>sparse LU decomposition) by factors of about x1.3 to x1.8, and maybe 
>> higher
>>for bigger cases.  It just involves invoking the MATLAB solver that is 
>> best
>>suited for power system matrices (more below).
>>- GPU computing, at least in its current state, is not a panacea for
>>the powerflow problem, because the kind of sparse linear algebra that's
>>needed in powerflow problems just doesn't offer enough opportunities to
>>exploit GPU architectures.
>>
>>
>> Try this simple change in *newtonpf.m*:
>>
>> marinjl@xwing:~/Library/matpower6.0b2$ diff newtonpf.m.ORIG newtonpf.m
>> 89c89,93
>> < dx = -(J \ F);
>> ---
>> > %dx = -(J \ F);
>> > % JLM: For performance, factorize using MATLAB's Gilbert-Peierls +
>> AMD (see T. Davis)
>> > Q = sparse(amd(J), 1:size(J,1), 1);  % AMD reordering
>> > [L, U, P] = lu(Q'*J*Q, 1.0);  % this form ensures that G-P is
>> called. Last arg: 0==no pivoting, 1.0==partial pivoting.
>> > dx = - Q * ( U \ ( L \ (P * Q' * F)));
>>
>>
>> My benchmarks (3GHz Intel Core 2 Duo) yield the following speedups:
>>
>> %UMFPACKG-Pspeedup
>> % ==
>> % case300:  1.689e-021.230e-02  1.37
>> % case1354pegase:   5.848e-023.652e-02  1.60
>> % case1888rte:  2.792e-022.130e-02  1.31
>> % case1951rte:  3.037e-022.169e-02  1.40
>> % case2383wp:   1.429e-018.711e-02  1.64
>> % case2736sp:   1.327e-017.512e-02  1.76
>> % case2737sop:  1.632e-019.074e-02  1.79
>> % case2746wop:  1.589e-018.725e-02  1.82
>> % case2746wp:   1.448e-018.636e-02  1.67
>> % case2848rte:  4.652e-023.301e-02  1.40
>> % case2868rte:  1.681e-019.661e-02  1.74
>> % case2869pegase:   2.402e-011.335e-01  1.79
>> % case3012wp:   1.125e-016.733e-02  1.67
>> % case3120sp:   2.155e-011.228e-01  1.75
>> % case3375wp:   1.598e-019.565e-02  1.67
>> % case6468rte:  2.155e-011.443e-01  1.49
>> % case6470rte:  2.316e-011.434e-01  1.61
>> % case6495rte:  2.118e-011.433e-01  1.47
>> % case6515rte:  2.244e-011.443e-01  1.55
>> % case9241pegase:   1.185e+006.392e-01  1.85
>> % case13659pegase:  1.484e+007.991e-01  1.85
>> %
>> %
>> % max speedup: 1.85
>> % min speedup: 1.31
>>
>> (note: I used pf.tol=1.0e-10 in order to trigger a few more iterations. I
>> ran each case 10 times, and kept the best timing)
>>
>>
>> *Longer explanation:*
>> The Jacobian matrices that one typically obtains in power networks are *very,
>> very* sparse.  They are very different than sparse matrices obtained in
>> other areas (such as, e.g., finite-element modeling), where the matrix
>> still has fairly dense blocks here and there (or at least it can be
>> re-arranged to have them).  Power system matri

Re: RunPF vs GPU calculations

[Jose Luis Marín]

Hello all,

I've been meaning to report the following findings to Ray and the MATPOWER
team for some time, but never got around doing it.  Well, here it is in
summarized form:

   - There's a very cheap way to speedup the J \ F step (actually, the
   sparse LU decomposition) by factors of about x1.3 to x1.8, and maybe higher
   for bigger cases.  It just involves invoking the MATLAB solver that is best
   suited for power system matrices (more below).
   - GPU computing, at least in its current state, is not a panacea for the
   powerflow problem, because the kind of sparse linear algebra that's needed
   in powerflow problems just doesn't offer enough opportunities to exploit
   GPU architectures.


Try this simple change in *newtonpf.m*:

marinjl@xwing:~/Library/matpower6.0b2$ diff newtonpf.m.ORIG newtonpf.m
89c89,93
< dx = -(J \ F);
---
> %dx = -(J \ F);
> % JLM: For performance, factorize using MATLAB's Gilbert-Peierls +
AMD (see T. Davis)
> Q = sparse(amd(J), 1:size(J,1), 1);  % AMD reordering
> [L, U, P] = lu(Q'*J*Q, 1.0);  % this form ensures that G-P is called.
Last arg: 0==no pivoting, 1.0==partial pivoting.
> dx = - Q * ( U \ ( L \ (P * Q' * F)));


My benchmarks (3GHz Intel Core 2 Duo) yield the following speedups:

%UMFPACKG-Pspeedup
% ==
% case300:  1.689e-021.230e-02  1.37
% case1354pegase:   5.848e-023.652e-02  1.60
% case1888rte:  2.792e-022.130e-02  1.31
% case1951rte:  3.037e-022.169e-02  1.40
% case2383wp:   1.429e-018.711e-02  1.64
% case2736sp:   1.327e-017.512e-02  1.76
% case2737sop:  1.632e-019.074e-02  1.79
% case2746wop:  1.589e-018.725e-02  1.82
% case2746wp:   1.448e-018.636e-02  1.67
% case2848rte:  4.652e-023.301e-02  1.40
% case2868rte:  1.681e-019.661e-02  1.74
% case2869pegase:   2.402e-011.335e-01  1.79
% case3012wp:   1.125e-016.733e-02  1.67
% case3120sp:   2.155e-011.228e-01  1.75
% case3375wp:   1.598e-019.565e-02  1.67
% case6468rte:  2.155e-011.443e-01  1.49
% case6470rte:  2.316e-011.434e-01  1.61
% case6495rte:  2.118e-011.433e-01  1.47
% case6515rte:  2.244e-011.443e-01  1.55
% case9241pegase:   1.185e+006.392e-01  1.85
% case13659pegase:  1.484e+007.991e-01  1.85
%
%
% max speedup: 1.85
% min speedup: 1.31

(note: I used pf.tol=1.0e-10 in order to trigger a few more iterations. I
ran each case 10 times, and kept the best timing)


*Longer explanation:*
The Jacobian matrices that one typically obtains in power networks are *very,
very* sparse.  They are very different than sparse matrices obtained in
other areas (such as, e.g., finite-element modeling), where the matrix
still has fairly dense blocks here and there (or at least it can be
re-arranged to have them).  Power system matrices, as well as many
electronic circuit matrices, have dense blocks of sizes 3x3 or 4x4 at most,
which is rather small.  Now, the vast majority of research going into
sparse LU solvers focuses on the former kind of matrices, i.e. those with
large dense blocks.  The methods go by names like "supernodal" and
"multifrontal", and basically consist in designing the algorithm so that it
can exploit the BLAS for processing those dense blocks.  The point is that
power system matrices are ill-suited to exploit the power of the BLAS, and
therefore those methods are actually slower than the more straight-forward
sparse solvers that don't try to use the BLAS at all.

MATLAB contains sparse solvers of both types, but the problem is that it's
not clearly documented how the backslash operator (or the lu call) uses
internally one algorithm or the other.  It turns out that the way you call
lu() influences the algorithm being used!!!  The default is UMFPACK, but if
you use the "three-output, two-input args" call as above then you'll be
invoking a Gilbert-Peierls method.  As you may have guessed by now, UMFPACK
is one of those methods using the BLAS, while G-P is not.  As shown in the
benchmarks above, UMFPACK clearly has an overhead.

Coming back to the issue of GPUs, one can see why the potential performance
gains for this type of sparse linear algebra are not as great as they seem.

If you are interested in more details, I can provide references.  See for
instance SuiteSparse by Tim Davies, which in terms of sparse solvers is
probably the best in class now  (some of it is already integrated within
MATLAB, btw).

-- 
Jose L. Marin
Grupo AIA




2016-11-30 11:02 GMT+01:00 davor sutic :

> Hello Jaroslaw,
> I find your analysis very interesting. Would you mind sharing your
> performance findings about a broader range of network sizes (e.g. 118, 300,
> perhaps even those of a few thousand)?
>
> Also, are the times you stated the mean 

Re: Help on radial LV distribution network

[Jose Luis Marín]

If you want to do things properly, I guess you would have to do *all* of
the following:

   - First, MATPOWER expects baseMVA to be in MW.  If you want to use a
   p.u. system based around Volts and 1 kW (instead of the traditional kV and
   100 MW), then you need baseMVA=0.001 at the beginning of your case file.
   - Expressing voltages in pu should be no problem, just divide the
   quantity in volts by your voltage levels (in your case I see it's only one,
   415 Volts)
   - Now, MATPOWER expects all power quantities (P, Q) to be expressed in
   MW.  So if your values are in kW, divide them by 1000.
   - Finally, to convert resistances and reactances (R, X), use your
   baseMVA and voltage base as follows: take the initial quantity in Ohms, and
   multiply it by:  baseMVA / Vbase^2 = 1000 Watts / (415 Volts)^2 =
   5.806357961968356e-03
   - You don't have Bshunt values in your case, but if you had, the
   conversion factor would be just the inverse of the one used for resistance
   and  reactance.

I hope I'm not missing anything, I think that's all you need in your case.

-- 
Jose L. Marin
Grupo AIA




2016-10-26 14:44 GMT+02:00 Nazurah Nasir :

> But does that means I should not divide my input power data by 1000 to
> make it in MW? If I do that, it won't converge. For example, these are my
> Power input for one time:
>
> Columns 1 through 6
> 1.3. 0 0 0 0
> 2.1.0.00120.0004 0 0
> 3.1.0.00190.0006 0 0
> 4.1.0.00060.0002 0 0
> 5.1.0.00240.0008 0 0
> 6.1.0.00120.0004 0 0
> 7.1.0.00100.0003 0 0
> 8.1.0.00230.0008 0 0
> 9.1.0.00050.0002 0 0
>10.1.0.00060.0002 0 0
>11.1.0.00120.0004 0 0
>12.1. 0 0 0 0
>   Columns 7 through 12
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
> 1.1. 00.41501.1.1000
>   Column 13
> 0.9400
> 0.9400
> 0.9400
> 0.9400
> 0.9400
> 0.9400
> 0.9400
> 0.9400
> 0.9400
> 0.9400
> 0.9400
> 0.9400
>
> Thank you very much for the help
>
> Yours sincerely,
> Nur
>
> On Wed, Oct 26, 2016 at 11:14 PM, Nazurah Nasir 
> wrote:
>
>> Aren't I supposed to make the R and X in p.u. if I want to use them in
>> MATPOWER? Regardless, your simulation seems to be more sensible. But, I
>> just curious, so we don't necessarily change the R and X into p.u. values?
>>
>> Thanks for the response.
>>
>> On Wed, Oct 26, 2016 at 3:26 PM, Saranya A  wrote:
>>
>>> Hi Nur,
>>> Dont divide R and X with the voltage. I get the following power flow
>>> without those two lines.
>>>
>>> 
>>> runpf('LV10')
>>>
>>> MATPOWER Version 6.0b1, 01-Jun-2016 -- AC Power Flow (Newton)
>>>
>>> Newton's method power flow converged in 5 iterations.
>>>
>>> Converged in 0.02 seconds
>>> 
>>> 
>>> | System Summary
>>>   |
>>> 
>>> 
>>>
>>> How many?How much?  P (MW)Q
>>> (MVAr)
>>> ----  -
>>>  -
>>> Buses 12 Total Gen Capacity 261.0-302.0 to
>>> 302.0
>>> Generators11 On-line Capacity   250.0-300.0 to
>>> 300.0
>>> Committed Gens 1 Generation (actual)  1.2   0.5
>>> Loads 10 Load 1.0   0.3
>>>   Fixed   10   Fixed  1.0   0.3
>>>   Dispatchable 0   Dispatchable  -0.0 of -0.0  -0.0
>>> Shunts 0 Shunt (inj) -0.0   0.0
>>> Branches  11 Losses (I^2 * Z) 0.23  0.16
>>> Transformers   0 Branch Charging (inj) -0.0

Re: Problem with runcpf

[Jose Luis Marín]

Just a guess: OPF observes limits QMIN and QMAX on generators, while a
powerflow (with the standard option pf.enforce_q_lims=0) does not.
Although not typical, it is certainly possible for a powerflow case to be
infeasible when solving without generator MVAR limits, and solvable when
enforcing them.

Take a look at your OPF solution and check if there are many generators
pegged at their Q limits or not.

-- 
Jose L. Marin
Grupo AIA



2016-09-08 14:45 GMT+02:00 Gamze Dogan :

> Thank you for your answer,
>
>
> I double check and my OPF does really converge. I also checked the buses
> and all my generators are PV and all my non-generators buses are PQ .
> I have zero island. If anyone else has an idea why I do not converge even
> if my opf converge?
>
> Thank you,
>
>
> Regards,
>
>
> Dogan
>
>
> --
> *De :* bounce-120759695-68616...@list.cornell.edu <
> bounce-120759695-68616...@list.cornell.edu> de la part de Abhyankar,
> Shrirang G. 
> *Envoyé :* mercredi 7 septembre 2016 14:16
>
> *À :* MATPOWER discussion forum
> *Objet :* Re: Problem with runcpf
>
> Are you sure the OPF converged? If it did, then I guess you are not
> setting the bus types correctly or there is still some error in your system
> data. Note that the first step of a continuation power flow is running a
> load flow. So if the load flow does not converge then continuation power
> flow will not either.
>
> On Sep 5, 2016, at 9:52 AM, Gamze Dogan  wrote:
>
> Thank you for your answer,
>
>
> I actually had an island, but the problem remains. Now I am sure that
> there are no island.
>
> The load flow does not converge.
> I tried to launch a fast decoupled, but it does not converge.
> The OPF converged and I used the results of Vm and Va as the initial guess
> to relauch my load flow, but it still does not converge.
>
>
>
>
> I tried to check the system load using runcpf, but it does not converge,
> here is the message:
>
> MATPOWER Version 5.1, 20-Mar-2015 -- AC Continuation Power Flow
> step   1 : lambda = 26502.011, corrector did not converge in 10 iterations
>
> >  Did NOT converge (0.20 seconds)  <
>
> Do you have any suggestion on how I can make it converge?
>
> Thank you,
>
> Dogan
>
>
> --
> *De :* bounce-120746902-68616...@list.cornell.edu <
> bounce-120746902-68616...@list.cornell.edu> de la part de Abhyankar,
> Shrirang G. 
> *Envoyé :* vendredi 2 septembre 2016 13:42
> *À :* MATPOWER discussion forum
> *Objet :* Re: Problem with runcpf
>
> I suspect you have isolated buses or islands that do not have a swing bus.
> Have you tried I-iv in MatPower faq on power flow divergence.
>
> http://www.pserc.cornell.edu/matpower/#pfconvergence
>
> Shri
>
> On Sep 2, 2016, at 6:51 AM, Gamze Dogan  wrote:
>
> Dear all,
>
> I have a problem when using runpf which does not converge. Therefore, I
> tried to check if the system load exceeded the steady-state loading limit
> (as mentionned on the website)
>
> *define_constants;
> mpcbase = loadcase('casefile');
> mpcbase.bus(:, PD) = 0;
> mpcbase.bus(:, QD) = 0;
> mpcbase.gen(:, PG) = 0;
> mpctarget = loadcase('casefile');
> results = runcpf(mpcbase, mpctarget);
> results.cpf.max_lam
> *But I have a problem using the function runcpf, it gives me the following 
> error message:
>
> Warning: Matrix is singular to working precision.
> > In cpf_corrector at 124
>   In runcpf at 257
> step   1 : lambda =NaN, corrector did not converge in 10 iterations
> ??? NaN's cannot be converted to logicals.
>
> Error in ==> printpf at 175
> nzld = find((bus(:, PD) | bus(:, QD)) & bus(:, BUS_TYPE) ~= NONE);
>
> Error in ==> runcpf at 392
> printpf(results, 1, mpopt);
>
> I do not understand what the problem is. Any suggestion is wellcome.
>
> Thank you,
>
> G. Dogan
>
>
>

Re: restrictions regarding the ratio of PQ/PV buses in test cases

[Jose Luis Marín]

Hello Davor,

This is actually a very interestng question from the point of view of the
numerical stability of the equations.  I asume that by "changing the bus
type PQ --> PV" you mean the following: start by solving a given case, then
take a particular PV bus and convert its generator injections PG,QG (values
taken from the solved case!) into a fixed PQ load, then remove the
generator, and finally switch the bus type to PQ.  Then you repeat this
progressively until you're left with no PV buses.

The resulting system obviously shares the same mathematical solution, but,
as you're guessing, trying to solve the converted system definitely shows
some numerical effects.  The way I understand it, there are differences
coming from these two sources:

   - The condition number of the Jacobian degrades: the Jacobian matrix of
   the NR method (and also FD methods) is related to the admittance matrix of
   the network, which is in turn related to the Laplacian matrix of the graph
   representing the transmission grid.  Swing buses remove the zero
   eigenvalue/eigenvector (the "uniform translation" mode in voltage), and PV
   buses sort of push the lowest non-zero eigenvalue (related to the so-called
   Fiedler vector) to higher values, thus improving the condition number of
   the resulting matrix.  It can be shown that a system with no PV buses would
   result in a condition number degrading linearly with network size.
   - The basins of attraction of the iterative methods will be different.
   This effect is in general much harder to analyze, so one has to resort to
   numerical experiments.


I experimented a while back with this and I found the second effect to be
much more dominant in practice.  My set up was the following: try to solve
cases always from a flat start, while progressively switching PV buses to
PQ.  What I found was that the loss of precision resulting from having
slighly worse condition numbers in the Jacobian was negligible, because one
would encounter NR non-convergence problems way before that could become
the dominant problem.

I'm attaching a test script I used for experimenting with all this.

By the way, there is also another important effect to take into account: PV
buses that are pushed over to "the other side" of their V-Q curve.  In
these buses, when you flip them to PQ, you may obtain a different voltage
solution (it depends on your initial seed). This is because in those cases
the solution as PV corresponds to a low-voltage branch when viewed as PQ.
Bus 191 in case IEEE 300 is a perfect example of this; you can try it with
my script.

-- 
Jose L. Marin
Grupo AIA


2016-07-27 9:13 GMT+02:00 davor sutic :

> Should I expect a convergent and (reasonably) accurate power flow
> solution, when I experiment with changing the bus type PQ->PV and
> vice versa?
>
> If so, is there a minimum number of each bus type required, e.g. consider
> a system where there are only PQ buses (apart from one slack, of course)?
>
> The data in the test cases seems sufficient for such changes, however I'm
> worried about the stability of the system of equations.
>
> Thanks a lot
>


convert_all_PV_to_PQ.m
Description: application/vnd.wolfram.mathematica.package

Re: How to use a Switch or Breaker in OPF?

[Jose Luis Marín]

MATPOWER, as most off-line planning tools, is "bus-branch" oriented, which
means that it simply does not have the concept of switches and breakers.
Each substation is typically represented with a single busbar at each
nominal voltage level.

However, you can simulate the action of breakers by setting the *status* of
generators (GEN_STATUS), lines (BR_STATUS), and buses (set BUS_TYPE=4).

-- 
Jose L. Marin
Grupo AIA



2016-07-19 21:26 GMT+02:00 RUBEN PEREZ GARCIA <100275...@alumnos.uc3m.es>:

> Hi everyone;
>
> I’m trying to solve a power network OPF and i have a question. In my grid
> i have 4 breakers. In one case 2 breakers are open and the others are
> closed, but in other case 3 are opened and one closed. I know that i can't
> switch the state during the solve, but i want to know how to indicate to
> the program that the breaker is opened or it's closed.
>
> Thank you very much.
>
>

Re: Replacing constant power with constant impedance loads

[Jose Luis Marín]

I think that you have to invert the sign of QD.  From
http://www.pserc.cornell.edu//matpower/docs/ref/matpower6.0b1/idx_bus.html:

columns 1-13 must be included in input matrix (in case file)
1  BUS_I   bus number (positive integer)
2  BUS_TYPEbus type (1 = PQ, 2 = PV, 3 = ref, 4 = isolated)
3  PD  Pd, real power *demand* (MW)
4  QD  Qd, reactive power *demand* (MVAr)
5  GS  Gs, shunt conductance (MW *demanded* at V = 1.0 p.u.)
6  BS  Bs, shunt susceptance (MVAr *injected* at V = 1.0 p.u.)
7  BUS_AREAarea number, (positive integer)
8  VM  Vm, voltage magnitude (p.u.)
9  VA  Va, voltage angle (degrees)
10 BASE_KV baseKV, base voltage (kV)
11 ZONEzone, loss zone (positive integer)
12 VMAXmaxVm, maximum voltage magnitude (p.u.)
13 VMINminVm, minimum voltage magnitude (p.u.)


-- 
Jose L. Marin
Grupo AIA


2016-07-06 12:07 GMT+02:00 Georgiadis Dionysios <
dionysios.georgia...@frs.ethz.ch>:

> Dear all,
>
>
> After some quick calculations one can see that it is possible to replace a
> constant power load with a constant impedance load of exactly equal power
> consumption.
>
> Here is an example, say bus 1 absorbs *S1* MVAs under *V1* volts. Bus 1
> is constant power.
> This would be equal to a constant impedance load, with impedance 
> *conjugate(Z1)
> = (Vn/V1)^2 * S1,* where Vn is the nominal voltage.
>
> Thus, we can set the PD and QD entries of the case to zero, and replace
> them with the right values for GS and BS.
>
> Specifically, after solving a case with all constant power loads we will
> obtain the Vm vector (voltage amplitudes) and the Pd and Qd vectors. Then,
> we can set all PD and QD values to zero and place* (1/Vm).^2 * Pd *and 
> *(1/Vm).^2
> * Qd* in the GS and BS columns respectively.
>
> From my understanding the output should be exactly the same. The solved
> case should remain a perfectly good solution to the PF equations. Yet, when
> I do so, the solver iterates.
>
> Can someone shed light on this? I apologise if I am asking the question in
> the wrong place.
>
>
> Best regards,
> Dionysios Georgiadis
>

Re: Performing N-1 contingency analysis

[Jose Luis Marín]

I'm attaching a structured script I wrote some time ago for carrying out
N-1 contingency analysis. You may use it as a guide, but I'm sure you will
need to make changes in order to adapt it to your specific needs (in
particular, the metrics I compute were very specific of the analysis I was
doing).  I only included branches, but it's easy to include other elements.

As Ray said, you've got to be careful with islanding (also beware of
possible changes in the swing bus, because of that).

-- 
Jose L. Marin
Grupo AIA



2016-07-05 18:36 GMT+02:00 Haroon Malik :

> Thank you very much for your response.
>
>
>
>
>
>
>  Original message 
> From: Ray Zimmerman 
> Date: 05/07/2016 5:10 p.m. (GMT+00:00)
> To: MATPOWER Discussion List 
> Subject: Re: Performing N-1 contingency analysis
>
> MATPOWER does not automatically handle contingency analysis, so you will
> have to handle it yourself by removing generators and branches from the
> network one at a time by setting GEN_STATUS and BR_STATUS to zero. Keep
> in mind that you can create islands or isolated buses by taking branches
> out of service. The case_info()
> 
>  and find_islands()
> 
>  can
> helpful with detecting such cases.
>
>Ray
>
> On Jul 3, 2016, at 7:57 AM, Haroon Malik 
> wrote:
>
>  Hello Dr. Zimmerman,
>
> I would want to perform a N-1 contingency analysis on the 118 bus system.
> What is the steps or procedure of applying such analysis using MATPOWER.
>
>
> Many thanks for your time,
>
>
> Mohammad Haroon Malik
>
>
>


rank_contingencies.m
Description: application/vnd.wolfram.mathematica.package

Re: Convergence problems in a Microgrid

[Jose Luis Marín]

I guess the answer is "maybe", as it depends on the reactive power needs of
both lines and loads.  The best way to know is to start solving just
regular powerflows before you attepmpt OPF.  Make all solar panel
generators run as PV-type (use VG=1 for all of them), and without enforcing
Q limits (use powerflow option pf.enforce_q_lims=0).  Look at the Q values
resulting from the solutions you obtain, and you will get a feeling as to
how much reactive power they will need to inject/absorb. That will help you
to bracket the [QMIN, QMAX] values that you need for feasibility.  Also the
values of the power factors that you would need to use when type-switching
them to PQ mode.

I would also make sure that you obtain good powerflows with all solar
panels running as PQ, before attempting OPF, just in case.

-- 
Jose L. Marin
Grupo AIA


2016-05-16 19:26 GMT+02:00 Zé Miguel :

> Hi Jose, thank you for your help.
>
> So, the only change that I need to implement in my .m file is the QMAX of
> each PV generator, according to the specified power factor. Right ?
>
> Cheers,
> José Miguel
>
> --
> From: mari...@aia.es
> Date: Fri, 13 May 2016 17:39:08 +0200
> Subject: Re: Convergence problems in a Microgrid
> To: matpowe...@list.cornell.edu
>
>
>
> You won't even have a mathematically solvable powerflow unless you allow
> at least one generator to inject/absorb reactive power (since you do have
> loads with a reactive component).  PV generation actually can do this,
> through their inverters.  The following article may help:
> http://solarprofessional.com/articles/design-installation/reactive-power-control-in-utility-scale-pv?v=disable_pagination
> Reactive capabilities of about -0.90 to +0.90 power factor are cited in
> there.
>
> --
> Jose L. Marin
> Grupo AIA
>
>
>
> 2016-05-13 17:19 GMT+02:00 Zé Miguel :
>
> Hello,
>
> Thank you for your help Ray. What do you suggest to solve this issue,
> since my test-microgrid only have PV generation ?
>
> Cheers,
> Zé Miguel
>
> --
> From: r...@cornell.edu
> Subject: Re: Convergence problems in a Microgrid
> Date: Thu, 12 May 2016 13:17:13 -0400
> To: matpowe...@list.cornell.edu
>
> It appears that you have no reactive generation capacity to meet any
> reactive demand.
>
> Ray
>
>
> >> case_info('rede_teste')
> Checking connectivity ... single fully connected network
> Elapsed time is 0.007586 seconds.
>
> 
> Full
>System
> Number of:   --
>   buses33
>   loads28
> on 28
> off -
> fixed  28
> dispatchable-
>   on-
>   off   -
>   generators   15
> on 15
> off -
>   shunt elements-
>   branches 32
> on 32
> off -
> ties (off)  -
>
> Load
>   active (MW)
> dispatched  0.2
>   fixed 0.2
>   dispatchable  -
> nominal 0.2
>   on0.2
>   off   -
>   fixed 0.2
>   dispatchable  -
> on  -
> off -
>   reactive (MVAr)
> dispatched  0.1
>   fixed 0.1
>   dispatchable  -
> nominal 0.1
>   on0.1
>   off   -
>   fixed 0.1
>   dispatchable  -
> on  -
> off -
>
> Generation
>   active (MW)
> dispatched  -
> max capacity0.1
>   on0.1
>   off   -
> min capacity-
>   on-
>   off   -
>   reactive (MVAr)
> dispatched  -
> max capacity-
>   on-
>   off   -
> min capacity-
>   on-
>   off   -
>
> Shunt Injections
> active (MW) -
> reactive (MVAr) -
>
> Branch Losses
> active (MW) -
> reactive (MVAr) -
>
> DC line
>   export (MW)
> dispatch-
> max capacity-
>   on-
>   off   -
> min capacity-
>   on-
>   off   -
>
> Reference Buses
>   num of ref buses  1
>   ref bus numbers   1
> >>
>
>
> On May 12, 2016, at 7:51 AM, Zé Miguel  wrote:
>
> Dear all,
>
> I'm trying to 

Re: Change 1 PQ to PV -> infeasible problem

[Jose Luis Marín]

Chris,

Iterative methods in general (whether they are N-R, FD, G-S, etc.) are
prone to failure; you can sometimes mitigate this problem, but you can
never eliminate it 100%. It sounds like you may have hit one of those cases
where NR iteration fails and it's hard to find "the right seed".  Our
company provides tools that help precisely in cases like this, but I don't
want to spam the list.  We can continue this conversation off-list.

-- 
Jose L. Marin
Grupo AIA


2016-05-10 10:40 GMT+02:00 Chris Prokop <christophprok...@gmail.com>:

> Dear Jose and Ray,
>
> thanks for your replies, I've tested another few things:
>
> I deactivated the other two PV-buses, so there is only one reference bus,
> many PQ-buses and the problematic PQ/PV-bus.
>
> With the PQ (result is ok):
> - Slack inserts: -60 Mvar
> - Problematic bus: 0.118 Mvar, V=1.0335 pu
>
> With PV without limits (result is ok)
> - Slack inserts: -61 Mvar
> - Problematic bus: 0.851 Mvar, V=1.05 pu (as set as Vg)
>
> With PV with limits (not ok):
> - Slack inserts: +84 Mvar (145 Mvar difference)
> - Problematic bus: 0.7 Mvar, V=0.15 pu (generator has reached its upper
> limit, hence it was converted to PQ)
>
> But the result of V of the problematic bus should be somewhere in between
> 1.0335 and 1.05 p.u.. All other P and Q of the generators and buses remain
> constant within all variants. For me it looks like the result is in the
> unstable area of the QV-curve.
>
> @ Jose: when setting Vg=1.0335 than the result is ok again. The problem
> starts at the point the generator reaches it's Q-Limit and gets converted
> to PQ again.
>
> I don't understand why the result is in the unstable area of the QV-curve,
> because the reactive power is capacitive and not inductive, hence I
> wouldn't expect a stability Problem (I'm increasing the voltage).
>
> Interesting point: with the fast-decoupled solver the power flow with
> Q-limit converges to V=1.0465 p.u. (result is ok), with Gauss-Seidel it
> doesn't converge. Maybe there is a problem with Newton's method?!
>
> Thanks for your help/explanations, nice regards,
> Chris
>
>
> 2016-05-06 7:33 GMT+02:00 Jose Luis Marín <mari...@aia.es>:
>
>>
>> Sorry for the confusion, Ray:  where I said, *"...solvable with the bus
>> running as PV, but unfeasible when the bus is running as PQ"*, I meant
>> to say:
>>
>>"*a valid powerflow solution* with the bus running as PV, but an
>> *unstable* one when the bus is running as PQ"
>>(it is one of the multiple low-voltage solutions of the full
>> mathematical problem--these solutions are unstable form the point of view
>> of real operations)
>>
>> So yes, in that situation the solution obtained with the bus type being
>> PV is (of course) still a mathematical solution to the same problem,
>> regardless of the bus type being reinterpreted as PQ.  But the point is
>> that, only if you look at the problem in this second case (bus type
>> switched to PQ) it is then possible to see that the solution is a low
>> voltage one, and therefore not valid.  One way to detect this is to look at
>> the slope dQ / dV of the Q-V curve, because this situation corresponds to
>> being "on the wrong side" of the curve.  It may happen either if the
>> setpoint V was too low, or if the P output is too high.
>>
>> What the Newton method may do in that case is hard to say, since N-R does
>> not "see" the kind of instability that we are talking about here.  It may
>> converge to this unstable solution, or to the operationally correct one, or
>> diverge.  It all depends on the complex landscape of the respective basins
>> of attraction for each solution, which as you know is fractal.  But unless
>> the case is right at the bifurcation point (here that would be the minimum
>> of the Q-V curve), the basin of attraction around this unstable solution
>> will be finite, so there's a good chance that if you provide the original
>> solution (obtained when the bus was type PV) and use it as the seed to the
>> new problem (bus type now switched to PQ), then NR will converge to the
>> same solution.  You'd probably need to start from a flat start in order to
>> find the high-voltage solution instead.
>>
>> --
>> Jose L. Marin
>> Grupo AIA
>>
>>
>> 2016-05-05 22:17 GMT+02:00 Ray Zimmerman <r...@cornell.edu>:
>>
>>> Thanks Jose. I’ve added that to the manual as you suggested.
>>>
>>> I also agree with your suggestions for Chris. However, I’m curious about
>>> your [*] note. In the situation you describe, it’s just that the Newton

Re: A question about transformer modeling with uk, ur

[Jose Luis Marín]

Certainly, testing with a small system is the only way to find out what the
internal model *really* is (since you can't inspect the code in this
case).  Anyway, I would add another suggestion: export from Netomac to
PSS/E RAW format, and inspect the transformer records very carefully.  Read
the specs of the RAW format (the specific version of the format used), and
have a look at the internal model used by PSS/E (I'd assume that Netomac
uses the same, but who knows).  You can find it here:


http://w3.usa.siemens.com/datapool/us/SmartGrid/docs/pti/2009July/PDFs/Modeling_of_two_winding_voltage_regulating_transformers.pdf

Compare this model with the one use by MATPOWER (page 22 of the PDF
manual).  You'll see they're not exactly the same.  You said you're not
using magnetizing branches (Bs), which simplifies things a lot.  But pay
attention to the tap ratios (both sides!), just in case.

-- 
Jose L. Marin
Grupo AIA



2016-05-10 13:50 GMT+02:00 Chris Prokop :

> Dear Shuo,
>
> I'd try to simulate a simply network with 1 reference bus with 1
> transformer and 1 subsequent load to test the behaviour between Matpower
> and netomac.
>
> I'm not used to netomac, so can't help here directly.
>
> Nice regards,
> Chris
>
> 2016-05-10 13:22 GMT+02:00 Shuo Chen :
>
>> Dear Chris,
>>
>> thanks for your reply, actually we are getting the same transformer data.
>> In my power system  uk is far more larger than ur (more than 50:1), so it
>> can be assumed that ux = uk.
>>
>> i'm wondering whether there is any little difference between netomac and
>> matpower by modeling the transformer, so that even if I give them the same
>> input, PF results could still differ.
>>
>> best regards
>> Shuo
>>
>> Zitat von Chris Prokop ：
>>
>> Dear Shuo Chen,
>>>
>>> I'm used to similar data and calculate r and x as (ur, uk in p.u.,
>>> base_MVA
>>> and S_transformer in MVA):
>>> - r = ur * base_MVA / S_transformer
>>> - x = sqrt(uk^2 - ur^2) * base_MVA / S_transformer
>>> - ratio = 1 (in your case, where rated = nominal voltage)
>>> - angle = 0 (in your case, without phase shifting)
>>>
>>> Maybe one of your problems is using uk instead of ux (uk^2 = ux^2 +
>>> ur^2)?
>>>
>>> Nice regards,
>>> Chris
>>>
>>> 2016-05-10 0:33 GMT+02:00 Shuo Chen :
>>>
>>>
 Dear Ray and matpower users,

 i'm writing my thesis about a data-converter for two power system
 simulation softwares: PSSE@Netomac from Siemens and matpower. Here is a
 problem with transformer modeling, i've read a lot in this archive but
 still can't solve it, so i decide to post my question here.

 My goal is to get the same or a simular power flow convergency using the
 Newton's method after converting the net topology from one to the other.
 However, I'm stucked when I try to convert a Netomac net into a matpower
 case. Comparing the PF results of the 2 softwares, there is always a big
 deviation (by bus voltage magnitude, bus voltage angle and branch P/Q
 injection). More specifically, the bus voltages in pu of matpower are
 all
 lower than those of Netomac, the deviation could be up to 5%, like 0.977
 instead of 0.997 in Netomac.

 I simply grab all the bus/gen/branch parameters from Netomac except for
 the transformer impedance, so I guess there might be a mistake when the
 transformer model is built in matpower.
 The transformers have 2-windings and no tap changer

 The parameters I could get from Netomac are:

 - Rated voltage of HV side winding (Un1)
 - Nominal network voltage HV side (UB1? the value is equal to Un1)
 - Rated voltage of LV side winding (Un2)
 - Nominal network voltage LV side (UB2? the value is equal to Un2)
 - Rated apparent power (equal to baseMVA)
 - ur in % (may not be zero)
 - uk in % (must be larger than ur)
 - P0 and I0 are not given
 - vector group YY0
 in matpower a transformer is treated as a transmission line, where
 r(p.u.)
 and x(p.u.) are needed for the power flow calculation. (b is omitted in
 my
 case)
 According to the formels:

 z_pu = uk% / 100
 x_pu = uk% *(Un/UB).^2 * (SB/Sn) / 100
 r_pu = sqrt(z_pu.^2 - x_pu.^2) (here r_pu = ur%/ 100)

 Since in my case Un/UB = 1, SB/Sn = 1, so x = uk%/100, r = ur%/100, b =
 0
 I set ratio = 1 and angle = 0, status = 1, angmin/max = +-360
 other values are set to 0

 Netomac can also export a .raw file for PSSE, the conversion above will
 generate equivalent r and x values as the exported .raw file.
 Theoretically, with the same net topology and P/Q accuracy, the PF
 convergency of the two simulation tools should be almost the same.

 Here is one more hint, as far as I know, the transformer impedance in
 matpower is modeled at "to" side, while in Netomac it's modeled at "HV"
 side 

Re: Change 1 PQ to PV -> infeasible problem

[Jose Luis Marín]

Sorry for the confusion, Ray:  where I said, *"...solvable with the bus
running as PV, but unfeasible when the bus is running as PQ"*, I meant to
say:

   "*a valid powerflow solution* with the bus running as PV, but an
*unstable* one when the bus is running as PQ"
   (it is one of the multiple low-voltage solutions of the full
mathematical problem--these solutions are unstable form the point of view
of real operations)

So yes, in that situation the solution obtained with the bus type being PV
is (of course) still a mathematical solution to the same problem,
regardless of the bus type being reinterpreted as PQ.  But the point is
that, only if you look at the problem in this second case (bus type
switched to PQ) it is then possible to see that the solution is a low
voltage one, and therefore not valid.  One way to detect this is to look at
the slope dQ / dV of the Q-V curve, because this situation corresponds to
being "on the wrong side" of the curve.  It may happen either if the
setpoint V was too low, or if the P output is too high.

What the Newton method may do in that case is hard to say, since N-R does
not "see" the kind of instability that we are talking about here.  It may
converge to this unstable solution, or to the operationally correct one, or
diverge.  It all depends on the complex landscape of the respective basins
of attraction for each solution, which as you know is fractal.  But unless
the case is right at the bifurcation point (here that would be the minimum
of the Q-V curve), the basin of attraction around this unstable solution
will be finite, so there's a good chance that if you provide the original
solution (obtained when the bus was type PV) and use it as the seed to the
new problem (bus type now switched to PQ), then NR will converge to the
same solution.  You'd probably need to start from a flat start in order to
find the high-voltage solution instead.

-- 
Jose L. Marin
Grupo AIA


2016-05-05 22:17 GMT+02:00 Ray Zimmerman :

> Thanks Jose. I’ve added that to the manual as you suggested.
>
> I also agree with your suggestions for Chris. However, I’m curious about
> your [*] note. In the situation you describe, it’s just that the Newton
> method will diverge, correct? The original solution will still be a
> solution of the power flow equations won’t they?
>
> Ray
>
>
> On Apr 30, 2016, at 2:01 PM, Jose Luis Marin 
> wrote:
>
>
> Interesting, I didn't know that enforcing Q-limits also affects generators
> when their bus-type is set to PQ.  Ray, I suggest documenting this behavior
> in the manual, probably at the end of the last paragraph in Section 4.1.
> Something to the effect of *"... Note also that this option affects
> generators even if the bus they are attached to is already of type PQ."*
>
> Going back to Chris's problem, I suggest you approach this as consisting
> of two orthogonal issues:
>
>- Solving with Q-limits enforced vs. solving without
>- Solving with a given generator running as PQ-type vs. running as
>PV-type (by manually switching it)
>
> With regards to the first issue, I suggest to start by analyzing the
> behavior of your case *without* enforcing Q limits: in particular, pay
> close attention to the QG injections obtained *in the solution* for those
> gens operating as PV-type (and also the specified QG for those operating as
> PQ), and *compare* those to your specified [QMIN, QMAX] values.  If you
> observe any large violation, then solving this case with Q limits enforced
> will yield a *very* different solution.
>
> Concerning the second issue, the main thing to check when manually
> switching a given bus from PQ-type to PV-type, is that the gen setpoint VG
> has to be set to the bus voltage VM obtained in the previous solution,
> otherwise your next solution could be again very different.  Conversely,
> when manually switching a bus from PV-type to PQ-type [*], you have to
> specify a gen injection QG equal to the value obtained in the previous
> solution.
>
> [*] As I mentioned before, the switch in this direction is not always
> guaranteed to give you the same solution.  If the setpoint VG is decreased
> too much, or if the real power PG is increased too much, you can end up
> with a case that is solvable with the bus running as PV, but unfeasible
> when the bus is running as PQ (under the switching procedure described
> here).  Technically: this happens when you're running the PV generator at
> the unstable branch of its Q-V curve.
>
> --
> Jose L. Marin
> Grupo AIA
>
>
>
> On Fri, Apr 29, 2016 at 10:12 AM, Chris Prokop  > wrote:
>
>> Hi,
>>
>> oh, I just realized that the Q-limits also hold on PQ-buses, hence the
>> generators are limited to it - I wouldn't have thought that. So without
>> Q-limits the results are the same - huch.
>>
>> Anyway, when changing the bus to a PV-bus I've got the problem that the Q
>> of the PV-buses (generators) are calculated way to