 # Automated bottleneck (critical line) analysis with Python in DIgSILENT PowerFactory

On the 8th of January 2021, the system separation incident in Continental Europe started with the tripping of a line in Croatia due to overcurrent protection. The subsequent shift of the power flow overloaded the neighboring lines and caused them to trip as well. These cascading events are undoubtedly very critical for power systems.

This shows how important it is to analyze the bottlenecks of the power system and under which conditions they appear.

In this post, I will provide Python code to conduct a basic bottleneck analysis with (probably) any DIgSILENT PowerFactory model. The basic analysis takes the actual generator and load setpoints and increases the total load (and so generation) until a pre-defined number of lines is overloaded. The distribution factors for the increase of the load and generator setpoints are determined by their initial setpoints.

If you want to do a more advanced bottleneck analysis, click on the banner below to download the advanced script. The advanced version includes a Monte Carlo simulation where the distribution of load and generator setpoints is varied to get a more complete picture of the possible critical lines aka bottlenecks. You can find more details about the advanced analysis at the end of the post.

If you are new to Python-PowerFactory simulation check this out: Run PowerFactory through Python

## Before we continue…

For the demonstration, 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.

## Basic bottleneck analysis in PowerFactory with Python

In the following section I will show you the different steps of the analysis. They include:

1. Import of libraries and definition of functions
2. The core of the bottleneck analysis
3. Print the summary of results
4. Restore the initial setpoints of generators and loads

### Import libraries and define functions

As always, at the beginning of the script definition of the needed libraries and functions. Here we need NumPy and pandas. Additionally, I have also defined two functions to query the load and generation setpoints.

Important: don’t forget to replace the path (shown in red) with the path to the Python directory of your PowerFactory installation as indicated.

``````"""
michael@thesmartinsights.com
www.thesmartinsights.com
"""
import numpy as np
import pandas as pd

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')
#------------------------------------------------------------------------------

#get active/reactive power setpoints of all loads in service

#get active/reactive power setpoints, operation mode and min/max capacity of all generators in service
def getGenSpecs():
Generators = app.GetCalcRelevantObjects('*.ElmSym')
name, P, Q, avr_mode, Pmin, Pmax = [], [], [], [], [], []
for gen in Generators: #loop through list
if getattr(gen, 'outserv') == 0: #only SGs in service
name.append(getattr(gen, 'loc_name')) # get name of the generator
P.append(getattr(gen, 'pgini')) #active power
Q.append(getattr(gen, 'qgini')) #reactive power
avr_mode.append(getattr(gen, 'av_mode')) #AVR mode
Pmin.append(getattr(gen, 'Pmin_uc'))
Pmax.append(getattr(gen, 'Pmax_uc'))
SpecsGen = pd.DataFrame(zip(name,P,Pmin,Pmax,Q,avr_mode), columns=['name','P','Pmin','Pmax','Q','avr_mode'])
return SpecsGen``````

### The core of the bottleneck analysis

After you have replaced the path to your PowerFactory installation (see above) you onle have to specify the following points and you are ready to run the analysis.

• project name and study case

The analysis runs until the line loadings exceed the defined number of permissible overloads, e.g. maxOverloads = 1 and OL_threshold = 100 % – the analysis stops when the loading of one line is above 100 %. More information on each step in the analysis can be directly found in the comments of the code (please see below).

``````#define project name and study case
projName =   '_TSI_39_Bus_New_England_System'
study_case = '01_StudyCase.IntCase'

#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 initial setpoints of generators and loads by load flow with actual settings
setattr(oLoadflow, 'iopt_net', 0) # AC LF balanced
setattr(oLoadflow, 'iopt_apdist', 0) # P as dispatched
setattr(oLoadflow, 'iopt_plim', 1) # consider active power limits
setattr(oLoadflow, 'iopt_lim', 1) # consider reactive power limits
# setattr(oLoadflow, 'iPbalancing', 1) # slack: reference machine

#get the generators and their active power operational limits
Generators = app.GetCalcRelevantObjects('*.ElmSym')
name, P, Q, avr_mode, Pmin, Pmax = [], [], [], [], [], []
for gen in Generators: #loop through list
if getattr(gen, 'outserv') == 0: #only SGs in service
name.append(getattr(gen, 'loc_name')) # get name of the generator
# P.append(getattr(gen, 'pgini'))
P.append(getattr(gen, 'm:P:bus1'))
Q.append(getattr(gen, 'qgini')) #reactive power
avr_mode.append(getattr(gen, 'av_mode')) #AVR mode
Pmin.append(getattr(gen, 'Pmin_uc'))
Pmax.append(getattr(gen, 'Pmax_uc'))
SpecsGen_init = pd.DataFrame(zip(name,P,Pmin,Pmax,Q,avr_mode), columns=['name','P','Pmin','Pmax','Q','avr_mode'])
# SpecsGen_init = pd.DataFrame(zip(name,P,Pmin,Pmax), columns=['name','P','Pmin','Pmax'])

