Open loop losses in three winding transformers

[WinfriedL]

The three winding transformer in pandapower models open loop losses/iron losses as the losses of one of the three two winding transformers. Is this a common approach? The models I know use a shunt element at the auxiliary bus. There is a difference between both methods depending on the values of pfe_kw and i0_percent and the voltage drop across the transformer.

[lthurner]

I don't remember exactly why we did it this way, but I think it was just easier to give it as i0_percent and pfe_kw to the HV transformer, because the functions were already implemented. But I agree that the star point would be the better place to model the open loop losses.

I looked at the values in the trafo3w test and there is actually a deviation of ~0.015kVA at the HV side while the deviation at the low voltage sides is <1e-11 kVA. This might be the effect of the open loop losses which are not modeled in the star point. The deviation is not big enough to let the tests fail in this case, but the effect could of course be larger for different transformers. Would definitely be interesting to see if modeling the open loop losses in the start point gets rid of this deviation!

[WinfriedL]

When we first calculated three winding transformers with pypower we tried splitting i0_percent and pfe_kw evenly to all three transformers. But I guess this is as questionable as the current approach. Now I think the only straightforward and unambiguous way is by putting it at the star point. I don't know which method is used in your test data, maybe you could try to increase pfe_kw and i0_percent? If the deviation gets larger it would be a clear sign that the current model is the problem?

Of course the current approach by reusing the two winding transformer functionality is easier at first sight, but maye in the end the code in _branch_df_from_trafo3w() gets cleaner since you don't have to distinguish between the three transformers anymore?

I don't know the code well enough to judge how complex adding a shunt element at the star point would be though.

[lthurner]

Of course the current approach by reusing the two winding transformer functionality is easier at first sight, but maye in the end the code in _branch_df_from_trafo3w() gets cleaner since you don't have to distinguish between the three transformers anymore?

This function is not currently used, it was an attempt at refactoring the 3W-transformer creation as discussed in #54, but it was never finished. I didn't realise these functions are still there, we should comment out or remove them as long as no one is working on this.

The function that is actually being used is _trafo_df_from_trafo3w. The losses are passed to the HV transformer here. Adding a shunt at the auxiliary bus instead should not be very difficult.

[lthurner]

I tried it out, but the voltages are now further from the powerfactory result than before. I added this code block in build_bus.py at line 400:

    trafo3w = net["trafo3w"]
    if len(trafo3w) > 0:
        vn_lv = get_values(ppc["bus"][:, BASE_KV], trafo3w["ad_bus"].values, bus_lookup)
        baseR = np.square(vn_lv) / net.sn_kva * 1e3
 
        vnl_squared = trafo3w["vn_hv_kv"].values ** 2
        b_real = trafo3w["pfe_kw"].values / (1000. * vnl_squared) * baseR
        i0 = trafo3w["i0_percent"].values
        pfe = trafo3w["pfe_kw"].values
        sn = trafo3w["sn_hv_kva"].values
        b_img = (i0 / 100. * sn / 1000.) ** 2 - (pfe / 1000.) ** 2
    
        b_img[b_img < 0] = 0
        b_img = np.sqrt(b_img) * baseR / vnl_squared
        q = np.hstack([q, b_img])
        p = np.hstack([p, b_real])
        b = np.hstack([b, trafo3w["ad_bus"].values])

which adds the susceptance and conductance at the auxiliary bus. Then I replaced ttab.pfe_kw and ttab.pfe_kw.i0_percent with 0 in build_branch.py. The open loop losses that the transformer then absorbs stay exactly the same, but the voltage difference to PowerFactory rises from ~1e-7 to >1-6. Feel free to try this in your case and see if results get better with this solution or if you can find a mistake in this implementation?

[WinfriedL]

Just a quick update, did not test much so far....

I have a testgrid with one 3 winding transformer and an identical grid where I did a manual repesentation by three 2winding transformers and an explicitly created shunt element at the star point.
The manual representation gives the "correct" results. Your implementation with the three winding transformer not. So there seems to be a mistake (either in your code or in my copy&pasting) but I did not have enough time so far for investigating. Will do it later.

[WinfriedL]

I am not sure, but I think you are calculating the wrong parameters...
It seems like you are calculating the shunt admittance, but pypower expects the consumed power in MW/Mvar.

Here is my implementation which works at least for my testcase:

