Concepts

This page explains the physical and numerical model that groundinsight implements. It is aimed at readers who want to understand why the code returns the numbers it returns — for example to validate results against measurements or analytical expressions.

Problem statement

Consider a power grid composed of substations (buses) and the cables and overhead lines (branches) that connect them. Each bus possesses a grounding grid tying it to remote earth; each branch carries a grounding conductor (cable shield, overhead-line earth wire) that is bonded to the grids of its two terminals. A single-phase-to-ground fault at one bus drives a fault current back towards the source. That current splits into two parallel paths:

• the local earth path through the grounding grid at the fault bus, and
• the metallic return path through the grounding conductor(s) of the connecting branches.

The split depends on the impedances of both paths and on the mutual coupling between the faulted phase and the grounding conductor that runs alongside it. The resulting EPR, the reduction factor and the grounding impedance are the quantities of interest.

Objects

groundinsight represents the grid with four Pydantic models:

• Bus — a node. Carries a grounding impedance \(Z_{\text{B}}(\rho_E, f)\) that couples the bus to remote earth.
• Branch — an edge between two buses with a self impedance \(Z_{\text{self}}(\rho_E, f, l)\) (series impedance of the grounding conductor) and a mutual impedance \(Z_{\text{mutual}}(\rho_E, f, l)\) (coupling between faulted phase and grounding conductor).
• Source — a current source anchored at a bus. Holds a dictionary that maps each frequency to the injected phasor current.
• Fault — a marker at a bus. Holds a dictionary of frequency-dependent scaling factors that reduce the source current for each harmonic.

Bus and branch properties are derived from reusable types (BusType, BranchType) whose impedances are parameterised by formula strings. The formulas are parsed with SymPy, turned into lambdify callables and evaluated for every network frequency.

Impedance formulas

Formula strings may contain the symbols

Symbol	Meaning
rho	specific earth resistance \(\rho_E\) in \(\Omega\,\text{m}\)
f	frequency \(f\) in Hz
l	branch length \(l\) in km
j	imaginary unit (internally 1j)

Expressions are evaluated symbolically, so any analytic expression supported by SymPy is admissible. The literal string "nan" maps to an effectively infinite impedance and can be used to model open ends (broken shields, isolators etc.).

Nodal-admittance formulation

All computations take place per frequency \(f\) in the phasor domain. groundinsight assembles an admittance matrix \(Y(f)\) of size \(N\times N\) (where \(N\) is the number of buses):

\[ Y_{ii} = \frac{1}{Z_{\text{B},i}(\rho_E, f)} + \sum_{k \in \mathcal{E}(i)} \frac{1}{Z_{\text{self},k}(\rho_E, f, l_k)}, \qquad Y_{ij} = -\frac{1}{Z_{\text{self},k}(\rho_E, f, l_k)} \quad (k \text{ connects } i,j). \]

Each entry can be switched off by setting grounding_conductor=False on the corresponding branch, which is useful for modelling insulated shields.

The right-hand-side vector \(\underline{i}(f)\) holds the source injections scaled by the active fault's frequency scaling and — crucially — the mutual-coupling contributions. For every branch on a path from source to fault a Norton equivalent is added: the phase current \(I_{\text{phase}}\) driving the branch induces a shield current of magnitude \(I_{\text{mut}} = I_{\text{phase}}\,Z_{\text{mutual}}/Z_{\text{self}}\), injected as \(-I_{\text{mut}}\) at the from bus and \(+I_{\text{mut}}\) at the to bus (signs follow the direction source → fault).

The EPR vector is then

\[ \underline{u}(f) = Y(f)^{-1}\,\underline{i}(f). \]

Numerically the system is solved with SciPy's sparse LU decomposition (scipy.sparse.csc_matrix + splu) — this scales well to meshed low-voltage networks with thousands of buses.

Path finding

Mutual-coupling injections require a direction. groundinsight derives that direction by enumerating every simple path from each source bus to the active fault bus via a depth-first search (PathFinder). Each path is stored as an ordered list of Branch objects; its injection signs follow the traversal order.

In ring or meshed topologies a single source–fault pair yields multiple paths. By default every path carries the full source current. The optional parallel_coefficient on a branch lets you pre-scale the current share of individual parallel legs; if you set auto_parallel_coefficients=True on run_fault, groundinsight solves a reduced phase-only network first and uses its current distribution as the per-path scaling.

Derived quantities

Once \(\underline{u}(f)\) is known for every frequency, three result families are computed:

Earth potential rise (EPR)

The per-frequency bus voltages are stored as ComplexNumber entries on ResultBus objects. RMS values across all frequencies are computed via

\[ U_{\text{RMS},i} = \sqrt{\sum_{f} |u_i(f)|^2}. \]

Branch currents

For every branch and frequency the current through the grounding conductor is

\[ I_{\text{branch}}(f) = \frac{u_{\text{from}}(f) - u_{\text{to}}(f)} {Z_{\text{self}}(f)} + I_{\text{mut}}(f), \]

where the second term accounts for the Norton source representing the mutual coupling (compute_branch_currents).

Reduction factor

The reduction factor \(r\) at the fault bus is defined as

\[ r(f) = \frac{|u_{\text{fault}}^{\text{(with mutual)}}(f)|} {|u_{\text{fault}}^{\text{(without mutual)}}(f)|}. \]

groundinsight obtains the denominator by re-solving the same network with all mutual-coupling Norton sources removed. For a single shielded line directly between source and fault with identical impedances the expression collapses to the familiar analytical form

