Cache for non-linear methods.

Fields

• x::Vector{Float64}: the current state vector.
• r::Vector{Float64}: the current residual.
• Δx_nr::Vector{Float64}: the step under the Newton-Raphson method.

The remainder of the fields are only used in the TrustRegionACPowerFlow:

• r_predict::Vector{Float64}: the predicted residual at x+Δx_proposed, under a linear approximation: i.e J_x⋅(x+Δx_proposed).
• Δx_proposed::Vector{Float64}: the suggested step Δx, selected among Δx_nr, Δx_cauchy, and the dogleg interpolation between the two. The first is chosen when x+Δx_nr is inside the trust region, the second when both x+Δx_cauchy and x+Δx_nr are outside the trust region, and the third when x+Δx_cauchy is inside and x+Δx_nr outside. The dogleg step selects the point where the line from x+Δx_cauchy to x+Δx_nr crosses the boundary of the trust region.
• Δx_cauchy::Vector{Float64}: the step to the Cauchy point if the Cauchy point lies within the trust region, otherwise a step in that direction.

Accept a trust region step: update cached residual and autoscale vector d. The caller is responsible for recomputing the Jacobian via J(time_step) before calling this, so that Jv reflects the new state.

Check error and do refinement.

Sets Δx_proposed equal to the Δx by which we should update x. Decides between the Cauchy step Δx_cauchy, Newton-Raphson step Δx_nr, and the dogleg interpolation between the two, based on which fall within the trust region.

Log final residual, report convergence, compute optional post-processing factors, and return true/false. Shared by all AC power flow drivers.

Attempt Iwamoto damping on a rejected trust region step.

Uses the already-evaluated trial-point residual to compute an optimal damped step. Returns true if the damped step was accepted, false if reverted.

Compute the optimal Iwamoto step multiplier μ ∈ [IWAMOTOMUMIN, IWAMOTOMUMAX] by minimizing g(μ) = (1-μ)²g₀ + 2μ²(1-μ)g₁ + μ⁴g₂.

The stationary points satisfy the cubic g'(μ)/2 = 2g₂μ³ - 3g₁μ² + (g₀+2g₁)μ - g₀ = 0. All real roots are found analytically via the depressed-cubic trigonometric/Cardano form, and the global minimizer of g over the domain is returned. O(1), zero-allocation.

Evaluate the Iwamoto objective g(μ) = ‖(1-μ)f₀ + μ²f₁‖² expanded as g(μ) = (1-μ)²g₀ + 2μ²(1-μ)g₁ + μ⁴g₂ where g₀=‖f₀‖², g₁=f₀ᵀf₁, g₂=‖f₁‖².

Does a single iteration of Newton-Raphson with Iwamoto step control. Computes the Newton step, takes a full trial step, and checks whether the residual norm decreased. If not, computes an optimal damping multiplier μ and applies a damped step instead. When the damped step also fails to reduce the residual, the step is reverted to avoid divergence.

Returns true if the step made progress (residual decreased), false if the step was reverted. Consecutive reverts signal stagnation and the caller should terminate early.

Runs the full NewtonRaphsonACPowerFlow.

Keyword arguments:

• maxIterations::Int: maximum iterations. Default: 50.
• tol::Float64: tolerance. The iterative search ends when norm(abs.(residual)) < tol. Default: 1.0e-9.
• refinement_threshold::Float64: If the solution to J_x Δx = r satisfies norm(J_x Δx - r, 1)/norm(r, 1) > refinement_threshold, do iterative refinement to improve the accuracy. Default: 0.05.
• refinement_eps::Float64: run iterative refinement on J_x Δx = r until norm(Δx_{i}-Δx_{i+1}, 1)/norm(r,1) < refinement_eps. Default: 1.0e-6

Runs the full TrustRegionNRMethod.

Keyword arguments:

• maxIterations::Int: maximum iterations. Default: 50.
• tol::Float64: tolerance. The iterative search ends when maximum(abs.(residual)) < tol. Default: 1.0e-9.
• factor::Float64: the trust region starts out with radius factor*norm(x_0, 1), where x_0 is our initial guess, taken from data. Default: 1.0.
• eta::Float64: improvement threshold. If the observed improvement in our residual exceeds eta times the predicted improvement, we accept the new x_i. Default: 0.0001.
• iwamoto_fallback::Bool: when a trust region step is rejected, attempt Iwamoto damping to salvage the step before reverting. Default: true.

Sets the Newton-Raphson step. Usually, this is just J.Jv \ stateVector.r, but J.Jv might be singular.

Does a single iteration of NewtonRaphsonACPowerFlow. Updates the r and x fields of the stateVector, and computes the Jacobian at the new x.

Returns a stand-in matrix for singular J's.

Solve for the Newton-Raphson step, given the factorization object for J.Jv (if non-singular) or its stand-in (if singular).

Does a single iteration of the TrustRegionNRMethod: updates the x and r fields of the stateVector and computes the value of the Jacobian at the new x, if needed. Unlike _simple_step, this has a return value, the updated value of delta`.

If μ ∈ [IWAMOTOMUMIN, IWAMOTOMUMAX] and g(μ) < bestg, return the improved (μ, g(μ)); otherwise return (bestμ, best_g) unchanged.

AdamConfig(; learning_rate=0.01, beta1=0.9, beta2=0.999, epsilon=1e-8)
AdamConfig(settings::Dict{Symbol, Any})

Configuration for the Adam optimizer used by GradientDescentACPowerFlow.

AdamState(n::Int)

Pre-allocated mutable state for the Adam optimizer. All working arrays are allocated once before the iteration loop to ensure a zero-allocation inner loop.

Fields

• g::Vector{Float64}: gradient vector (length n)
• m::Vector{Float64}: 1st moment estimate (length n)
• v::Vector{Float64}: 2nd moment estimate (length n)
• t::Int: step counter for bias correction

Interpolate x between x_save (α=0) and current x (α=1) at fraction α. Computes x .= x_save .+ α .* (x .- x_save) in-place without allocation.

Driver for the GradientDescentACPowerFlow method: sets up the data structures, runs the Adam-based power flow method with backtracking line search, then handles post-processing.

adam_step!(x::Vector{Float64}, state::AdamState, cfg::AdamConfig)

Applies one Adam update to x in-place. All intermediate quantities are written directly into state.m and state.v; bias-corrected values are computed as scalars.

compute_gradient!(state::AdamState, J::ACPowerFlowJacobian, R::ACPowerFlowResidual)

Overwrites state.g with Jᵀ·F in-place. Uses mul! with the Transpose wrapper so the CSC column structure is traversed without materializing a new sparse matrix.

Compute value of gradient and Hessian at x.

Does A += B' * B, in a way that preserves the sparse structure of A, if possible. A workaround for the fact that Julia seems to run dropzeros!(A) automatically if I just do A .+= B' * B.

_update_hessian_matrix_values!(
    Hv::SparseMatrixCSC{Float64, Int32},
    F_value::Vector{Float64},
    data::ACPowerFlowData,
    time_step::Int64
)

Update the Hessian matrix values for the robust homotopy power flow solver.

Description

This function sets Hv equal to:

\[\sum_{k=1}^{2n} F_k(x) H_{F_k}(x)\]

where $F_k$ denotes the $k$th power balance equation and $H_{F_k}$ denotes its Hessian matrix.

This computes only the terms in the Hessian that come from the second derivatives of the power balance equations. The full Hessian of the objective function also includes a $J^T J$ term, which is computed separately.

Sparse Structure

The Hessian is organized into 2×2 blocks, each corresponding to a pair of buses. For a pair of buses $i$ and $k$ connected by a branch, the sparse structure of their block depends on the bus types:

\[\begin{array}{c|cc|cc|cc} & \text{REF} & & \text{PV} & & \text{PQ} & \\ & P_i & Q_i & Q_i & V_i & V_i & \theta_i \\ \hline \text{REF: } P_k & & & & & & \\ Q_k & & & & & & \\ \hline \text{PV: } Q_k & & & & & & \\ V_k & & & & \bullet & \bullet & \bullet \\ \hline \text{PQ: } V_k & & & & \bullet & \bullet & \bullet \\ \theta_k & & & & \bullet & \bullet & \bullet \end{array}\]

where $\bullet$ represents a potentially non-zero entry.

Diagonal blocks (where $i = k$) follow the same pattern as if each bus is its own neighbor. Off-diagonal blocks for pairs of buses not connected by a branch are structurally zero.

Arguments

• Hv::SparseMatrixCSC{Float64, Int32}: The Hessian matrix to be updated (modified in-place).
• F_value::Vector{Float64}: Current values of the power balance residuals.
• data::ACPowerFlowData: The power flow data containing bus and network information.
• time_step::Int64: The time step for which to compute the Hessian.

Driver for the LevenbergMarquardtACPowerFlow method: sets up the data structures (e.g. residual), runs the power flow method via calling _run_power_flow_method on them, then handles post-processing (e.g. loss factors).

In order to in-place modify the numeric values of a CHOLMOD matrix, we need to write our own wrapper around CHOLMOD.Sparse.

set_values!(mat::FixedStructureCHOLMOD, new_vals::AbstractVector{Float64})

In-place update of the numeric values in the CHOLMOD matrix.

A cached linear solver using KLU. Carefully written so as to minimize allocations: solve! and numeric_refactor! are completely non-allocating.

Fields:

• K: the underlying KLU object.
• reuse_symbolic::Bool: reuse the symbolic factorization. Defaults to true.
• check_pattern::Bool: if true, numeric_refactor! verifies that the new matrix has the same sparsity structure. Defaults to true.
• rf_common, rf_symbolic, rf_numeric: internal usage. Stored to avoid allocations.

Constructor.

Frees numeric factorization stored by cache, if non-null. If cache.check_pattern is true and the sparse matrix structure of A doesn't match the cached one, throws an error. Finally, computes the numeric factorization of A and stores that to cache.

symbolic_factor!(cache::KLULinSolveCache{T},  A::SparseMatrixCSC{Float64, T})

Frees up the current symbolic and numeric factorizations stored by cache, if non-null. Then computes the symbolic factorization of A and stores that to cache.

Symbolic refactor. Behavior depends on the values of cache.reuse_symbol and cache.check_pattern. There are 3 cases:

• !reuse_symbol: always refactor. Just calls symbolic_factor(cache, A).
• reuse_symbol && check_pattern: checks if the symbolic structure of A matches the cached one, and throws an error if it doesn't. This is to prevent bad input: we expected the structure to be the same, but it isn't.
• reuse_symbol && !check pattern: do nothing. Assume the structure of A matches the cached one.

Abstract supertype for all cached linear solvers. Subtypes must implement: symbolic_factor!, symbolic_refactor!, numeric_refactor! (which doubles as numeric_factor!), and solve!.

Build the full generator list including optional sources, storages, condensers, and HVDC synthetics.

Build v33 switched shunt step data (N, B pairs padded to 8).

Build v35 switched shunt step data (S, N, B triplets padded to 8).

Build all transformer-related metadata mappings (names, control objectives, winding groups, impedance, taps).

Calculate the STAT field for a 3-winding transformer based on per-winding availability.

Collect winding data for a 3-winding transformer.

Compute active power limits considering HVDC scaling.

Compute common DC line fields (record 1) for Two-Terminal DC export.

Compute inverter-side fields for Two-Terminal DC export.

Compute rectifier-side fields for Two-Terminal DC export.

Compute generator active and reactive power considering HVDC scaling.

Compute reactive power limits considering HVDC scaling.

Compute VSC converter fields for one side (from or to) of a VSC DC line.

Try to make an informative one or two character name for the load/generator/etc.

• generator-1234-AB -> AB
• 123_CT_7 -> 7
• load1234 -> 34

Setting a value of zero 0.0 when having a value greater than or equal to INFINITE_BOUND reverses the operation done in the PSY parsing side, according to PSSE Manual.

Load transformer components and create circuit ID mappings.

Returns a tuple of:

• transformerswithnumbers: 2-winding transformers with their bus numbers
• transformers3wwith_numbers: 3-winding transformers with their bus numbers
• transformercktmapping: Circuit ID mapping for 2-winding transformers
• transformer3wckt_mapping: Circuit ID mapping for 3-winding transformers

Create a synthetic generator (PSY.ThermalStandard) representing one end of a TwoTerminalGenericHVDCLine for export purposes. The generator is initialized with parameters reflecting the HVDC line's state.

Notes

- The generator's name is constructed as "<hvdc_line_name>_<suffix>".
- The `ext` field includes `"HVDC_END"` to indicate the end ("FR"/"TO").

Validate that the Sienna area/zone names parse as PSS/E-compatible area/zone numbers, output a mapping

Given a vector of Sienna bus names, create a dictionary from Sienna bus name to PSS/E-compatible bus name. Guarantees determinism and minimal changes.

WRITTEN TO SPEC: PSS/E 33.3 POM 5.2.1 Bus Data

Given a vector of Sienna bus numbers, create a dictionary from Sienna bus number to PSS/E-compatible bus number. Assumes that the Sienna bus numbers are positive and unique. Guarantees determinism: if the input contains the same bus numbers in the same order, the output will. Guarantees minimal changes: that if an existing bus number is compliant, it will not be changed.

WRITTEN TO SPEC: PSS/E 33.3 POM 5.2.1 Bus Data

Given a vector of Sienna transformer names, create a dictionary from Sienna transformer name to PSS/E-compatible transformer name. Guarantees determinism and minimal changes.

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Transformer Data

Parse a value to Int64, handling strings, floats, and PSSEDEFAULT. Returns PSSEDEFAULT for non-whole numbers or invalid values.

Update the parameters of synthetic generators created from HVDC lines, so they reflect the current setpoints and limits of the HVDC devices in the system.

Write the first record line for a 2-winding transformer.

Write the second record line (impedance data) for a 2-winding transformer.

Write the third record line (winding 1 data) for a 2-winding transformer.

Write the fourth record line (winding 2 data) for a 2-winding transformer.

Write the second record line (impedance data) for a 3-winding transformer.

Write winding records for a 3-winding transformer.

Write a DiscreteControlledACBranch record as a non-transformer branch (v33 path).

Write generator record for PSS/E v33 format.

Write generator record for PSS/E v35 format.

Write impedance correction table points in v33 format (T, F pairs).

Write impedance correction table points in v35 format (T, Re(F), Im(F) triplets).

Write a regular (Line/MonitoredLine) branch record to the buffer.

Write a TapTransformer with UNDEFINED control objective as a non-transformer branch record.

Temporary, very specialized proof of concept patch for https://github.com/JuliaLang/julia/issues/55835

Throw a NotImplementedError if the psse_version is not supported

If val is empty, returns T(); if not, asserts that val isa T and returns val. Has nice type checker semantics.

Examples

convert_empty(Vector{String}, [])  # -> String[]
convert_empty(Vector{String}, ["a"])  # -> ["a"]
convert_empty(Vector{String}, [2])  # -> TypeError: in typeassert, expected Vector{String}, got a value of type Vector{Int64}
Base.return_types(Base.Fix1(convert_empty, Vector{String}))  # -> [Vector{String}]

Given a vector of component names and a corresponding vector of container IDs (e.g., bus numbers), create unique-per-container PSS/E-compatible IDs, output a dictionary from (container ID, component name) to PSS/E-compatible component ID. The "singlesto1" flag detects components that are the only one on their bus and gives them the name "1".

Expand a single PowerFlowEvaluationModel into its possibly multiple parts for separate evaluation. Namely, if pfem contains a non-nothing exporter, return [pfem, exporter], else return [pfem].

Collects all AC branches (Line, MonitoredLine, DiscreteControlledACBranch) from the system, sorts them by their bus numbers, and returns a vector of tuples (branch, bus_numbers).

Arguments

• exporter::PSSEExporter: The exporter containing the system.

Returns

• Vector{Tuple{<:PSY.Branch, Tuple}}: Each tuple contains a branch and its associated bus numbers.

Force all cached information (serialized metadata, component lists, etc.) to be regenerated

Take the output of create_component_ids and make it more suitable for JSON serialization

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Bus Data. Sienna voltage limits treated as PSS/E normal voltage limits; PSSE emergency voltage limits left as default.

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Case Identification Data

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 FACTS Device Data

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Fixed Bus Shunt Data

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Generator Data

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Load Data

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Non-Transformer Branch Data

WRITTEN TO SPEC: PSS/E 33.3 POM 5.2.1 Q Record

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Switched Shunt Data

WRITTEN TO SPEC: PSS/E 35.4 POM 5.2.1 System Switching Device Data

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Transformer Data

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Transformer Impedance Correction Tables

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Two-Terminal DC Transmission Line Data

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Voltage Source Converter (VSC) DC Transmission Line Data

WRITTEN TO SPEC: PSS/E 33.3/35.4 POM 5.2.1 Zone Data

Write v35 header comments for a given section if applicable.

Apply BranchFlowEntry results to branch objects. Entries and branches are matched by iterating in the same order they were generated. Returns the number of entries consumed.

Non-AC: distribute pre-computed arc-level flows to individual branches. Losses are computed per-segment inside _distribute_arc_flows.

AC: recompute per-segment flows from solved voltages using _compute_segment_flows.

Returns a dictionary of bus index to power contribution at that bus from FixedAdmittance components, as a tuple of (active power, reactive power).

_compute_segment_flows(arc_entry, data, arc, time_step) -> Vector{BranchFlowEntry}

Compute per-segment branch flow entries from arc-level data and endpoint voltages. Dispatches on the arc entry type (direct, 3WT, parallel, series).

Distribute pre-computed arc-level flows to individual branches for non-AC power flow. Returns a Vector{BranchFlowEntry}, analogous to _compute_segment_flows for AC. Losses are computed per-segment as R * flow^2.

_get_arc_endpoint_voltages(data, arc, time_step)

Look up the complex voltages at the two endpoints of an arc from the solved power flow data.

_segment_flow_entry(segment, V_from, V_to)

Compute a BranchFlowEntry for a single segment given its endpoint voltages. Returns the from-to and to-from complex power flows, plus losses.

_set_series_interior_voltages!(sys, segment_sequence, equivalent_arc, V_endpoints, temp_bus_map)

Set the voltages at interior buses of a series chain from the solved interior voltages.

Method

Number the nodes in the series segment 0, 1, ..., n. Number the segments by their concluding node: 1, 2, ... n. The currents in the segments are given by:

\[\begin{bmatrix} y^i_{ff} & y^i_{ft} \\ y^i_{tf} & y^i_{tt} \end{bmatrix} \begin{bmatrix} V_{i-1} \\ V_i \end{bmatrix} = \begin{bmatrix} I_{i-1, i} \\ I_{i, i-1} \end{bmatrix}\]

where upper indices denote the segment number.

There are no loads or generators at the internal nodes, so $I_{i, i+1} + I_{i, i-1} = 0$. Substitute the above expressions for the currents and group by $V_i$:

\[y^i_{tf} V_{i-1} + (y_{tt}^i + y_{ff}^{i+1}) V_i + y_{ft}^{i+1} V_{i+1} = 0\]

For $i = 1$ and $i = n-1$, move the terms involving $V_0$ and $V_n$ (known) to the other side. This gives a tridiagonal system for $x = [V_1, \ldots, V_{n-1}]$:

\[A x = [-y^1_{tf} V_0, 0, \ldots, 0, -y^{n}_{ft} V_n]\]

where $A$ has diagonal entries $y_{tt}^i + y_{ff}^{i+1}$, subdiagonal entries $y_{tf}^{i+1}$, and superdiagonal entries $y_{ft}^i$.

In the implementation, $y_{11}$ is used instead of $y_{ff}$, $y_{12}$ instead of $y_{ft}$, etc.

_solve_series_interior_voltages(segment_sequence, equivalent_arc, V_endpoints)

Solve for the complex voltages at the interior nodes of a series chain by constructing and solving the tridiagonal system. Returns a Vector{ComplexF64} of length n-1 where n is the number of segments (i.e. one voltage per interior node). See the docstring of _set_series_interior_voltages! for the mathematical derivation.

Return the names of the arcs in the power flow data. Each arc is named by its from-to bus number pair, e.g. "123-456".

 update_system!(sys::PSY.System, data::PowerFlowData; time_step = 1)

Modify the values in the given System to correspond to the given PowerFlowData such that if a new PowerFlowData is constructed from the resulting system it is the same as data. See also write_power_flow_solution!. NOTE this assumes that data was initialized from sys and then solved with no further modifications.

Updates system voltages and powers with power flow results

write_results(
    ::ACPowerFlow{<:ACPowerFlowSolverType},
    sys::PSY.System,
    data::ACPowerFlowData,
    time_step::Int64,
) -> Dict{String, DataFrames.DataFrame}

Returns a dictionary containing the AC power flow results.

Only single-period evaluation is supported at the moment for AC Power flows. The resulting dictionary will therefore feature just one key linked to one DataFrame.

Arguments:

• ::ACPowerFlow: use ACPowerFlow() storing AC power flow results.
• sys::PSY.System: container storing the system information.
• result::Vector{Float64}: vector containing the results for one single time-period.

write_results(
    data::Union{PTDFPowerFlowData, vPTDFPowerFlowData, ABAPowerFlowData},
    sys::PSY.System,
)

Returns a dictionary containing the DC power flow results. Each key corresponds to the name of the considered time periods, storing a DataFrame with the power flow results.

Arguments:

• data::Union{PTDFPowerFlowData, vPTDFPowerFlowData, ABAPowerFlowData}: PowerFlowData structure containing power flows and bus angles.
• sys::PSY.System: A PowerSystems.System object storing the system information.

Return set of all bus numbers that can be PV: i.e. have an available generator, or certain voltage regulation devices.

Return the active and reactive power generation from a generator component. It's pg=0 as default for synchronous condensers since there's no field in the component for active power.

Return the active power limits that should be used in power flow calculations and PSS/E exports. Redirects to PSY.get_active_power_limits in all but special cases.

Return the reactive power limits that should be used in power flow calculations and PSS/E exports. Redirects to PSY.get_reactive_power_limits in all but special cases.

Return set of all bus numbers that must be PV: i.e. have an available generator.

Compute per-bus active power range Rk = sum(Pmax - P_setpoint) for generators at REF/PV buses. Used for headroom-proportional distributed slack.

Only writes to column 1: PF uses single-value active power limits and setpoints (not time-varying). The caller copies column 1 to all time steps. For time-varying headroom, compute slack weights in PSI.

Similar to above: A*X where X is a matrix.

Matrix multiplication A*x. Written this way because a VirtualPTDF matrix does not store all of its entries: instead, it calculates them (or retrieves them from cache), one element or one row at a time.

For pretty printing floats in debugging messages.

Weighted dot product of two vectors.

Weighted norm of two vectors.

Check if a device has attribute 'active_power' for active power consumption or generation.