trafo3w = net["trafo3w"]
if len(trafo3w) > 0:
    pfe_kw = trafo3w["pfe_kw"].values
    i0 = trafo3w["i0_percent"].values
    sn_kv = trafo3w["sn_hv_kva"].values
    
    q_kvar= np.sqrt( (sn_kv * i0 / 100.) ** 2 - pfe_kw **2 ) 
    
    q = np.hstack([q, q_kvar / np.float64(1000.)])
    p = np.hstack([p, pfe_kw / np.float64(1000.)])
    b = np.hstack([b, trafo3w["ad_bus"].values])

This should work as long as the auxiliary bus voltage equals vn_hv_kv , but I have not tested that in depth.

[WinfriedL]

While it worked for my very simple testcase, it did not work for a real three winding transformer. I probably did a simplification along the way which I did not take into account. I will look further...

[WinfriedL]

I think I found the problem. I assumed that the auxiliary bus gets the vn_hv_kv from the transformer, but in fact the voltage is equal to to bus.vn_kv where the transformer is connected to. So I have to rescale the shunt parameters. My current implementation:

if len(trafo3w) > 0:
    pfe_kw = trafo3w["pfe_kw"].values
    i0 = trafo3w["i0_percent"].values
    sn_kv = trafo3w["sn_hv_kva"].values
    
    q_kvar= np.sqrt( (sn_kv * i0 / 100.) ** 2 - pfe_kw **2 ) 
    
    vn_hv_trafo=trafo3w["vn_hv_kv"].values
    
    hv_buses = net.bus.loc[net.trafo3w.hv_bus.values]
    vn_hv_bus=hv_buses.vn_kv.values

    v_ratio = ( vn_hv_bus / vn_hv_trafo) ** 2
    
    q = np.hstack([q, q_kvar / np.float64(1000.) * v_ratio])
    p = np.hstack([p, pfe_kw / np.float64(1000.) * v_ratio])
    b = np.hstack([b, trafo3w["ad_bus"].values])

This seems to work in my testcases.

Few comments:
Is there a reaseon why the hv bus voltage is used for the auxiliary bus? It would seem to me that the vn_hv_kv from the transformer would be more straightforward...

In _calc_shunts_and_add_on_ppc, it might be usefull to add a comment explaining those lines:
https://github.com/lthurner/pandapower/blob/fbdf93fcd56d7709076130e67e469b3b0b722370/pandapower/build_bus.py#L379

It took me some time to understand what you are doing there...

[lthurner]

I am not sure, but I think you are calculating the wrong parameters...
It seems like you are calculating the shunt admittance, but pypower expects the consumed power in MW/Mvar.

In per unit values, that is actually the same parameter. From s=u * i and i=u * y it follows that s=u^2*y. That means that at rated voltage (u=1), it is s=y. So S can be calculated through y.

This seems to work in my testcases.

How do you determine what is working? Are you testing against a different tool, if so which tool?
Here is the test code that I am using with the grid from the pandapower test suite:

from pandapower.test.conftest import result_test_network
net = result_test_network()
v_tol=1e-6
i_tol=1e-6
s_tol=2e-2
l_tol=1e-3
buses = net.bus[net.bus.zone == "test_trafo3w"]
trafos = [x for x in net.trafo3w.index if net.trafo3w.hv_bus[x] in buses.index]
b2 = buses.index[1]
b3 = buses.index[2]
b4 = buses.index[3]
t3 = trafos[0]

uhv = 1.010117166
umv = 0.955501331
ulv = 0.940630980

load = 37.21
qhv = 1.64375
qmv = 0
qlv = 0

ihv = 0.00858590198
imv = 0.20141269123
ilv = 0.15344761586

phv = 300.43
pmv = -200.00
plv = -100.00


assert abs((net.res_bus.vm_pu.at[b2] - uhv)) < v_tol
assert abs((net.res_bus.vm_pu.at[b3] - umv)) < v_tol
assert abs((net.res_bus.vm_pu.at[b4] - ulv)) < v_tol

assert abs((net.res_trafo3w.loading_percent.at[t3] - load)) < l_tol

assert abs((net.res_trafo3w.p_hv_kw.at[t3] - phv)) < s_tol
assert abs((net.res_trafo3w.p_mv_kw.at[t3] - pmv)) < s_tol
assert abs((net.res_trafo3w.p_lv_kw.at[t3] - plv)) < s_tol

assert abs((net.res_trafo3w.q_hv_kvar.at[t3] - qhv)) < s_tol
assert abs((net.res_trafo3w.q_mv_kvar.at[t3] - qmv)) < s_tol
assert abs((net.res_trafo3w.q_lv_kvar.at[t3] - qlv)) < s_tol

