Polar formulation
Generic templates
ExaPF.AbstractVariable
— TypeAbstractVariable
Variables corresponding to a particular formulation.
ExaPF.AbstractFormulation
— TypeAbstractFormulation
Interface between the data and the mathemical formulation.
ExaPF.State
— TypeState <: AbstractVariable
All variables $x$ depending on the variables Control
$u$ through the non-linear equation $g(x, u) = 0$.
ExaPF.Control
— TypeControl <: AbstractVariable
Independent variables $u$ used in the reduced-space formulation.
Structure and variables
ExaPF.PolarForm
— TypePolarForm{T, IT, VT, MT} <: AbstractPolarFormulation
Implement the polar formulation associated to the network's equations.
Wrap a PS.PowerNetwork
network to load the data on the target device (CPU()
and CUDABackend()
are currently supported).
Example
julia> const PS = ExaPF.PowerSystem;
julia> network_data = PS.load_case("case9.m");
julia> polar = PolarForm(network_data, ExaPF.CPU())
Polar formulation (instantiated on device CPU(false))
Network characteristics:
#buses: 9 (#slack: 1 #PV: 2 #PQ: 6)
#generators: 3
#lines: 9
giving a mathematical formulation with:
#controls: 5
#states : 14
ExaPF.BlockPolarForm
— TypeBlockPolarForm{T, IT, VT, MT} <: AbstractFormulation
Block polar formulation: duplicates k
different polar models to evaluate them in parallel.
ExaPF.load_polar
— Functionload_polar(case, device=CPU(); dir=PS.EXADATA)
Load a PolarForm
instance from the specified benchmark library dir
on the target device
(default is CPU
). ExaPF uses two different benchmark libraries: MATPOWER (dir=EXADATA
) and PGLIB-OPF (dir=PGLIB
).
Examples
julia> polar = ExaPF.load_polar("case9")
Polar formulation (instantiated on device CPU(false))
Network characteristics:
#buses: 9 (#slack: 1 #PV: 2 #PQ: 6)
#generators: 3
#lines: 9
giving a mathematical formulation with:
#controls: 5
#states : 14
ExaPF.NetworkStack
— TypeNetworkStack{VT,VD,MT} <: AbstractNetworkStack{VT}
NetworkStack(polar::PolarForm)
NetworkStack(nbus::Int, ngen::Int, nlines::Int, VT::Type)
Store the variables associated to the polar formulation. The variables are stored in the field input
, ordered as follows
input = [vmag ; vang ; pgen]
The object stores also intermediate variables needed in the expression tree, such as the LKMR basis ψ
.
Notes
The NetworkStack can be instantiated on the host or on the target device.
Examples
julia> polar = ExaPF.load_polar("case9");
julia> stack = ExaPF.NetworkStack(polar)
21-elements NetworkStack{Vector{Float64}}
julia> stack.vmag
9-element Vector{Float64}:
1.0
1.0
1.0
1.0
1.0
1.0
1.0
1.0
1.0
ExaPF.init!
— Functioninit!(polar::PolarForm, stack::NetworkStack)
Set stack.input
with the initial values specified in the base PS.PowerNetwork
object.
The state and the control are defined as mapping:
ExaPF.mapping
— Functionmapping(polar::PolarForm, ::Control)
Return the mapping associated to the Control()
in NetworkStack
according to the polar formulation PolarForm
.
Examples
julia> polar = ExaPF.load_polar("case9");
julia> mapu = ExaPF.mapping(polar, Control())
5-element Vector{Int64}:
1
2
3
20
21
mapping(polar::PolarForm, ::State)
Return the mapping associated to the State()
in NetworkStack
according to the polar formulation PolarForm
.
Examples
julia> polar = ExaPF.load_polar("case9");
julia> mapu = ExaPF.mapping(polar, State())
14-element Vector{Int64}:
11
12
13
14
15
16
17
18
4
5
6
7
8
9
Powerflow solver
ExaPF.run_pf
— Functionrun_pf(
polar::PolarForm, stack::NetworkStack;
rtol=1e-8, max_iter=20, verbose=0,
)
Solve the power flow equations $g(x, u) = 0$ w.r.t. the stack $x$, using the (NewtonRaphson
algorithm. The initial state $x$ is specified implicitly inside stack
, with the mapping mapping
associated to the polar formulation. The object stack
is modified inplace in the function.
The algorithm stops when a tolerance rtol
or a maximum number of iterations maxiter
is reached.
Arguments
polar::AbstractFormulation
: formulation of the power flow equationstack::NetworkStack
: initial values in the network
Examples
julia> polar = ExaPF.load_polar("case9");
julia> stack = ExaPF.NetworkStack(polar);
julia> conv = run_pf(polar, stack);
julia> conv.has_converged
true
julia> conv.n_iterations
4
ExaPF.nlsolve!
— Functionnlsolve!(
algo::NewtonRaphson,
jac::Jacobian,
stack::NetworkStack;
linear_solver=DirectSolver(jac.J),
nl_buffer=NLBuffer(size(jac, 2)),
)
Solve the nonlinear system of equations $g(x) = 0$ with a NewtonRaphson
algorithm. At each iteration, we update the variable $x$ as
\[ x_{k+1} = x_{k} - (∇g_k)^{-1} g(x_k) \]
till $\| g(x_k) \| < ε_{tol}$
In the implementation,
- the function $g$ is specified in
jac.func
, - the initial variable $x_0$ in
stack::NetworkStack
(with mappingjac.map
), - the Jacobian $∇g$ is computed automatically in
jac
, with automatic differentiation.
Note that stack
is modified inplace during the iterations of algorithm.
The Jacobian jac
should be instantied before calling this function. By default, the linear system $(∇g_k)^{-1} g(x_k)$ is solved using a LU factorization. You can specify a different linear solver by changing the optional argument linear_solver
.
Arguments
algo::NewtonRaphon
: Newton-Raphson object, storing the options of the algorithmjac::Jacobian
: Stores the function $g$ and its Jacobian $∇g$. The Jacobian is updated with automatic differentiation.stack::NetworkStack
: initial valueslinear_solver::AbstractLinearSolver
: linear solver used to compute the Newton stepnl_buffer::NLBuffer
: buffer storing the residual vector and the descent directionΔx
. Can be reused to avoid unecessary allocations.
Examples
julia> polar = ExaPF.load_polar("case9");
julia> powerflow = ExaPF.PowerFlowBalance(polar) ∘ ExaPF.PolarBasis(polar);
julia> jx = ExaPF.Jacobian(polar, powerflow, State());
julia> stack = ExaPF.NetworkStack(polar);
julia> conv = ExaPF.nlsolve!(NewtonRaphson(), jx, stack);
julia> conv.has_converged
true
julia> conv.n_iterations
4
ExaPF.NewtonRaphson
— TypeNewtonRaphson <: AbstractNonLinearSolver
Newton-Raphson algorithm.
Attributes
maxiter::Int
(default 20): maximum number of iterationstol::Float64
(default1e-8
): tolerance of the algorithmverbose::Int
(default0
): verbosity level
Constraints
The different parts of the polar formulation are implemented in the following AbstractExpression
:
ExaPF.PolarBasis
— TypePolarBasis{VI, MT} <: AbstractExpression
PolarBasis(polar::AbstractPolarFormulation)
Implement the LKMR nonlinear basis. Takes as input the voltage magnitudes vmag
and the voltage angles vang
and returns
\[ \begin{aligned} & \psi_\ell^C(v, \theta) = v^f v^t \cos(\theta_f - \theta_t) \quad \forall \ell = 1, \cdots, n_\ell \\ & \psi_\ell^S(v, \theta) = v^f v^t \sin(\theta_f - \theta_t) \quad \forall \ell = 1, \cdots, n_\ell \\ & \psi_k(v, \theta) = v_k^2 \quad \forall k = 1, \cdots, n_b \end{aligned}\]
Dimension: 2 * n_lines + n_bus
Complexity
3 n_lines + n_bus
mul, n_lines
cos
and n_lines
sin
Examples
julia> polar = ExaPF.load_polar("case9");
julia> stack = ExaPF.NetworkStack(polar);
julia> basis = ExaPF.PolarBasis(polar)
PolarBasis (AbstractExpression)
julia> basis(stack)
27-element Vector{Float64}:
1.0
1.0
1.0
1.0
1.0
1.0
1.0
1.0
1.0
0.0
⋮
1.0
1.0
1.0
1.0
1.0
1.0
1.0
1.0
1.0
ExaPF.PowerFlowBalance
— TypePowerFlowBalance{VT, MT}
PowerFlowBalance(polar)
Implement a subset of the power injection corresponding to $(p_{inj}^{pv}, p_{inj}^{pq}, q_{inj}^{pq})$. The function encodes the active balance equations at PV and PQ nodes, and the reactive balance equations at PQ nodes:
\[\begin{aligned} p_i &= v_i \sum_{j}^{n} v_j (g_{ij}\cos{(\theta_i - \theta_j)} + b_{ij}\sin{(\theta_i - \theta_j})) \,, & ∀ i ∈ \{PV, PQ\} \\ q_i &= v_i \sum_{j}^{n} v_j (g_{ij}\sin{(\theta_i - \theta_j)} - b_{ij}\cos{(\theta_i - \theta_j})) \,. & ∀ i ∈ \{PQ\} \end{aligned}\]
Require composition with PolarBasis
.
Dimension: n_pv + 2 * n_pq
Complexity
2
SpMV
Examples
julia> polar = ExaPF.load_polar("case9");
julia> stack = ExaPF.NetworkStack(polar);
julia> powerflow = ExaPF.PowerFlowBalance(polar) ∘ ExaPF.PolarBasis(polar);
julia> round.(powerflow(stack); digits=6)
14-element Vector{Float64}:
-1.63
-0.85
0.0
0.9
0.0
1.0
0.0
1.25
-0.167
0.042
-0.2835
0.171
-0.2275
0.259
julia> run_pf(polar, stack); # solve powerflow equations
julia> isapprox(powerflow(stack), zeros(14); atol=1e-8)
true
ExaPF.VoltageMagnitudeBounds
— TypeVoltageMagnitudeBounds
Implement the bounds on voltage magnitudes not taken into account in the bound constraints. In the reduced space, this is associated to the the voltage magnitudes at PQ nodes:
\[v_{pq}^♭ ≤ v_{pq} ≤ v_{pq}^♯ .\]
Dimension: n_pq
Complexity
1
copyto
Note
In the reduced space, the constraints on the voltage magnitudes at PV nodes $v_{pv}$ are taken into account when bounding the control $u$.
Examples
julia> polar = ExaPF.load_polar("case9");
julia> stack = ExaPF.NetworkStack(polar);
julia> voltage_pq = ExaPF.VoltageMagnitudeBounds(polar);
julia> voltage_pq(stack)
6-element Vector{Float64}:
1.0
1.0
1.0
1.0
1.0
1.0
ExaPF.PowerGenerationBounds
— TypePowerGenerationBounds{VT, MT}
PowerGenerationBounds(polar)
Constraints on the active power productions and on the reactive power productions that are not already taken into account in the bound constraints. In the reduced space, that amounts to
\[p_{g,ref}^♭ ≤ p_{g,ref} ≤ p_{g,ref}^♯ ; C_g q_g^♭ ≤ C_g q_g ≤ C_g q_g^♯ .\]
Require composition with PolarBasis
.
Dimension: n_pv + 2 n_ref
Complexity
1
copyto, 1
SpMV
Examples
julia> polar = ExaPF.load_polar("case9");
julia> stack = ExaPF.NetworkStack(polar);
julia> run_pf(polar, stack); # solve powerflow equations
julia> power_generators = ExaPF.PowerGenerationBounds(polar) ∘ ExaPF.PolarBasis(polar);
julia> round.(power_generators(stack); digits=6)
4-element Vector{Float64}:
0.719547
0.24069
0.144601
-0.03649
ExaPF.LineFlows
— TypeLineFlows{VT, MT}
LineFlows(polar)
Implement thermal limit constraints on the lines of the network.
Require composition with PolarBasis
.
Dimension: 2 * n_lines
Complexity
4
SpMV, 4 * n_lines
quadratic, 2 * n_lines
add
Examples
julia> polar = ExaPF.load_polar("case9");
julia> stack = ExaPF.NetworkStack(polar);
julia> run_pf(polar, stack); # solve powerflow equations
julia> line_flows = ExaPF.LineFlows(polar) ∘ ExaPF.PolarBasis(polar);
julia> round.(line_flows(stack); digits=6)
18-element Vector{Float64}:
0.575679
0.094457
0.379983
0.723832
0.060169
0.588673
2.657418
0.748943
0.295351
0.560817
0.112095
0.38625
0.728726
0.117191
0.585164
2.67781
0.726668
0.215497
Objective
The production costs is given in the AbstractExpression
CostFunction
:
ExaPF.CostFunction
— TypeCostFunction{VT, MT} <: AutoDiff.AbstractExpression
CostFunction(polar)
Implement the quadratic cost function for OPF
\[ ∑_{g=1}^{n_g} c_{2,g} p_g^2 + c_{1,g} p_g + c_{0,g}\]
Require composition with PolarBasis
to evaluate the cost of the reference generator.
Dimension: 1
Complexity
1
SpMV, 1
sum
Examples
julia> polar = ExaPF.load_polar("case9");
julia> stack = ExaPF.NetworkStack(polar);
julia> cost = ExaPF.CostFunction(polar) ∘ ExaPF.PolarBasis(polar);
julia> cost(stack)
1-element Vector{Float64}:
4509.0275
Composition of expressions
The different expressions can be combined together in several different ways.
ExaPF.MultiExpressions
— TypeMultiExpressions <: AbstractExpression
Implement expressions concatenation. Takes as input a vector of expressions [expr1,...,exprN]
and concatenate them in a single expression mexpr
, such that
mexpr(x) = [expr1(x) ; expr2(x) ; ... ; exprN(x)]
ExaPF.ComposedExpressions
— TypeComposedExpressions{Expr1<:PolarBasis, Expr2} <: AbstractExpression
Implement expression composition. Takes as input two expressions expr1
and expr2
and returns a composed expression cexpr
such that ``` cexpr(x) = expr2 ∘ expr1(x)
Notes
Currently, only PolarBasis
is supported for expr1
.