Electrical network
Numerical core of groundinsight. ElectricalNetwork assembles
the nodal admittance matrix \(Y(f)\) per frequency, fills the
right-hand side with source and mutual-coupling Norton currents,
solves the linear system via sparse LU and derives branch
currents, reduction factors and grounding impedances.
Physical / modelling context
For each frequency \(f\) in the network's frequency list the model
reduces to a sparse complex linear system
\[
Y(f)\,\underline{u}(f) \;=\; \underline{i}(f),
\qquad
\underline{u}(f) \;=\; Y(f)^{-1}\,\underline{i}(f),
\]
with
- \(Y(f)\) the nodal admittance matrix. The diagonal entries are
the bus grounding admittances \(1/Z_E(f, \rho)\); the
off-diagonal entries are the branch self-admittances
\(1/Z_b(f, l)\) when
grounding_conductor=True.
- \(\underline{u}(f)\) the bus-EPR vector.
- \(\underline{i}(f)\) the source currents (scaled by the
fault scaling at \(f\)) plus the Norton equivalents of the
mutual coupling between phase and grounding conductor along
the source-to-fault path.
After the per-frequency solve the network derives:
- Branch currents from \(\Delta u\) and the impressed mutual
current.
- Reduction factors as the ratio
\(|\underline{u}_\text{with}|/|\underline{u}_\text{without}|\) at
the fault bus, where the second solve omits the mutual
Norton sources.
- Grounding impedance
\(Z_G = u_\text{EPR}/(r \cdot I_\text{fault})\) per frequency.
The ElectricalNetwork is held as a PrivateAttr of Network,
not exposed in the JSON / DB serialisation; it carries the raw
NumPy / SciPy working arrays used by the LU solver.
Example
import groundinsight as gi
# Assume `net` was built via gi.create_* factories with a fault `f1`.
gi.run_fault(network=net, fault_name="f1") # populates net._enet under the hood
# Inspect derived quantities through the public Network accessors:
df_bus = net.res_buses() # EPR per bus
df_branch = net.res_branches() # branch currents
df_zg = net.res_all_impedances() # grounding impedance + reduction factor
Direct use of ElectricalNetwork is rarely needed; access is
mostly via gi.run_fault and the Network.res_* DataFrame
accessors.
API reference
electrical_network
module for creating an electrical network based on the core models
The network is used to perform calculations based on the matrix form of the network:
Y * u = i
u = Y^-1 * i
where:
Y - admittance matrix, branches and buses are used to build this matrix
u - vector of bus voltages (earth potential rise per bus)
i - vector of source injections plus the Norton-equivalent injections that
represent the inductive coupling between phase conductors and shield
(grounding) conductors along each branch.
Sign convention for the mutual coupling (see _add_mutual_currents):
U_from - U_to = Z_self * I_s - Z_mutual * I_p
where I_p is the phase current through the branch in from->to direction and
I_s is the shield current in the same direction. The induced EMF enters the
nodal form as a Norton current i_mut = (Z_mutual / Z_self) * I_p which is
injected as +i_mut out of the from bus (nodal: i_vector[from] -= i_mut) and
+i_mut into the to bus (nodal: i_vector[to] += i_mut).
Two strategies are available to determine the phase current I_p per branch:
1) Path-based (default). For every simple path from source to fault, the
branch direction is derived from the actual path traversal (no index
heuristic). The user may scale the contribution per branch via
Branch.parallel_coefficient, which is the legacy knob for splitting
current between parallel paths.
2) Automatic distribution (auto_phase_currents=True). A reduced phase-only
network is solved per source with +I at the source bus and the fault bus
as reference. The resulting branch currents are used directly. In this
mode parallel_coefficient is ignored. This mode is topology-agnostic and
is the intended integration point for an external phase-current source
(e.g. pandapower single-phase short-circuit results).
ElectricalNetwork
ElectricalNetwork(
network: Network, auto_phase_currents: bool = False
)
Represents the electrical properties of a network, enabling calculations.
This class handles the construction of admittance matrices, voltage and current vectors,
and performs network analysis to compute results such as bus voltages, branch currents,
reduction factors, and grounding impedances.
Initialize the ElectricalNetwork with a given Network model.
Sets up necessary data structures and initializes the network calculations.
Parameters:
| Name |
Type |
Description |
Default |
network
|
Network
|
The Network instance containing buses, branches, sources, and faults.
|
required
|
auto_phase_currents
|
bool
|
If True the phase current through each
branch is computed by solving a reduced phase-only network (Variant B).
If False (default) the phase current is derived from the enumerated
source-to-fault paths using each branch's parallel_coefficient
(Variant A). Defaults to False.
|
False
|
Source code in src/groundinsight/electrical_network.py
| def __init__(self, network: Network, auto_phase_currents: bool = False):
"""
Initialize the ElectricalNetwork with a given Network model.
Sets up necessary data structures and initializes the network calculations.
Parameters
----------
network : Network
The Network instance containing buses, branches, sources, and faults.
auto_phase_currents : bool, optional
If True the phase current through each
branch is computed by solving a reduced phase-only network (Variant B).
If False (default) the phase current is derived from the enumerated
source-to-fault paths using each branch's ``parallel_coefficient``
(Variant A). Defaults to False.
"""
self.network = network
self.auto_phase_currents = auto_phase_currents
self.bus_indices = {}
self.Y_matrices = {} # Admittance matrices for each frequency
self.u_vectors = {} # Voltage vectors for each frequency
self.u_vectors_no_mutual = {} # Voltage vectors without mutual currents
self.i_vectors_no_mutual = {} # Current vectors without mutual currents
self.i_vectors = {} # Current vectors for each frequency
self.results: Result = Result() # Stores the calculation results
self.i_mutuals = {} # Store mutual currents per frequency per branch
self.total_source_currents = {} # Store total source currents per frequency
# Source-only nodal injection per frequency, i.e. the ``i`` vector *before*
# the mutual Norton equivalents are added on top. This is the current that
# physically enters the grounding system through a lumped connection at the
# bus -- the earthing conductor (EN 50522 "Erdungsleiter") -- as opposed to
# ``ResultBus.ia``, which is the share dissipated into the soil through the
# earth electrode ("Erder"). The mutual terms are a distributed modelling
# artefact of the line/shield coupling and are deliberately excluded here:
# no lumped conductor carries them into the node.
self.source_injections = {}
self.phase_currents = {} # Signed phase current per branch per frequency
# One-shot latch so the phase-impedance proxy warning is emitted once
# per network, not once per frequency of the sweep.
self._phase_proxy_warned = False
# Substitute impedance per frequency for elements that are a short
# circuit there. Only 0 Hz can have an entry, and only when the network
# actually contains such an element. See _dc_substitute_at.
self._dc_substitutes: Dict[float, float] = {}
self._initialize()
|
compute_branch_currents
compute_branch_currents()
Compute branch currents for each frequency and store them in the Result object.
This method calculates the current flowing through each branch based on the bus voltages
and branch impedances. The results are stored as ResultBranch instances within the
network's results object.
Source code in src/groundinsight/electrical_network.py
| def compute_branch_currents(self):
"""
Compute branch currents for each frequency and store them in the Result object.
This method calculates the current flowing through each branch based on the bus voltages
and branch impedances. The results are stored as `ResultBranch` instances within the
network's results object.
"""
fault_name = self.network.active_fault
if fault_name is None:
raise ValueError("No active fault set in the network.")
if fault_name not in self.network.results:
raise ValueError(f"No results available for fault '{fault_name}'.")
result = self.network.results[fault_name]
for branch in self.network.branches.values():
i_s_freq = {}
# Inactive branches and branches with at least one inactive endpoint
# are open-circuited: their shield current is zero by construction.
is_open = (
not branch.active
or branch.from_bus not in self.bus_indices
or branch.to_bus not in self.bus_indices
)
if is_open:
for freq in self.network.frequencies:
i_s_freq[freq] = ComplexNumber(real=0.0, imag=0.0)
rms_current = 0.0
result_branch = ResultBranch(
name=branch.name, i_s=rms_current, i_s_freq=i_s_freq
)
result.branches.append(result_branch)
continue
from_idx = self.bus_indices[branch.from_bus]
to_idx = self.bus_indices[branch.to_bus]
for freq in self.network.frequencies:
from_voltage = self.u_vectors[freq][from_idx]
to_voltage = self.u_vectors[freq][to_idx]
Z_self_complex = self._resolved_impedance(
branch.self_impedance.get(freq), freq
)
if (
Z_self_complex is not None
and branch.type.grounding_conductor
and not _is_open_circuit(Z_self_complex)
):
Y_self_complex = 1 / Z_self_complex
delta_voltage = to_voltage - from_voltage
# Stored mutual term already carries the sign to combine with
# delta_voltage = u_to - u_from. See _add_mutual_currents.
i_mutual = self.i_mutuals.get(freq, {}).get(branch.name, 0)
current = delta_voltage * Y_self_complex + i_mutual
i_s_freq[freq] = ComplexNumber(real=current.real, imag=current.imag)
else:
i_s_freq[freq] = ComplexNumber(real=0.0, imag=0.0)
# Calculate RMS current
rms_current = self._calculate_rms(i_s_freq)
result_branch = ResultBranch(
name=branch.name, i_s=rms_current, i_s_freq=i_s_freq
)
result.branches.append(result_branch)
# Update the result in the network's results dictionary
self.network.results[fault_name] = result
self.results = result # Update self.results
|
compute_grounding_impedance
compute_grounding_impedance()
Compute the grounding impedance for the fault bus.
This method calculates the grounding impedance using the formula::
where I_fault is the (signed) sum of source injections at the
active fault. For current-mode sources this is Σ scaling * I_src
as before; for Thevenin (voltage-mode) sources it is the corresponding
Norton injection Σ scaling * U_emf / Z_src. In Thevenin mode the
resulting Z_G is therefore the EPR per Norton ampere, which
depends on both the grounding network and Z_src; it recovers the
classic grounding impedance in the limit Z_src -> ∞ with
U_emf = I_src * Z_src held constant.
The results are stored in the network's results object.
Source code in src/groundinsight/electrical_network.py
| def compute_grounding_impedance(self):
"""
Compute the grounding impedance for the fault bus.
This method calculates the grounding impedance using the formula::
Z_G = Z_E = U_E / I_E
where ``I_fault`` is the (signed) sum of source injections at the
active fault. For current-mode sources this is ``Σ scaling * I_src``
as before; for Thevenin (voltage-mode) sources it is the corresponding
Norton injection ``Σ scaling * U_emf / Z_src``. In Thevenin mode the
resulting ``Z_G`` is therefore the EPR per Norton ampere, which
depends on both the grounding network and ``Z_src``; it recovers the
classic grounding impedance in the limit ``Z_src -> ∞`` with
``U_emf = I_src * Z_src`` held constant.
The results are stored in the network's results object.
"""
fault_name = self.network.active_fault
if fault_name is None:
raise ValueError("No active fault set in the network.")
fault_bus = self.network.faults[fault_name].bus
fault_bus_idx = self.bus_indices[fault_bus]
grounding_impedances = (
{}
) # Dictionary to store grounding impedance per frequency
frequencies = self.network.frequencies
result = self.network.results[fault_name]
# Ensure that reduction factors are computed
if not result.reduction_factor:
raise ValueError(
"Reduction factors not computed. Please compute reduction factors before grounding impedance."
)
for freq in frequencies:
# Z_E is the earthing impedance of the bonded earthing system in the
# EN 50522 sense: the earthing voltage of that system divided by the
# current it passes into the soil. It is computed alongside the
# current-based reduction factor, where the group is determined, and
# it closes the chain U_E = 3*I_0 * Z_E * r by construction.
z_earthing = result.reduction_factor.z_earthing.get(freq)
if z_earthing is not None:
grounding_impedances[freq] = ComplexNumber(
real=complex(z_earthing).real, imag=complex(z_earthing).imag
)
continue
# No earth-return current at this frequency -- nothing flows into
# the soil, so there is no earthing impedance to report. This is a
# different statement from an impedance of zero.
grounding_impedances[freq] = None
# Store the grounding impedance in the result
result_grounding_impedance = ResultGroundingImpedance(
fault_bus=fault_bus, value=grounding_impedances
)
result.grounding_impedance = result_grounding_impedance
# Update the result in the network's results dictionary
self.network.results[fault_name] = result
self.results = result # Update self.results
|
compute_reduction_factors
compute_reduction_factors()
Compute the reduction factors by solving the network with and without mutual currents.
This method calculates how much the presence of mutual currents affects the Earth Potential Rise (EPR).
The reduction factors are stored in the network's results object.
Source code in src/groundinsight/electrical_network.py
| def compute_reduction_factors(self):
"""
Compute the reduction factors by solving the network with and without mutual currents.
This method calculates how much the presence of mutual currents affects the Earth Potential Rise (EPR).
The reduction factors are stored in the network's results object.
"""
fault_name = self.network.active_fault
if fault_name is None:
raise ValueError("No active fault set in the network.")
fault_bus = self.network.faults[fault_name].bus
fault_bus_idx = self.bus_indices[fault_bus]
reduction_factors = {}
uepr_with_mutual = {}
uepr_without_mutual = {}
frequencies = self.network.frequencies
# Step 1: Solve network with mutual currents (already done)
# Voltages are stored in self.u_vectors
# Store uepr with mutual currents
for freq in frequencies:
voltage = self.u_vectors[freq][fault_bus_idx]
uepr_with_mutual[freq] = voltage
# Step 2: Create i_vectors without mutual currents
self._construct_vectors_no_mutual()
# Step 3: Solve network without mutual currents
self.u_vectors_no_mutual = {}
for freq in frequencies:
Y_matrix = self.Y_matrices[freq]
i_vector = self.i_vectors_no_mutual[freq]
try:
# Solve for u_vector without mutual currents
u_vector = np.linalg.solve(Y_matrix, i_vector)
self.u_vectors_no_mutual[freq] = u_vector
except np.linalg.LinAlgError as e:
logger.error(
"Error solving network equations at frequency %s without mutual currents: %s",
freq,
e,
exc_info=True,
)
continue
# Store uepr without mutual currents. A frequency whose no-mutual solve
# failed above is absent from ``u_vectors_no_mutual``; it is skipped here
# and reported as ``None`` rather than raising a bare KeyError out of the
# solver core.
for freq in frequencies:
if freq not in self.u_vectors_no_mutual:
continue
voltage = self.u_vectors_no_mutual[freq][fault_bus_idx]
uepr_without_mutual[freq] = voltage
# Step 4: Compute reduction factors
for freq in frequencies:
if freq not in uepr_without_mutual:
reduction_factors[freq] = None
continue
v_with = uepr_with_mutual[freq]
v_without = uepr_without_mutual[freq]
# Compute magnitudes
mag_with = abs(v_with)
mag_without = abs(v_without)
if mag_without != 0:
reduction_factor = mag_with / mag_without
else:
reduction_factor = None # Handle division by zero
reduction_factors[freq] = reduction_factor
# Step 5: the current-based reduction factor -- see the method docstring.
(
reduction_factors_current,
earth_currents,
earth_buses,
earthing_voltages,
earthing_impedances,
) = self._compute_current_reduction_factors(fault_bus_idx)
# The two definitions answer different questions and only coincide in
# the limit of ideally bonded stations. Saying so where they part
# company beats leaving a factor of forty to be discovered later.
self._warn_on_reduction_factor_divergence(
reduction_factors, reduction_factors_current
)
# Store the reduction factors in the result
result = self.network.results[fault_name]
result_reduction_factor = ResultReductionFactor(
fault_bus=fault_bus,
value=reduction_factors,
value_current=reduction_factors_current,
i_earth=earth_currents,
earth_buses=earth_buses,
u_earthing=earthing_voltages,
z_earthing=earthing_impedances,
)
result.reduction_factor = result_reduction_factor
# Update the result in the network's results dictionary
self.network.results[fault_name] = result
self.results = result # Update self.results
|
solve_network
Solve the network equations Y * u = i for each frequency.
This method computes the bus voltages by solving the admittance matrix equations for each frequency.
The results are stored in the network's results object.
It uses the csc_matrix and splu functions from scipy, assuming the Y-Matrix is a sparse matrix.
.. warning::
This **replaces** ``network.results[fault]`` with a fresh
:class:`Result` that carries bus rows only. The branch rows are
filled in by :meth:`compute_branch_currents`, which must be
called afterwards -- :func:`groundinsight.run_fault` does exactly
that. Calling ``solve_network`` on its own, e.g. on a hand-built
:class:`ElectricalNetwork` used to inspect ``Y``, ``i`` or ``u``,
therefore discards the branch results of a previous ``run_fault``.
Clearing them is deliberate -- stale branch currents next to fresh
bus voltages would be silently inconsistent -- and the thermal
checks raise on the resulting gap rather than reporting an
incomplete result as free of violations.
To inspect the nodal system without touching the stored results,
build the :class:`ElectricalNetwork` (its constructor has no side
effects on ``network.results``) and read ``u`` back from
``network.results[fault].buses`` instead of re-solving.
Source code in src/groundinsight/electrical_network.py
| def solve_network(self):
"""
Solve the network equations Y * u = i for each frequency.
This method computes the bus voltages by solving the admittance matrix equations for each frequency.
The results are stored in the network's results object.
It uses the csc_matrix and splu functions from scipy, assuming the Y-Matrix is a sparse matrix.
.. warning::
This **replaces** ``network.results[fault]`` with a fresh
:class:`Result` that carries bus rows only. The branch rows are
filled in by :meth:`compute_branch_currents`, which must be
called afterwards -- :func:`groundinsight.run_fault` does exactly
that. Calling ``solve_network`` on its own, e.g. on a hand-built
:class:`ElectricalNetwork` used to inspect ``Y``, ``i`` or ``u``,
therefore discards the branch results of a previous ``run_fault``.
Clearing them is deliberate -- stale branch currents next to fresh
bus voltages would be silently inconsistent -- and the thermal
checks raise on the resulting gap rather than reporting an
incomplete result as free of violations.
To inspect the nodal system without touching the stored results,
build the :class:`ElectricalNetwork` (its constructor has no side
effects on ``network.results``) and read ``u`` back from
``network.results[fault].buses`` instead of re-solving.
"""
fault_name = self.network.active_fault
if fault_name is None:
raise ValueError("No active fault set in the network.")
result = Result(buses=[], branches=[], fault=fault_name)
for freq in self.network.frequencies:
Y_matrix = self.Y_matrices[freq]
i_vector = self.i_vectors[freq]
singular_msg = (
f"Singular admittance matrix at f={freq} Hz: the network has no "
"path to reference earth. Ensure at least one bus has a finite "
"grounding impedance and that the fault bus is connected to it."
)
# A NaN anywhere in the nodal system poisons the factorisation, and
# scipy reports that as a singular matrix -- which reads like a
# topology error and sends the engineer looking for a missing
# earth connection that is not missing. Catch it first and name
# what is actually NaN.
self._assert_finite_system(freq, Y_matrix, i_vector)
# Exact, solver-independent floating-network guard. If no active bus
# is referenced to earth at this frequency (every grounding impedance
# is zero, infinite or missing) the admittance matrix is singular.
# scipy's sparse ``splu`` handles a singular Y inconsistently across
# versions -- it may raise ``RuntimeError``, return a non-finite
# solution, or return an arbitrary finite one (a floating network has
# no unique EPR) -- so the common floating case is caught structurally
# here, before the solve, with no numerical tolerance.
grouped = self._classify_bus_grounding(freq)
if not grouped["referenced"]:
raise ValueError(self._no_ground_reference_message(freq, grouped))
try:
# Backstops for any remaining singular case (e.g. a disconnected
# ungrounded island): scipy may raise, or return a non-finite
# solution.
Y_matrix_sparse = csc_matrix(Y_matrix)
lu = splu(Y_matrix_sparse)
u_vector = lu.solve(i_vector)
except (np.linalg.LinAlgError, RuntimeError) as e:
raise ValueError(f"{singular_msg} Original solver error: {e}") from e
if not np.all(np.isfinite(u_vector)):
raise ValueError(singular_msg)
self.u_vectors[freq] = u_vector
# Create ResultBus instances. Inactive buses are not in ``bus_indices``
# and therefore do not appear in the result -- they are physically
# disconnected from the nodal system.
for bus_name, idx in self.bus_indices.items():
uepr_freq = {}
ia_freq = {}
i_inj_freq = {}
for freq in self.network.frequencies:
voltage = self.u_vectors[freq][idx]
bus = self.network.buses.get(bus_name)
Z_self_complex = self._resolved_impedance(
bus.impedance.get(freq), freq
)
if Z_self_complex is None:
current = 0
else:
if _is_open_circuit(Z_self_complex):
# Open end (Z = inf): no electrode, so no current leaves
# into the soil here. Computing the quotient would give
# NaN and store it in the I_a column.
current = 0
else:
current = voltage / Z_self_complex
uepr_freq[freq] = ComplexNumber(real=voltage.real, imag=voltage.imag)
ia_freq[freq] = ComplexNumber(
real=complex(current).real, imag=complex(current).imag
)
# Source-only injection at this bus (0 at every bus that is
# neither a source bus nor the fault bus).
inj_vector = self.source_injections.get(freq)
injection = 0j if inj_vector is None else complex(inj_vector[idx])
i_inj_freq[freq] = ComplexNumber(
real=injection.real, imag=injection.imag
)
# Calculate RMS values
rms_voltage = self._calculate_rms(uepr_freq)
rms_current = self._calculate_rms(ia_freq)
rms_injection = self._calculate_rms(i_inj_freq)
result_bus = ResultBus(
name=bus_name,
uepr=rms_voltage,
ia=rms_current,
i_inj=rms_injection,
uepr_freq=uepr_freq,
ia_freq=ia_freq,
i_inj_freq=i_inj_freq,
)
result.buses.append(result_bus)
# Store the result in the network's results dictionary
self.network.results[fault_name] = result
self.results = result # Also keep a reference in self.results
|