assert abs((net.res_trafo3w.i_hv_ka.at[t3] - ihv)) < i_tol
assert abs((net.res_trafo3w.i_mv_ka.at[t3] - imv)) < i_tol
assert abs((net.res_trafo3w.i_lv_ka.at[t3] - ilv)) < i_tol

This is testing the results of a 3W-transformer against results from PowerFactory. This test for me only passes with the current pandapower implementation. Every other implementation that I tried (the code that I posted yesterday and the code that you posted now) gives deviations from the PowerFactory values outside of the tolerances. Note that this transformer includes a tap changer with position outside of the neutral point and rated voltage values that are different from the rated bus voltages. I tried to define the most complicated possible case so that the test covers all parameters.

Is there a reaseon why the hv bus voltage is used for the auxiliary bus? It would seem to me that the vn_hv_kv from the transformer would be more straightforward...

These details are very difficult in the modeling process, because they are rarely described anywhere in literature with enough detail to actually implement such a model. In this case what I did is just try both possibilities and see what leads to the same result as in other (commercial) tools. Of course other tools are also not perfect, but there would have to be a good reason to deliberately implement a model that does not match up with different tools. In this case, using vn_hv_kv from the transformer would be more straightforward, but in my tests it does lead to different voltage results than PowerFactory, so that is why I implemented it this way. I don't think I tested it against other tools, so if other tools have different results, it would be a valid discussion what the correct implementation is that we should follow.

In _calc_shunts_and_add_on_ppc, it might be usefull to add a comment explaining those lines:

That is really a little unclear, could be refactored for better readability...

[WinfriedL]

In per unit values, that is actually the same parameter.

Right. It looked overly complicated and in that complexity there seemed to be an error somewhere....

I test against Integral, using a simple three bus grid with one three winding transformer. I have another Integral test grid where I replaced that three winding transformer with three two winding transformers and a shunt element at the star point.
I get the same results with both grids in Integral and I get the same results for both grids with pandapower and my implementation.

I have another larger grid with several three winding transformers, lines, impedances and the power deviations between Integral and pandapower are less than 0.05kw with my implementation and much larger with the original implementation.

That's what I meant with "work in my testcases."

Those transformers are all 380/220/30 kV Transformers without a tap changer.

These details are very difficult in the modeling process, because they are rarely described anywhere in literature with enough detail to actually implement such a model.

Indeed....

In this case what I did is just try both possibilities and see what leads to the same result as in other (commercial) tools.

I think it does not really matter which voltage level is used, as long as all parameters are correctly rescaled.... and that should only matter for the open loop losses. And maybe tap changers, I have not yet considered those for three winding transformers. Maybe that is the reason why the pandapower testcase fails (for me too).

I don't think I tested it against other tools, so if other tools have different results, it would be a valid discussion what the correct implementation is that we should follow.

Now you have an Integral test result ;)

How difficult would it be to generate data for your testcase with different tools? I could try to do it with Integral...

[lthurner]

Right. It looked overly complicated and in that complexity there seemed to be an error somewhere....

Thats because I copy pasted it from _calc_y_from_dataframe, which is already a tested version of calculating y from i0_percent and pfe_kw, so it was easier to just use this than to implement a new function calculating S directly.

I have another larger grid with several three winding transformers, lines, impedances and the power deviations between Integral and pandapower are less than 0.05kw with my implementation and much larger with the original implementation.

Interesting. I would have expected there are some differences between the different tools, but not that they would be so significant. But if we now have an implementation that matches up with PowerFactory and an implementation that matches up with Integral, we can of course also introduce a parameter in runpp that allows to switch between the different implementations. We have already done this for the 2W-transformer, where it is possible to switch between the pi and the T model. But if we want to include the model, the new version should also be tested with tap-changer, because that is usually where it really gets tricky.

How difficult would it be to generate data for your testcase with different tools? I could try to do it with Integral...

If you can replicate the test case in Integral and supply the results, that would be helpful for a comparison. Or propose a different use case if you think the one we have does not cover everything?

I could try to replicate the test case in NEPLAN and Sincal. However, Sincal uses a simple pi model for 2W-transformers even though transformers are T-shaped, so I am not sure how relevant their transformer modeling is. It would also be interesting to look into the documentation of those tools to see if there are differences in the models. I can't promise I will get around to this soon though...

[WinfriedL]

If you can replicate the test case in Integral and supply the results, that would be helpful for a comparison.