Lines=app.GetCalcRelevantObjects('*.ElmLne')
for line in Lines: #loop through list
if getattr(line, 'outserv') == 0: #only line is in service
LineName.append(getattr(line, 'loc_name')) # get name of the line

#total generation, generation capacity and load
TotGen_init = np.sum(SpecsGen_init['P'])
TotGenCap = np.sum(SpecsGen_init['Pmax'])

#tangent phi = P / Q (assumed constant)

#generation shift key / load shift key - based on total initial generaton and load, respectively
GSK = np.array(SpecsGen_init['P']) / TotGen_init

#load variation: starts at initial load and increases by 1 % after each step until permissible overloads are reached

i = i + 1
#set new SG setpoints according to GSK
Generators = app.GetCalcRelevantObjects('*.ElmSym')
for ig, gen in enumerate(Generators): #loop through list
if getattr(gen, 'outserv') == 0: #only SGs in service
setattr(gen, 'pgini', TotLoad_MinToMax[i] * GSK[ig]) # generator setpoints acc. to GSK

#set new load setpoints according to LSK
#save scenario
if app.GetActiveScenario() != None:
actScen = app.GetActiveScenario()
actScen.Save()

#load flow settings / now distributed slack by SGs to avoid too high loading of reference machine
setattr(oLoadflow, 'iopt_net', 0) # AC LF balanced
setattr(oLoadflow, 'iopt_apdist', 0) # P as dispatched
setattr(oLoadflow, 'iopt_plim', 1) # consider active power limits
setattr(oLoadflow, 'iopt_lim', 1) # consider reactive power limits
setattr(oLoadflow, 'iPbalancing', 4) # distributed slack by SGs

Lines=app.GetCalcRelevantObjects('*.ElmLne')
for line in Lines: #loop through list
if getattr(line, 'outserv') == 0: #only line is in service
LineName.append(getattr(line, 'loc_name')) # get name of the line

#get SG and load setpoints after calculation
SpecsGen_after = getGenSpecs()

### Print summary of results

The code below prints a summary of the results in the Python console. The summary includes the active power of the generators, active/reactive power of the loads, and the line loadings before and after the analysis. The lines are sorted according to their loading from highest to lowest. Additionally, the total load initially and at the end of the analysis is given as well as the number of overloaded lines.

``````Results_SG = pd.DataFrame()
Results_SG['name'] = SpecsGen_init['name']
Results_SG['P_SG_init'] = SpecsGen_init['P']
Results_SG['P_SG_after'] = SpecsGen_after['P']
Results_SG['Delta Generation'] = SpecsGen_after['P'] - SpecsGen_init['P']

Results_Lines = pd.DataFrame()

print('\n')
print('Summary------------------------------------------------')
print('\n')
print('Synchronous Generators----------------------')
print(Results_SG.to_string(index=False))
print('\n')
print('\n')
print('Lines---------------------------------------')
print(Results_Lines.to_string(index=False))
print('\n')

### Restore initial setpoints of generators and loads

If you want to restore the initial setpoints of the generators and loads, you can use the following code. If you don’t use it, the final setpoints at the end of the analysis will remain in the PowerFactory model.

``````for ig, gen in enumerate(Generators): #loop through list
if getattr(gen, 'outserv') == 0: #only SGs in service
setattr(gen, 'pgini', SpecsGen_init['P'][ig])

#save scenario
if app.GetActiveScenario() != None:
actScen = app.GetActiveScenario()
actScen.Save()``````

## Advanced bottleneck analysis using Monte Carlo simulations

In the following, I will provide a brief description of the capabilities of the advanced bottleneck analysis. In order to get a more complete picture of the possible bottlenecks of a power system, different load flow situations need to be analyzed. Here, I am doing that by generating load samples within a defined minimum and maximum total load. For each load sample, I generate different generator setpoints that meet the total load to achieve various load flow scenarios. That basically means different load and generation combinations.

### User-specific inputs for the advanced analysis

Basically, you can get started with the advanced analysis by specifying the PowerFactory project name, study case, and the total load range. The rest is already pre-defined as shown in the screenshot below.

The advanced bottleneck analysis has the following user-specific inputs:

• generation samples: number of generation samples analyzed for each load sample
• standard deviation of individual loads: defines how much the individual loads are varied around their mean (the initial load setpoints serve as the mean for each load individually)

### Results of the advanced analysis

Below you will find an example of the results from the advanced analysis of the New England system. The most important output of the analysis is the ranking of the overloaded lines by the number of occurrences in the load flow simulations. By default, the 10 most overloaded lines are shown.

However, the advanced analysis also allows to dig deeper into the load flow situations. For example, if you want to know how were the setpoints of a specific generator distributed when a specific line was overloaded. Below an example of the generator G 01 setpoints when Line 02 – 03 was overloaded.

I hope you find this article insightful. Consider leaving a comment to let me know your thoughts on the topic.

More insights within energy are to follow in the next article.

Btw, have you also checked how to do peak shaving with energy storage? If not, find it here: Peak shaving with energy storage

Michael