Automatized short-circuit analysis with Python in DIgSILENT PowerFactory
Sometimes you want to know the short-circuit currents and short-circuit powers at different locations in the power system. This could be of importance for different reasons. For example, design of protection systems or identification of weak grid areas in terms of short-circuit power. Therefore, in this post, I present an approach for automatized short-circuit analysis with Python in DIgSILENT PowerFactory. The code works on any PowerFactory model and conveniently calculates the short-circuit current/power at all buses of the network.
To showcase the functionality of the analysis I use the good old 39-bus New England system as shown below. The New England system is available in the example projects in PowerFactory.
Automatized short-circuit analysis in DIgSILENT PowerFactory
The approach for the automatized short-circuit analysis makes use of the ‘ComShc’ object in PowerFactory. It is the short-circuit calculation module that offers many options in terms of the calculation method, fault type, fault impedance, etc.
In the provided script, you can change the following short-circuit parameters:
- calculation method (I have only tested the method according to VDE 0102 Part 0 / DIN EN 60909-0, no guarantee that the other methods work as well)
- fault type: 3-phase-fault, 2-phase-fault, etc.
- fault impedance in terms of resistance and reactance
Import libraries and PowerFactory API
As always, at the beginning of the script definition of the needed libraries and functions. Here we need pandas and matplotlib. Additionally, we load the PowerFactory API.
Important: don’t forget to replace the path (shown in red) with the path to the Python directory of your PowerFactory installation as indicated.
""" firstname.lastname@example.org www.thesmartinsights.com """ import pandas as pd import matplotlib.pyplot as plt import sys sys.path.append(r'C:\Program Files\DIgSILENT\PowerFactory 2021 SP1\Python\3.8') if __name__ == "__main__": import powerfactory as pf app = pf.GetApplication() if app is None: raise Exception('getting Powerfactory application failed')
Core of the automatized short-circuit analyis
Below you can find the core of the analysis. First, you need to define the project name and the study case. Then you have the possibility to set a voltage range for the evaluation buses. Only buses within this voltage range (and in-service) are considered.
The code below contains detailed comments for each step. So, if you want to change something, it should be easy to find.
#%% main #define project name and study case projName = '_TSI_39_Bus_New_England_System_SC' study_case = '01_StudyCase.IntCase' #evaluate buses within line-to-line voltage range Vbusmin, Vbusmax = 0, 1000 # (kV) #activate project project = app.ActivateProject(projName) proj = app.GetActiveProject() #get the study case folder and activate project oFolder_studycase = app.GetProjectFolder('study') oCase = oFolder_studycase.GetContents(study_case) oCase.Activate() #get all the buses oBuses = app.GetCalcRelevantObjects('.ElmTerm') #select buses that are within specified voltage range and in-service oBuseseval =  for oBus in oBuses: Vbus = getattr(oBus, 'uknom') #bus voltage if (Vbusmin <= Vbus <= Vbusmax) & (getattr(oBus, 'outserv') == 0): oBuseseval.append(oBus) numbuseval = len(oBuseseval) print('#%i buses to evaluate' %numbuseval) #variable pre-allocation Bus_col =  Skss_col =  Ikss_col =  Vbus_col =  #loop through all buses (that are within specified voltage range) for oBus in oBuseseval: # Short-circuit settings oSC = app.GetFromStudyCase('ComShc') #get SC object setattr(oSC, 'iopt_mde', 0) # SC calculation method: 0 = VDE 0102 Part 0 / DIN EN 60909-0 (six other methods avail. -> check in SC calculation object in PF) setattr(oSC, 'iopt_shc', '3psc') # SC type: 3psc = 3-phase, 2psc = 2-phase, spgf = single-phase-to-ground, 2pgf = 2-phase-to-ground (more options -> check in SC calculation object in PF) setattr(oSC, 'Rf', 0) # (Ohm) fault resistance setattr(oSC, 'Xf', 0) # (Ohm) fault reactance setattr(oSC, 'iopt_allbus', 0) # fault location option: 0 = User Selection, 1 = Busbars and Junction Nodes, 2 = All Busbars setattr(oSC, 'shcobj', oBus) # fault location -> bus/terminal object oSC.Execute() #execute load flow #get variables of interest Skss = getattr(oBus, 'm:Skss') # (MVA) short-circuit power Ikss = getattr(oBus, 'm:Ikss') # (kA) initial symmetrical short-circuit current Vbus = getattr(oBus, 'uknom') # (kV) voltage level of bus #put together important data Bus_col.append(getattr(oBus, 'loc_name')) # bus name Skss_col.append(Skss) # (MVA) short-circuit power Ikss_col.append(Ikss) # (kA) initial symmetrical short-circuit current Vbus_col.append(Vbus) # (kV) voltage level of bus #Results: bus names, SC power, SC current, bus voltage RES = pd.DataFrame() RES['name'] = Bus_col RES['Skss'] = Skss_col RES['Ikss'] = Ikss_col RES['Vbus'] = Vbus_col #sorted results RES_hightolow_Sk = RES.sort_values('Skss', ascending=False, ignore_index=True) #by SC power - high to low RES_lowtohigh_Sk = RES.sort_values('Skss', ascending=True, ignore_index=True) #by SC power - low to high RES_hightolow_Ik = RES.sort_values('Ikss', ascending=False, ignore_index=True) #by SC current - high to low RES_lowtohigh_Ik = RES.sort_values('Ikss', ascending=True, ignore_index=True) #by SC current - low to high
Plot the results
To visualize the results, I have prepared some bar plots that show the results of the analysis in sorted order (find the results of the New England system below the script).
There are two plots, each containing two subplots. The first one shows the buses sorted from highest to lowest short-circuit power and vice versa. The second one shows the buses sorted from highest to lowest short-circuit current and vice versa.
#%% plot: buses with highest and lowest SC power #plot specifications left = 0.01 # the left side of the subplots of the figure right = 0.99 # the right side of the subplots of the figure bottom = 0.01 # the bottom of the subplots of the figure top = 0.99 # the top of the subplots of the figure wspace = 0.2 # the amount of width reserved for space between subplots, # expressed as a fraction of the average axis width hspace = 0.7 # the amount of height reserved for space between subplots, # expressed as a fraction of the average axis height numhigh_Sk = 10 if len(RES_hightolow_Sk) < numhigh_Sk: numhigh_Sk = len(RES_hightolow_Sk) numlow_Sk = 10 if len(RES_lowtohigh_Sk) < numlow_Sk: numlow_Sk = len(RES_lowtohigh_Sk) #plot highest SC powers plt.figure() plt.subplots_adjust(left=left, bottom=bottom, right=right, top=top, wspace=wspace, hspace=hspace) plt.subplot(2,1,1) plt.bar(range(numhigh_Sk), RES_hightolow_Sk['Skss'][0:numhigh_Sk] , tick_label=RES_hightolow_Sk['name'][0:numhigh_Sk], color='darkred') plt.xticks(range(numhigh_Sk), rotation='vertical') plt.ylabel('Skss (MVA)') plt.grid(linestyle='--') plt.title('Buses with highest short-circuit power') #plot lowest SC powers plt.subplot(2,1,2) plt.bar(range(numlow_Sk), RES_lowtohigh_Sk['Skss'][0:numlow_Sk] , tick_label=RES_lowtohigh_Sk['name'][0:numlow_Sk]) plt.xticks(range(numlow_Sk), rotation='vertical') plt.ylabel('Skss (MVA)') plt.grid(linestyle='--') plt.title('Buses with lowest short-circuit power') #%% plot: buses with highest and lowest SC currents numhigh_Ik = 10 if len(RES_hightolow_Ik) < numhigh_Ik: numhigh_Ik = len(RES_hightolow_Ik) numlow_Ik = 10 if len(RES_lowtohigh_Ik) < numlow_Ik: numlow_Ik = len(RES_lowtohigh_Sk) #plot highest SC currents plt.figure() plt.subplots_adjust(left=left, bottom=bottom, right=right, top=top, wspace=wspace, hspace=hspace) plt.subplot(2,1,1) plt.bar(range(numhigh_Ik), RES_hightolow_Ik['Ikss'][0:numhigh_Ik] , tick_label=RES_hightolow_Ik['name'][0:numhigh_Ik], color='red') plt.xticks(range(numhigh_Ik), rotation='vertical') plt.ylabel('Ikss (kA)') plt.grid(linestyle='--') plt.title('Highest SC currents') #plot lowest SC currents plt.subplot(2,1,2) plt.bar(range(numlow_Ik), RES_lowtohigh_Ik['Ikss'][0:numlow_Ik] , tick_label=RES_lowtohigh_Ik['name'][0:numlow_Ik], color='blue') plt.xticks(range(numlow_Ik), rotation='vertical') plt.ylabel('Ikss (kA)') plt.grid(linestyle='--') plt.title('Lowest SC currents')
Below you can find the results of the anaysis of the New England system.
The first plot shows the buses with the highest and lowest short-circuit power. If you want to show more/fewer buses, you can change the parameter ‘numhigh_Sk’ or ‘numlow_Sk’ in the code above.
Clearly, Bus 39 has the highest short-circuit power within the network. This makes sense because a large generator (10000 MVA) is connected to Bus 39. The generator represents the external network – the rest of the USA and Canada.
The second plot shows the buses with the highest and lowest short-circuit current. If you want to show more/fewer buses, you can change the parameter ‘numhigh_Ik’ or ‘numlow_Ik’ in the code above.
Here, we see a different distribution. Why? Because the voltage level of the bus plays a significant role. Bus 39 is connected at a 345-kV-bus and the short-circuit power and current are related by Skss = √3 * Vbus * Ikss, where Vbus…line-to-line bus voltage. On the other hand, the impedance of the grid seen from the bus highly influences the short-circuit current, i.e. the smaller the impedance, the higher the short-circuit current and vice versa. That is why it is important to look at both, the short-circuit power and current.
If you find this post insightful, you might also be interested in the following: Automated bottleneck analysis with Python in PowerFactory
More insights within power systems are to follow in the next article.
Anytype of video tutorials/courses for python and DIgSILENT?
Hi Muzamil, unfortunately I do not have video tutorials, but I can recommend you this to start with: https://www.academia.edu/35045732/DIgSILENT_PowerFactory_Application_Guide_Python_Tutorial_DIgSILENT_Technical_Documentation