It took me much longer than expected since I did some stupid mistakes implementing the testcase in Integral, but here it is:

#Test results Integral:
uhv = 1.01011711678
umv = 0.95550024145
ulv = 0.94062989256

load = 37.209
qhv = 1.660
qmv = 0
qlv = 0

ihv = 0.00858591110
imv = 0.20141290445
ilv = 0.15344776975

phv = 300.43
pmv = -200.00
plv = -100.00

When I put this data in test_results.py all tests pass with my implementation with a shunt at the star point.

So we have indeed two different approaches from commercial tools and we could either match one with the current implementation or the other with the new implementation.

BTW the parameters for the out of service transformer are not plausible, since pfe_kw is too large in comparison to i0_percent. This leads to a negative value in the squareroot of my implementation...

I see several options at the moment:

• ignore the Integral results (not good for me, since then I would have to maintain a fork from pandapower)
• implement both approaches and let the user select. (Fine with me, but adds complexity to the code)
• adopt the new implementation, maybe with some more tests, ideally a third commercial software which backs that up.

That is of course all up to you.

From all discussions we had internally we find the Integral approach more intuitive and from a physical standpoint more in line with the T-model for a two winding transformer. And while the testcase fails against the PowerFactory values, the deviations are still extrem small in my opinion, 1.08e-6 > 1e-6 if recall correctly.

On a side note: I am not convinced PowerFactory has the same implementation as you have at the moment. I don't have access to PowerFactory, but from a technical report I have seen in the past somewhere there were several shunt elements, neither only at the star point, as Integral does it, nor only at the hv transformer as you do it.

[lthurner]

Thanks for the values thats really helpful! So we have now discussed this in the pandapower developers meeting, and there was consensus that if we now have two versions that are validated, we prefer keeping both of them as options, even if it makes the code a little more complex.

As you mentioned, the deviation in the results is pretty small, so I don't think the differences matters much for the electric results of studies. It does however matter for validation purposes - even if the differences are small, one can never be sure if the difference really comes from the transformer modelling or from some other error in the conversion. That is why we prefer maintaining both implementations.

To do that, we have to

• introduce a new parameter in runpp and ensure it is stored in net._options
• check the parameter in build_bus and build_branch to switch between the models
• include the Integral parameters in the test and run two powerflows with the different parameters checking agains each set of results

This shouldn't be much work to implement if you already have the transformer implementation ready. Do you want to do this?

BTW the parameters for the out of service transformer are not plausible, since pfe_kw is too large in comparison to i0_percent. This leads to a negative value in the squareroot of my implementation...

I think this was intentional to test the backstop we have against this (although an explicit seperate test would have been better). The imaginary part of the impedance is limited to zero here, so that the square root canot be negative.

From all discussions we had internally we find the Integral approach more intuitive and from a physical standpoint more in line with the T-model for a two winding transformer. And while the testcase fails against the PowerFactory values, the deviations are still extrem small in my opinion, 1.08e-6 > 1e-6 if recall correctly.

On a side note: I am not convinced PowerFactory has the same implementation as you have at the moment. I don't have access to PowerFactory, but from a technical report I have seen in the past somewhere there were several shunt elements, neither only at the star point, as Integral does it, nor only at the hv transformer as you do it.

As far as I know, PowerFactory does not use the equivalent branch modeling with three 2W-transformers, but uses a 3-pole equivalent circuit model. Modeling the 3W-transformer as an actual 3-pole should generally be more accurate than the equivalent branch modeling that is used by Integral, pandapower or Sincal. So any modeling we use with three 2W-transformers will always only be an approximation of the PowerFactory model. Even though I agree the Integral approach looks more straightforward, our current equivalent branch implementation seems to come closer than any other one to the PowerFactory 3-pole modeling. Always keeping in mind however, that the differences are so small that it probably doesn't really matter except for validation with the different tools.

[lthurner]

Update: we have discussed about this issue again today, and checked PowerFactory again. We discovered that PowerFactory provides an option where the losses are supposed to be modeled:

I was not aware this a parameter when I made the test cases, and apparently the losses were modeled at the HV side by default. So it seems like we matched the implementation with the PowerFactory HV option, but it could also be modelled differently.

With this new information I would propose making a parameter "trafo3w_losses" that can be "HV", "MV", "LV" or "star". Then the loss parameters would either be passed to the respective equivalent transformer or be modeled as a shunt element in the start point. This should cover all implementation variants?

[WinfriedL]

This shouldn't be much work to implement if you already have the transformer implementation ready. Do you want to do this?

