Object-Oriented Programming and DC Power Flow#
Information |
Details |
|---|---|
Lead Author |
Hantao Cui |
Learning Objectives |
• Understand OOP basics through power system components |
Prerequisites |
• Python and NumPy basics |
Estimated Time |
90-150 minutes |
Topics |
• OOP in power systems |
Required Packages |
• NumPy |
Development Status |
Completed |
Last Updated on |
April 24, 2025 |
OOP Introduction in Power Systems Context#
Object-oriented programming (OOP) is a programming paradigm that uses “objects” – instances of classes – to structure code. In power systems analysis, OOP offers significant advantages by allowing us to model physical components (buses, generators, lines) as software objects that correspond naturally to real-world entities.
Power systems are inherently complex networks of interconnected components with specific behaviors. OOP helps manage this complexity by:
Encapsulating data and functionality together in classes
Abstracting away details so we can focus on the important aspects
Enabling inheritance to share common functionality
Providing polymorphism for specialized implementations of common interfaces
For instance, a transmission line in a power grid has properties (resistance, reactance) and behaviors (power flow). Using OOP, we can create a Line class that contains both these properties and methods that calculate relevant values.
Practical Considerations
When modeling power systems with OOP, try to match your class structures to physical components. This makes your code more intuitive and easier to maintain.
Designing Power System Component Classes#
Let’s create a set of classes to model a power system. We’ll start with a BaseModel abstract class that our other classes will inherit from:
BaseModel Class#
import numpy as np
from abc import ABC, abstractmethod
class BaseModel(ABC):
"""
Base class for all power system model components.
Provides common functionality and enforces interface requirements.
"""
def __init__(self):
# Basic tracking of object properties
self.addresses = {}
self.values = {}
@abstractmethod
def get_count(self):
"""Return the number of components in the model."""
raise NotImplementedError
def str(self):
"""Return a string representation of the model component."""
return f"{self.__class__.__name__} object"
def repr(self):
"""Return a detailed string representation of the model component."""
return f"{self.__class__.__name__}"
The BaseModel class serves as a foundation for all power system components in our code. Think of it as a blueprint that defines what every power system component (like buses, generators, or transmission lines) should be able to do. It uses Python’s Abstract Base Class (ABC) feature, which is like creating a template that other classes must follow.
Let’s break it down with a real-world analogy: Imagine you’re designing different types of vehicles. Every vehicle (car, truck, motorcycle) needs certain basic features - they all need to know how many wheels they have (like our get_count() method), and they all need to be describable (like our str() and repr() methods). The @abstractmethod on get_count() is like saying “every vehicle MUST have a way to count its wheels, but each type will count them differently.”
The addresses and values dictionaries in __init__() are like a vehicle’s registration book. They keep track of important information about the component. When we create specific components (like a Bus or Generator class), they will inherit from this BaseModel, just like how a Car class inherits from a Vehicle class - they get all these basic features and must implement the required methods.
Required abstract methods
When you create a new class that inherits from BaseModel, you must implement the get_count() method, or Python will give you an error. This helps ensure that all your power system components behave consistently.
Now, let’s implement the Bus class to store bus data:
Bus Class#
import numpy as np
class Bus(BaseModel):
"""
Bus class to store bus data.
Attributes
----------
bus_i : int
Bus number
bus_type : int
Bus type (1=PQ, 2=PV, 3=slack, 4=isolated)
Pd : float
Active power demand (MW)
Qd : float
Reactive power demand (MVAr)
Gs : float
Shunt conductance (MW demanded at V = 1.0 p.u.)
Bs : float
Shunt susceptance (MVAr injected at V = 1.0 p.u.)
area : int
Area number
Vm : float
Voltage magnitude (p.u.)
Va : float
Voltage angle (degrees)
baseKV : float
Base voltage (kV)
zone : int
Loss zone
Vmax : float
Maximum voltage magnitude (p.u.)
Vmin : float
Minimum voltage magnitude (p.u.)
"""
def __init__(self, bus_i, type, Pd, Qd, Gs, Bs,
area, Vm, Va, baseKV, zone, Vmax, Vmin):
"""Initialize the Bus class."""
super().__init__()
self.bus_i = bus_i
self.bus_type = np.array(type)
self.Pd = np.array(Pd)
self.Qd = np.array(Qd)
self.Gs = np.array(Gs)
self.Bs = np.array(Bs)
self.area = area
self.Vm = np.array(Vm)
self.Va = np.array(Va)
self.baseKV = np.array(baseKV)
self.zone = zone
self.Vmax = np.array(Vmax)
self.Vmin = np.array(Vmin)
def get_count(self):
"""Return the number of buses."""
if isinstance(self.bus_i, list):
return len(self.bus_i)
else:
return 1
def is_slack(self):
"""Check if the bus is a slack bus (type 3)."""
return self.bus_type == 3
def is_pv(self):
"""Check if the bus is a PV bus (type 2)."""
return self.bus_type == 2
def is_pq(self):
"""Check if the bus is a PQ bus (type 1)."""
return self.bus_type == 1
The Bus class represents an electrical bus (node) in a power system, inheriting from our BaseModel blueprint. Each bus has specific characteristics like its number (bus_i), type (PQ, PV, or slack), power demands (Pd, Qd), and voltage parameters (Vm, Va). When we create a new bus, the __init__ method sets up all these properties, using NumPy arrays for values that need mathematical operations. The super().__init__() call is crucial.It’s like telling the bus “before you set up your own properties, make sure you have all the basic features from BaseModel.” Think of it as building a specialized electric vehicle: first you need the basic vehicle frame (called by super()), then you add the electric-specific components.
The class also provides helper methods to identify the bus type (is_slack(), is_pv(), is_pq()), which is essential for power flow calculations. These methods make the code more readable and maintainable: instead of checking bus_type == 3 everywhere in your code, you can use the more descriptive is_slack(). The get_count() method, which we inherited from BaseModel, is implemented to handle both single buses and lists of buses, making the class flexible for different use cases.
Understanding super()
The super().__init__() call is like building a house: you need the foundation (BaseModel) before adding the rooms and furniture (Bus-specific features). Without this call, you’d lose all the basic functionality that BaseModel provides. Always call super().__init__() when inheriting from another class!
Power System Bus Types
PQ Bus (Type 1): Known real and reactive power
PV Bus (Type 2): Known real power and voltage magnitude (typically generators)
Slack Bus (Type 3): Reference bus with known voltage magnitude and angle
Isolated Bus (Type 4): Not connected to the system
These are the definitions from the MATPOWER format, which may have been inspired by earlier data formats.
Branch Class#
import numpy as np
class Branch(BaseModel):
"""
Branch class to store branch data (transmission lines and transformers).
Attributes
----------
fbus : int
"From" bus number
tbus : int
"To" bus number
r : float
Resistance (p.u.)
x : float
Reactance (p.u.)
b : float
Total line charging susceptance (p.u.)
rateA : float
MVA rating A (long term rating)
rateB : float
MVA rating B (short term rating)
rateC : float
MVA rating C (emergency rating)
ratio : float
Transformer off nominal turns ratio
angle : float
Transformer phase shift angle (degrees)
status : int
Initial branch status, 1 = in-service, 0 = out-of-service
angmin : float
Minimum angle difference (degrees)
angmax : float
Maximum angle difference (degrees)
"""
def __init__(self, fbus, tbus, r, x, b, rateA, rateB, rateC, ratio, angle, status, angmin, angmax):
"""Initialize the Branch class."""
super().__init__()
self.fbus = fbus
self.tbus = tbus
self.r = np.array(r)
self.x = np.array(x)
self.b = np.array(b)
self.rateA = rateA
self.rateB = rateB
self.rateC = rateC
self.ratio = ratio
self.angle = angle
self.status = status
self.angmin = angmin
self.angmax = angmax
def get_count(self):
"""Return the number of branches."""
return len(self.fbus) if isinstance(self.fbus, list) else 1
Generator Class#
import numpy as np
class Generator(BaseModel):
"""
Generator class to store generator data.
Attributes
----------
bus : np.ndarray
Bus number
Pg : np.ndarray
Active power output (MW)
Qg : np.ndarray
Reactive power output (MVAr)
Qmax : np.ndarray
Maximum reactive power output (MVAr)
Qmin : np.ndarray
Minimum reactive power output (MVAr)
Vg : np.ndarray
Voltage magnitude setpoint (p.u.)
mBase : np.ndarray
Total MVA base of machine (MVA)
status : np.ndarray
Machine status, 1 = in-service, 0 = out-of-service
Pmax : np.ndarray
Maximum active power output (MW)
Pmin : np.ndarray
Minimum active power output (MW)
"""
def __init__(self, bus, Pg, Qg, Qmax, Qmin, Vg, mBase, status, Pmax, Pmin,
Pc1=0, Pc2=0, Qc1min=0, Qc1max=0, Qc2min=0, Qc2max=0,
ramp_agc=0, ramp_10=0, ramp_30=0, ramp_q=0, apf=0):
"""Initialize the Generator class with all fields as NumPy arrays."""
super().__init__()
self.bus = np.array(bus)
self.Pg = np.array(Pg)
self.Qg = np.array(Qg)
self.Qmax = np.array(Qmax)
self.Qmin = np.array(Qmin)
self.Vg = np.array(Vg)
self.mBase = np.array(mBase)
self.status = np.array(status)
self.Pmax = np.array(Pmax)
self.Pmin = np.array(Pmin)
# Optional cost and ramp rate parameters (scalar or vector)
self.Pc1 = np.array(Pc1)
self.Pc2 = np.array(Pc2)
self.Qc1min = np.array(Qc1min)
self.Qc1max = np.array(Qc1max)
self.Qc2min = np.array(Qc2min)
self.Qc2max = np.array(Qc2max)
self.ramp_agc = np.array(ramp_agc)
self.ramp_10 = np.array(ramp_10)
self.ramp_30 = np.array(ramp_30)
self.ramp_q = np.array(ramp_q)
self.apf = np.array(apf)
def get_count(self):
"""Return the number of generators."""
return len(self.bus)
def is_active(self):
"""
Return a boolean array indicating which generators are active.
Returns
-------
np.ndarray
Boolean array of active status flags.
"""
return self.status == 1
Format
The classes above follow the MATPOEWR convention for power system data. MATPOWER is a widely used open-source package for power system analysis in MATLAB.
Now we’ll create a PowerSystem class to contain these components and provide system-level functionality:
PowerSystem Class#
import numpy as np
class PowerSystem:
"""
PowerSystem class to store and manage power system data.
This class contains buses, branches (lines/transformers), and generators.
It provides methods for system-wide operations like power flow analysis.
"""
def __init__(self, bus=None, branch=None, gen=None):
"""Initialize the PowerSystem class."""
self.bus = bus
self.branch = branch
self.gen = gen
self.baseMVA = 100.0 # Default base MVA
def get_count(self):
"""
Return the number of components in the system.
Returns
-------
dict
Dictionary with counts of buses, branches, and generators.
"""
return {
'buses': self.bus.get_count() if self.bus else 0,
'branches': self.branch.get_count() if self.branch else 0,
'generators': self.gen.get_count() if self.gen else 0
}
def get_adjacency_matrix(self):
"""
Create an adjacency matrix for the network topology.
Returns
-------
np.ndarray
Adjacency matrix representing connectivity between buses.
"""
n_bus = self.bus.get_count()
adj_matrix = np.zeros((n_bus, n_bus))
for i in range(self.branch.get_count()):
from_bus = self.branch.fbus[i] - 1 # Convert to 0-indexing
to_bus = self.branch.tbus[i] - 1
if self.branch.status[i] == 1: # Only in-service branches
adj_matrix[from_bus, to_bus] = 1
adj_matrix[to_bus, from_bus] = 1 # Undirected graph
return adj_matrix
def __str__(self):
"""
Return a string representation of the power system.
Returns
-------
str
Summary of the power system components.
"""
counts = self.get_count()
return (
f"PowerSystem with {counts['buses']} buses, "
f"{counts['branches']} branches, and "
f"{counts['generators']} generators"
)
The PowerSystem class serves as the central orchestrator in our power system analysis framework, bringing together all the individual components (buses, branches, and generators) into a cohesive system. Think of it as a conductor of an orchestra - while each musician (component) knows how to play their instrument, the conductor ensures they work together harmoniously. The class maintains these components and provides a common interface to access and manipulate them. The baseMVA attribute sets the system’s base power, which is crucial for per-unit calculations, defaulting to the industry-standard value of 100 MVA.
Beyond just storing components, the PowerSystem class provides system-wide analysis capabilities. Its get_adjacency_matrix() method creates a mathematical representation of how buses are connected through branches, which is fundamental for network analysis. This matrix representation is like a map showing which buses are connected, with 1’s indicating connections and 0’s indicating no direct connection.
The class’s __str__() method provides a human-readable summary of the system’s size, making it easy to quickly understand the scale of the power system being analyzed or other general information. This orchestration role makes PowerSystem the primary interface through which most power system analyses (like power flow studies) will be conducted.
Adjacency Matrix
The adjacency matrix is a fundamental tool in network analysis:
# For a 3-bus system with branches between buses 1-2 and 2-3:
[[0, 1, 0], # Bus 1's connections
[1, 0, 1], # Bus 2's connections
[0, 1, 0]] # Bus 3's connections
Each row and column represents a bus, and 1’s indicate direct connections.
Data Parsing from XLSX (MATPOWER Format)#
To use our classes, we need to load data from Excel files. Let’s implement functions to read data from Excel and convert it to our class objects:
Data Parsing Functions#
import pandas as pd
def read_from_excel(path):
"""
Read power system data from Excel file.
Parameters
----------
path : str
Path to the Excel file
Returns
-------
dict
Dictionary containing pandas DataFrames for bus, gen, and branch data
"""
# Read all sheets from the Excel file
data_dict = pd.read_excel(path, sheet_name=None)
# Verify that the required sheets exist
required_sheets = ['bus', 'gen', 'branch']
for sheet in required_sheets:
if sheet not in data_dict:
raise ValueError(f"Required sheet '{sheet}' not found in Excel file")
return data_dict
def column_to_obj(df, cls):
"""
Convert a DataFrame to an object.
Parameters
----------
df : pd.DataFrame
DataFrame to be converted
cls : class
Class to be instantiated
Returns
-------
object
Object of the specified class
"""
# Create a dictionary of column names and lists of values
obj_dict = {col: df[col].tolist() for col in df.columns}
return cls(**obj_dict)
# Map sheet names to their corresponding classes
name_to_class = {
'bus': Bus,
'branch': Branch,
'gen': Generator
}
def create_system_from_data(data):
"""
Create a PowerSystem object from a dictionary of data.
Parameters
----------
data : dict
Dictionary of DataFrames
Returns
-------
PowerSystem
Instantiated PowerSystem object
"""
model_objects = {key: column_to_obj(data[key], cls) for key, cls in name_to_class.items()}
return PowerSystem(**model_objects)
def load_data(path):
"""
Load data from an Excel file and create a PowerSystem object.
Parameters
----------
path : str
Path to the Excel file
Returns
-------
PowerSystem
PowerSystem object
"""
data = read_from_excel(path)
return create_system_from_data(data)
Let’s see how this works with a simple example:
# Load a power system from an Excel file
system = load_data('case3.xlsx')
print(system)
# Access system components
print(f"Voltage at bus 1: {system.bus.Vm[0]:.4f} p.u.")
print(f"First branch connects bus {int(system.branch.fbus[0])} to bus {int(system.branch.tbus[0])}")
PowerSystem with 3 buses, 3 branches, and 1 generators
Voltage at bus 1: 1.0200 p.u.
First branch connects bus 1 to bus 2
Building the Admittance Matrix#
The admittance matrix (Y-bus) is fundamental in power system analysis. It represents the network’s topology and electrical characteristics. Let’s implement a method to build this matrix in our PowerSystem class:
Y-Bus Matrix Construction#
def build_ybus(self):
"""
Build the bus admittance matrix (Y-bus) for the power system.
Returns
-------
np.ndarray
Complex bus admittance matrix (Y-bus)
"""
n_bus = self.bus.get_count()
y_bus = np.zeros((n_bus, n_bus), dtype=complex)
# Add branch admittances
for i in range(self.branch.get_count()):
if self.branch.status[i] == 0: # Skip out-of-service branches
continue
# Get bus indices (convert from 1-indexing to 0-indexing)
f_bus = int(self.branch.fbus[i]) - 1
t_bus = int(self.branch.tbus[i]) - 1
# Branch parameters
r = self.branch.r[i]
x = self.branch.x[i]
b = self.branch.b[i]
# Series admittance
y_series = 1.0 / complex(r, x)
# Shunt admittance (half at each end)
y_shunt = complex(0, b / 2)
# Off-diagonal elements (mutual admittances)
y_bus[f_bus, t_bus] -= y_series
y_bus[t_bus, f_bus] -= y_series
# Diagonal elements (self admittances)
y_bus[f_bus, f_bus] += y_series + y_shunt
y_bus[t_bus, t_bus] += y_series + y_shunt
# Add bus shunt admittances (from bus data)
for i in range(n_bus):
gs = self.bus.Gs[i]
bs = self.bus.Bs[i]
y_bus[i, i] += complex(gs, bs) / self.baseMVA
return y_bus
Add this method to the PowerSystem class.
PowerSystem.build_ybus = build_ybus
The line above assigns the method to the PowerSystem class. This is allowed in Python.
Unit
The Y-bus matrix construction assumes that all values are in per-unit. Make sure your input data is properly converted before using this method.
Let’s look at the Y-bus matrix in more detail:
The diagonal elements of Y-bus (Y_ii) represent the sum of all admittances connected to bus i.
The off-diagonal elements (Y_ij) represent the negative of the admittance between buses i and j.
The Y-bus matrix is used in both power flow and short circuit analysis. For a system with n buses, the Y-bus matrix is an n×n complex-valued matrix.
More Notes on Y-bus Matrix
For large power systems, the Y-bus matrix is typically sparse (most entries are zero), since most buses are connected to only a few other buses. Special sparse matrix data structures should be used for efficient computation.
This admittance matrix does not consider transformers.
Assumption on bus index
This code also makes an assumption about the bus indices. It assumes the indices start at 1 and are contiguous. Otherwise, the code will fail.
A general implementation will convert the user-named bus index to a contiguous 0-based internal index. To focus on conceptual discussion, we will skip that practical part.
DC Power Flow Analysis#
Now, let’s implement DC power flow analysis, which is a linearized approximation of the AC power flow. It’s useful for quick analysis and is computationally efficient. DC power flow makes the following assumptions:
Flat voltage profile (all voltages are 1.0 p.u.)
Line resistances are negligible (r << x)
Voltage angle differences are small (sin(θ) ≈ θ)
The resulting DC power flow equation becomes:
where:
\(P\) is the vector of real power injections
\(B\) is the susceptance matrix
\(\theta\) is the vector of voltage angles
Solving for Voltage Angles#
To solve for voltage angles, we need to:
Remove the slack bus row and column from B (creating \(B'\))
Remove the slack bus element from P (creating \(P'\))
Solve the linear system:
Branch Power Flows#
Once voltage angles are known, the power flow on a branch between buses \(i\) and \(j\) can be calculated as:
This linear relationship between power flows and angle differences is what makes DC power flow analysis computationally efficient and suitable for many planning and operational studies.
Implementation#
The function implementation is broken down into several parts, each led by an explanation. To follow this training, you are encouraged to type up the code in your local Jupyter Notebook.
def run_dc_power_flow(self):
"""
Perform DC power flow analysis.
DC power flow is a linearized approximation of the AC power flow equations,
based on the following assumptions:
- Flat voltage profile (all voltages = 1.0 p.u.)
- Line resistances are negligible (r << x)
- Voltage angle differences are small (sin(θ) ≈ θ)
Returns
-------
dict
Dictionary containing results of the DC power flow:
- voltage_angles : np.ndarray
- branch_flows : np.ndarray
- slack_bus_injection : float
"""
n_bus = self.bus.get_count()
In DC power flow analysis, we derive a simplified linear relationship between power injections and voltage angles. The B matrix, also known as the susceptance matrix, plays a central role in this formulation.
From Y Matrix to B Matrix#
The complex admittance matrix (Y matrix) for a power system is given by:
where:
\(G_{ij}\) is the conductance between buses \(i\) and \(j\)
\(B_{ij}\) is the susceptance between buses \(i\) and \(j\)
Under the three assumptions, the B matrix is constructed from the imaginary part of the Y matrix:
where \(x_{ij}\) is the reactance of the branch between buses \(i\) and \(j\).
# Step 1: Build the B matrix (imaginary part of Y-bus)
y_bus = self.build_ybus()
b_matrix = np.imag(y_bus)
# Step 2: Identify the slack bus and remove it from the matrices
slack_bus_idx = None
for i in range(n_bus):
if self.bus.bus_type[i] == 3: # Slack bus
slack_bus_idx = i
break
if slack_bus_idx is None:
raise ValueError("No slack bus found in the system")
# Remove slack bus row and column from B matrix
b_reduced = np.delete(b_matrix, slack_bus_idx, axis=0)
b_reduced = np.delete(b_reduced, slack_bus_idx, axis=1)
Next, we will form the power injections. Note that the slack needs to be removed, because the full system is rank deficient (all power injections add up to zero).
# Step 3: Calculate net power injections (P = Pg - Pd)
p_injection = np.zeros(n_bus)
# Add generation
for i in range(self.gen.get_count()):
if self.gen.status[i] == 1: # In-service generators only
bus_idx = int(self.gen.bus[i]) - 1 # Convert to 0-indexing
p_injection[bus_idx] += self.gen.Pg[i]
# Subtract load
for i in range(n_bus):
p_injection[i] -= self.bus.Pd[i]
# Convert to per-unit
p_injection /= self.baseMVA
# Remove slack bus injection
p_reduced = np.delete(p_injection, slack_bus_idx)
With the injections calculated, we can solve the linear equation for the phase angles.
# Step 4: Solve for voltage angles (Bθ = P)
try:
theta_reduced = np.linalg.solve(b_reduced, p_reduced)
except np.linalg.LinAlgError as e:
raise ValueError(f"Failed to solve DC power flow: {str(e)}")
# Step 5: Re-insert the slack bus angle (0)
theta = np.zeros(n_bus)
j = 0
for i in range(n_bus):
if i != slack_bus_idx:
theta[i] = theta_reduced[j]
j += 1
Then, we can use the DC power flow equation \(P_{ij} = (\theta_i - \theta_j)/X\) to calculate the branch flows:
# Step 6: Calculate branch flows
branch_flows = np.zeros(self.branch.get_count())
for i in range(self.branch.get_count()):
if self.branch.status[i] == 1:
f_bus = int(self.branch.fbus[i]) - 1
t_bus = int(self.branch.tbus[i]) - 1
x = self.branch.x[i]
branch_flows[i] = (theta[f_bus] - theta[t_bus]) / x
# Convert branch flows back to MW
branch_flows *= self.baseMVA
There are multiple ways to calculate the slack injections. The following shows a branch-based approach. One can also calculate the total power injection using the \(B\) matrix without the slack bus to obtain the missing active power, and then flip the sign to get the slack injection.
# Step 7: Compute slack bus injection
slack_injection = 0
for i in range(self.branch.get_count()):
if self.branch.status[i] == 1:
f_bus = int(self.branch.fbus[i]) - 1
t_bus = int(self.branch.tbus[i]) - 1
if f_bus == slack_bus_idx:
slack_injection += branch_flows[i]
elif t_bus == slack_bus_idx:
slack_injection -= branch_flows[i]
# Return results
results = {
'voltage_angles': theta,
'branch_flows': branch_flows,
'slack_bus_injection': slack_injection
}
return results
PowerSystem.run_dc_power_flow = run_dc_power_flow
system = load_data('case3.xlsx')
dc_results = system.run_dc_power_flow()
Let’s use the DC power flow method on our example system:
# Run DC power flow analysis
system = load_data('case3.xlsx')
dc_results = system.run_dc_power_flow()
# Print results
print("\nVoltage Angles (radians):")
for i, angle in enumerate(dc_results['voltage_angles']):
print(f"Bus {i+1}: {angle:.6f}")
print("\nBranch Flows (MW):")
for i, flow in enumerate(dc_results['branch_flows']):
from_bus = int(system.branch.fbus[i])
to_bus = int(system.branch.tbus[i])
print(f"Branch {i+1} (Bus {from_bus} to Bus {to_bus}): {flow:.2f} MW")
print(f"\nSlack Bus Injection: {dc_results['slack_bus_injection']:.2f} MW")
Voltage Angles (radians):
Bus 1: 0.000000
Bus 2: 0.012692
Bus 3: 0.067286
Branch Flows (MW):
Branch 1 (Bus 1 to Bus 2): -4.23 MW
Branch 2 (Bus 2 to Bus 3): -54.59 MW
Branch 3 (Bus 3 to Bus 1): 67.29 MW
Slack Bus Injection: -71.52 MW
Notes
DC power flow is a great tool for: Initial screening of potential system violations Transmission planning studies Generator dispatch optimization Contingency analysis where speed is more important than accuracy However, it doesn’t consider reactive power flows or voltage magnitudes, which are important for a complete system analysis.
Advanced OOP Concepts in Practice#
Now that we’ve implemented the basic classes and functionality, let’s explore some advanced OOP concepts that can enhance our power system analysis toolkit.
Inheritance and Polymorphism#
We’ve already seen inheritance in action with our classes inheriting from BaseModel. We can extend this concept by creating specialized versions of our analysis tools:
Advanced OOP Example - Contingency Analysis#
import copy
class ContingencyAnalysis:
"""
Class to perform contingency analysis on power systems.
Contingency analysis simulates the outage of system components
(lines, generators, etc.) to assess system reliability.
"""
def __init__(self, power_system):
"""Initialize with a power system."""
self.original_system = power_system
self.results = {}
def run_n_1_analysis(self, contingency_type='branch'):
"""
Perform N-1 contingency analysis.
Parameters
----------
contingency_type : str
Type of contingency to analyze ('branch' or 'gen')
Returns
-------
dict
Dictionary containing results of each contingency
"""
self.results = {}
if contingency_type == 'branch':
# Analyze branch (line) outages
for i in range(self.original_system.branch.get_count()):
if self.original_system.branch.status[i] == 0:
continue # Skip already out-of-service branches
system_copy = copy.deepcopy(self.original_system)
system_copy.branch.status[i] = 0
try:
results = system_copy.run_dc_power_flow()
overloads = []
for j, flow in enumerate(results['branch_flows']):
if system_copy.branch.status[j] == 1:
limit = system_copy.branch.rateA[j]
if limit > 0 and abs(flow) > limit:
overloads.append((j, flow, limit))
from_bus = int(self.original_system.branch.fbus[i])
to_bus = int(self.original_system.branch.tbus[i])
key = f"Branch_{i+1}_{from_bus}_{to_bus}"
self.results[key] = {
'flow_results': results,
'overloads': overloads,
'is_converged': True
}
except Exception as e:
from_bus = int(self.original_system.branch.fbus[i])
to_bus = int(self.original_system.branch.tbus[i])
key = f"Branch_{i+1}_{from_bus}_{to_bus}"
self.results[key] = {
'error': str(e),
'is_converged': False
}
elif contingency_type == 'gen':
# Analyze generator outages
for i in range(self.original_system.gen.get_count()):
if self.original_system.gen.status[i] == 0:
continue # Skip already out-of-service generators
system_copy = copy.deepcopy(self.original_system)
system_copy.gen.status[i] = 0
try:
results = system_copy.run_dc_power_flow()
overloads = []
for j, flow in enumerate(results['branch_flows']):
limit = system_copy.branch.rateA[j]
if limit > 0 and abs(flow) > limit:
overloads.append((j, flow, limit))
bus_id = int(self.original_system.gen.bus[i])
key = f"Gen_{i+1}_Bus{bus_id}"
self.results[key] = {
'flow_results': results,
'overloads': overloads,
'is_converged': True
}
except Exception as e:
bus_id = int(self.original_system.gen.bus[i])
key = f"Gen_{i+1}_Bus{bus_id}"
self.results[key] = {
'error': str(e),
'is_converged': False
}
return self.results
def get_critical_contingencies(self):
"""
Identify critical contingencies that result in overloads.
Returns
-------
list
List of tuples: (contingency name, list of overloads)
"""
if not self.results:
raise ValueError("No contingency analysis results available")
critical = []
for contingency, result in self.results.items():
if result.get('is_converged') and result.get('overloads'):
critical.append((contingency, result['overloads']))
return critical
This ContingencyAnalysis class demonstrates several advanced OOP concepts:
Composition - This class contains a
PowerSystemobject rather than inheriting from itEncapsulation - The analysis logic is contained within its own class
Polymorphism - The analysis can handle different types of contingencies (branch or generator)
“Has-a” versus “Is-a”
When designing OOP systems, consider “has-a” relationships (composition) versus “is-a” relationships (inheritance). In this case, a contingency analysis has a power system, it isn’t a power system itself.
Decorator Pattern#
Another advanced OOP concept is the decorator pattern, which allows adding new behavior to objects dynamically. Let’s implement a simple monitor decorator:
def add_monitoring(system):
"""
Decorator function that adds monitoring capabilities to a PowerSystem.
This wraps the `run_dc_power_flow` method with additional logging and
branch flow monitoring.
Parameters
----------
system : PowerSystem
The power system to enhance.
Returns
-------
PowerSystem
The enhanced power system with monitoring.
"""
import time
original_run_dc_power_flow = system.run_dc_power_flow
def monitored_dc_power_flow(*args, **kwargs):
"""
Enhanced DC power flow with timing and branch loading warnings.
Returns
-------
dict
DC power flow results with monitoring.
"""
print("Starting DC power flow computation...")
start_time = time.time()
# Call the original method
results = original_run_dc_power_flow(*args, **kwargs)
elapsed_time = time.time() - start_time
print(f"DC power flow completed in {elapsed_time:.4f} seconds")
# Add monitoring for branch loading
branch_flows = results['branch_flows']
for i, flow in enumerate(branch_flows):
if system.branch.status[i] == 1: # Only in-service branches
limit = system.branch.rateA[i]
if limit > 0: # Only branches with defined limits
loading_percent = abs(flow) / limit * 100
if loading_percent > 90:
print(f"WARNING: Branch {i + 1} loaded at {loading_percent:.1f}% of capacity")
return results
# Replace the original method with the enhanced version
system.run_dc_power_flow = monitored_dc_power_flow
return system
Let’s see the decorator in action:
# Load a power system
system = load_data('case3.xlsx')
# Add monitoring capabilities
monitored_system = add_monitoring(system)
# Run DC power flow with monitoring
results = monitored_system.run_dc_power_flow()
Starting DC power flow computation...
DC power flow completed in 0.0001 seconds
Conclusion#
In this tutorial, we’ve developed a comprehensive OOP framework for power systems analysis. We’ve covered:
Basic power system component classes (buses, branches, generators)
Data parsing from Excel files in MatPower format
Building the admittance matrix for network representation
DC power flow analysis for linear power flow solutions
Advanced OOP concepts like inheritance, composition, and decorators
This framework provides a solid foundation for power system analysis using Python and object-oriented programming. You can extend it further with AC power flow, optimal power flow, stability analysis, and other advanced power system studies.