\[ r = \left| 1 - \frac{Z_{\text{mutual}}}{Z_{\text{self}}} \right|. \]

The same closed form applies to a fully symmetric ring with the fault diametrically opposite the source — the Norton injections in the two ring halves are then perfectly anti-parallel and superpose to the same expression as the single-line case.

Frequency dependence of the reduction factor

For a shielded cable the impedances split into a resistive and an inductive part,

\[ Z_{\text{self}}(f) = R + j\,\omega L, \qquad Z_{\text{mutual}}(f) = j\,\omega M, \qquad \omega = 2\pi f. \]

With full coupling between the faulted phase and the grounding conductor (the geometric ideal \(M = L\)) the closed form simplifies to

\[ r(f) \;=\; \left| 1 - \frac{j\,\omega L}{R + j\,\omega L} \right| \;=\; \frac{R}{\sqrt{R^{2} + (\omega L)^{2}}}. \]

Two limits follow directly:

• \(f \to 0\): \(\omega L \to 0\), so \(Z_{\text{mutual}}/Z_{\text{self}} \to 0\) and \(r \to 1\). At DC the shield carries no induced current, the entire fault current flows through the local earth path and the reduction factor equals 1.
• \(f \to \infty\): the imaginary parts dominate, so \(Z_{\text{mutual}}/Z_{\text{self}} \to M/L = 1\) and \(r \to 0\). At high frequency the shield short-circuits the inductive coupling, almost the entire fault current returns metallically and the EPR collapses.

In a real cable the resistive part is small compared with \(\omega L\) already at power frequency, which is why MV cables typically reach \(r \approx 0.3 \ldots 0.4\) at 50 Hz and the reduction factor decays quickly above one or two hundred Hz. This convergence is exercised explicitly by tests/test_topology_and_reduction.py::test_reduction_factor_sweep_*, which sweeps a single MV cable section and a 20-bus symmetric ring from 50 Hz to 5 kHz and asserts both the closed form above and the monotonic decay towards zero.

Grounding impedance

The effective grounding impedance seen at the fault bus is

\[ Z_G(f) = \frac{u_{\text{EPR}}(f)}{r(f)\,I_{\text{fault}}(f)}. \]

It is exposed per frequency and as RMS-scalar through net.res_all_impedances().

Active flag and outage studies

Both Bus and Branch carry a boolean active field (default True). The flag has a clean physical interpretation:

• An inactive Bus is removed from the nodal system entirely. Its row and column drop from \(Y(f)\) and the bus contributes nothing to the right-hand-side vector \(\underline{i}(f)\).
• An inactive Branch behaves as an open circuit: no contribution to the admittance matrix, no Norton-equivalent injection from the mutual coupling, and a zero current in the result.

PathFinder skips inactive elements when enumerating source-to-fault paths, so the per-path Norton bookkeeping stays consistent. The flag is round-tripped through SQLite and JSON; existing payloads load with active=True for every element, which keeps backwards compatibility intact.

That makes maintenance scenarios, planned outages, broken shields and N-1 contingencies expressible without rebuilding the network. The groundinsight.simulation.outage sub-package wraps this into two convenience entry points:

• outage_context(network, outage) — context manager that flips the listed elements to active=False for the duration of a with block and restores the previous state (including the cached path list) afterwards.
• run_outage_study(network, fault, scenarios=[...]) — executes the base case plus one fault calculation per Outage scenario and returns an OutageStudyResult whose compare_buses(...) / compare_branches(...) accessors yield long-format Polars DataFrames with absolute and relative deltas against a chosen reference scenario.

Inverse rho analysis

The forward solve answers "given \(\rho_E\), what is the EPR?" The sister inverse question — "how large can \(\rho_E\) become before the EPR at the fault bus exceeds a touch-voltage limit \(u_{\max}\)?" — is answered by groundinsight.analysis.find_max_rho_scaling. It log-bisects a uniform scaling factor \(c\) of Bus.specific_earth_resistance on a user-selected bus set, re-evaluates each BusType.impedance_formula through the existing SymPy machinery and triggers run_fault at every trial value. The output is the largest \(c\) that still satisfies \(|U_\text{EPR}(f)|_{\text{RMS}} \le u_{\max}\), the EPR at that point, and the per-bus \(\rho_{\max} = c\,\rho_0\). The original \(\rho\) values are restored via a finally block, so the network is unchanged after the call. A frequency-dependent variant (find_max_rho_f_scaling) extends the same idea to two-parameter rho-f curves.

Summary of the calculation pipeline

run_fault(network, fault_name) executes the following steps (see network_operations.run_fault and ElectricalNetwork for details):

1. set_active_fault — select the target fault.
2. define_paths — if no paths are set yet, call PathFinder to enumerate them with a DFS.
3. build_electrical_network — create an ElectricalNetwork helper that holds the numerical arrays.
4. solve_network — build \(Y(f)\) and \(\underline{i}(f)\) and solve the linear system per frequency (sparse LU).
5. compute_branch_currents — derive branch currents from \(\Delta u\) plus Norton contributions.
6. compute_reduction_factors — re-solve without mutual Norton sources and take the EPR ratio at the fault bus.
7. compute_grounding_impedance — evaluate \(Z_G\) per frequency and as RMS.

The persistent state of the calculation (sparse matrices, vectors) lives on a private ElectricalNetwork attribute of the Network instance, so re-solving after topology changes only requires a fresh run_fault call.