In principle yes, but not sure when I can do it. My list of todos is rather long and I have a working version for my testcases...

I think this was intentional to test the backstop we have against this

This test is missing in the new shunt-at-starpoint implementation, so I have to add it.

We discovered that PowerFactory provides an option where the losses are supposed to be modeled:

Actually, Integral has the same option. But at least for Integral it does not change the actual equivalent model of the transformer, but only on which bus (and thereby in which area of the grid) the losses are accounted for.

[lthurner]

In principle yes, but not sure when I can do it. My list of todos is rather long and I have a working version for my testcases...

Ok, it's actually not a big deal to do this for me, but I would prefer not copying in your code and taking credit for it. I have now introduced the parameter trafo3w_losses in runpp. The loss parameters are now handed to the HV, MV or LV transformer depending on this parameter value in build_branch. HV is the default value, so that the behaviour does not change. If the parameter is defined as "star", a NotImplemented Error is raised in build_bus. Can you copy in your current implementation there and the Integral results in the test that is checked against a power flow with trafo3w_losses="star"?

Actually, Integral has the same option. But at least for Integral it does not change the actual equivalent model of the transformer, but only on which bus (and thereby in which area of the grid) the losses are accounted for.

Weirdly enough, I just tested it in PowerFactory and the power flow result doesn't change at all when I select a different location for the losses, not even in the 8th or 9th decimal point. Not sure if some other option overrides this or what is happening, but there should be some notable difference...

[WinfriedL]

Weirdly enough, I just tested it in PowerFactory and the power flow result doesn't change at all when I select a different location for the losses, not even in the 8th or 9th decimal point. Not sure if some other option overrides this or what is happening, but there should be some notable difference...

No, that is exactly what I meant. That parameter does not change anything in the calculation. But after the powerflow is calculated, this parameter determines in which area the losses are accounted for.

[lthurner]

Oh ok, now I understand what you meant. That is not what the PowerFactory Technical Reference indicates though:

[WinfriedL]

That's the model I meant earlier with shunt elements at all buses. I don't really get how this is connected to the real implementation. I tried to implement this kind of model in pandapower but either I got the wrong distribution of r/x-values for the single shunt elements or that diagram does not represent the real implementation at all.

While it is possible that the position of the shunts is selectable by the mentioned PowerFactory setting "Verlustzuordnung", from the similarity to Integral and your observations I don't think that those are connected.

[lthurner]

Ok I found the right parameter:

This one actually changes the electric behaviour. I am going to analyse a little more if I can replicate this. Can you include your single-shunt element implementation in build_bus and put in a pull request?

[WinfriedL]

Can you include your single-shunt element implementation in build_bus and put in a pull request?

Yes, hopefully I will do that today.

That PowerFactory parameter raises the question, if there should only a shunt element at one bus. If I recall correctly, you are using a pi-model (at least in the testcase) for the hv transformer, this would add a shunt at the hv bus and the auxiliary bus, right?

Out of curiosity, what are the possible options for that parameter?

[lthurner]

That PowerFactory parameter raises the question, if there should only a shunt element at one bus. If I recall correctly, you are using a pi-model (at least in the testcase) for the hv transformer, this would add a shunt at the hv bus and the auxiliary bus, right?

Yes, the option that was chosen by default is "star point", and the best fit for that option seems to be to give the losses to the HV transformer, which splits them between star point and HV side. No idea why actually putting them at the star point gives "worse" result. I am guessing it just fits the 3-pole circuit better in equivalent 2-pole branches by chance?

Out of curiosity, what are the possible options for that parameter?

HV, MV, LV or star point

[KemperJan]

Please find attached some results, a brief summary and some hints by using the test case in NEPLAN and PSS Sincal.
In my opinion the pi-model of the three winding transformer at the HV-side is applicable in both cases, as shown in the results.

Results of the test case by comparing PSS Sincal with Pandapower:
Vergleich_Ergebnisse_Sincalmodel_PandaPow_UPDATE.xlsx
Results of the test case by comparing NEPLAN with Pandapower:
Vergleich_NEPLANmodel_PandaPow_UPDATE.xlsx
Hints:
Beachtung bei der Eingabe von Leistungsangaben aus.docx
Comments/Summary on the 3W-Trafomodell:
Kommentare_zu_den_Anpassungen_der_Trafomodelle.docx
The the test case implemented in NEPLAN & PSS Sincal:
Testfiles_NEPLAN_SINCAL.zip

.