FullSpaceEvaluator

FullSpaceEvaluator models the original OPF problem in the full-space.

Initialization

A FullSpaceEvaluator can be instantiated both from a MATPOWER file:

datafile = joinpath(INSTANCES_DIR, "case9.m")
flp = Argos.FullSpaceEvaluator(datafile)

A FullSpaceEvaluator object
    * device: KernelAbstractions.CPU(false)
    * #vars: 19
    * #cons: 36

or equivalently, from a ExaPF.PolarForm object:

polar = ExaPF.PolarForm(datafile)
flp = Argos.FullSpaceEvaluator(polar)

A FullSpaceEvaluator object
    * device: KernelAbstractions.CPU(false)
    * #vars: 19
    * #cons: 36

Note

One can remove the line-flow constraints in the model simply by adding the keyword argument flp = Argos.FullSpaceEvaluator(polar; line_constraints=false).

Attributes

FullSpaceEvaluator is just a thin wrapper on top of ExaPF, and comes with just a few attributes. One can query the original ExaPF model with

model = Argos.model(flp)

Polar formulation (instantiated on device KernelAbstractions.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

The initial variable:

x = Argos.initial(flp)

19-element Vector{Float64}:
 0.0
 0.0
 0.0
 0.0
 0.0
 0.0
 0.0
 0.0
 1.0
 1.0
 1.0
 1.0
 1.0
 1.0
 1.0
 1.0
 1.0
 1.63
 0.85

The dimensions of the problem:

(n, m) = Argos.n_variables(flp), Argos.n_constraints(flp)

(19, 36)

The bounds on the variables (x, u)

xlb, xub = Argos.bounds(flp, Argos.Variables())
[xlb xub]

19×2 Matrix{Float64}:
 -Inf   Inf
 -Inf   Inf
 -Inf   Inf
 -Inf   Inf
 -Inf   Inf
 -Inf   Inf
 -Inf   Inf
 -Inf   Inf
   0.9   1.1
   0.9   1.1
   0.9   1.1
   0.9   1.1
   0.9   1.1
   0.9   1.1
   0.9   1.1
   0.9   1.1
   0.9   1.1
   0.1   3.0
   0.1   2.7

The bounds on the constraints:

clb, cub = Argos.bounds(flp, Argos.Constraints())
[clb cub]

36×2 Matrix{Float64}:
   0.0  0.0
   0.0  0.0
   0.0  0.0
   0.0  0.0
   0.0  0.0
   0.0  0.0
   0.0  0.0
   0.0  0.0
   0.0  0.0
   0.0  0.0
   ⋮    
 -Inf   6.25
 -Inf   6.25
 -Inf   2.25
 -Inf   9.0
 -Inf   2.25
 -Inf   6.25
 -Inf   6.25
 -Inf   6.25
 -Inf   6.25

Sparse Jacobian and sparse Hessian

More importantly, FullSpaceEvaluator stores two AD backends to evaluate the Jacobian and the Hessian in sparse format using ExaPF.

For the Hessian, the AD backend is

flp.hess

ExaPF.FullHessian{ExaPF.PolarForm{Float64, Vector{Int64}, Vector{Float64}, Matrix{Float64}}, ExaPF.ComposedExpressions{ExaPF.PolarBasis{Vector{Int64}, SparseArrays.SparseMatrixCSC{Float64, Int64}}, ExaPF.MultiExpressions}, ExaPF.NetworkStack{Vector{Float64}, Vector{ForwardDiff.Dual{Nothing, Float64, 8}}, NamedTuple{(:c, :sfp, :sfq, :stp, :stq, :∂edge_vm_fr, :∂edge_vm_to, :∂edge_va_fr, :∂edge_va_to), NTuple{9, Vector{ForwardDiff.Dual{Nothing, Float64, 8}}}}}, Vector{ForwardDiff.Dual{Nothing, Float64, 8}}, SparseArrays.SparseMatrixCSC{Float64, Int64}, Vector{Int64}}(Polar formulation (instantiated on device KernelAbstractions.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.ComposedExpressions{ExaPF.PolarBasis{Vector{Int64}, SparseArrays.SparseMatrixCSC{Float64, Int64}}, ExaPF.MultiExpressions}(PolarBasis (AbstractExpression), ExaPF.MultiExpressions(ExaPF.AutoDiff.AbstractExpression[CostFunction (AbstractExpression), PowerFlowBalance (AbstractExpression), PowerGenerationBounds (AbstractExpression), LineFlows (AbstractExpression)])), [11, 12, 13, 14, 15, 16, 17, 18, 4, 5, 6, 7, 8, 9, 1, 2, 3, 20, 21], 21-elements NetworkStack{Vector{ForwardDiff.Dual{Nothing, Float64, 8}}}, 21-elements NetworkStack{Vector{ForwardDiff.Dual{Nothing, Float64, 8}}}, [1, 1, 1, 2, 3, 4, 2, 3, 4, 5, 6, 7, 5, 6, 7, 8, 8, 1, 1], 8, ForwardDiff.Dual{Nothing, Float64, 8}[Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0)  …  Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0)], ForwardDiff.Dual{Nothing, Float64, 8}[Dual{Nothing}(0.0,0.0,0.0,0.0,2.121995791e-314,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,2.121995791e-314,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0)  …  Dual{Nothing}(0.0,0.0,0.0,0.0,2.2754302785599115e-305,0.0,0.0,-6.798821520133666e303,0.0), Dual{Nothing}(0.0,6.463054385729356e243,0.0,0.0,-1.7141420012175656e210,0.0,0.0,8.72835702500226e-219,0.0), Dual{Nothing}(0.0,2.121995791e-314,0.0,0.0,0.0,0.0,0.0,9.07443846736522e166,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,2.121995791e-314,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,3.273434022483496e-172,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0), Dual{Nothing}(0.0,0.0,0.0,0.0,0.0,0.0,0.0,-5.662090190527122e-166,0.0)], sparse([1, 7, 13, 16, 2, 5, 11, 17, 3, 4  …  1, 7, 13, 16, 2, 5, 11, 17, 18, 19], [1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5, 6, 6, 6, 6, 6, 6, 7, 7, 7, 7, 7, 7, 7, 7, 8, 8, 8, 8, 8, 8, 9, 9, 9, 9, 9, 9, 9, 10, 10, 10, 10, 10, 10, 11, 11, 11, 11, 11, 11, 11, 11, 12, 12, 12, 12, 12, 12, 13, 13, 13, 13, 13, 13, 13, 13, 14, 14, 14, 14, 14, 14, 15, 15, 15, 16, 16, 16, 16, 17, 17, 17, 17, 18, 19], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0  …  0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], 19, 19))

generating the matrix

flp.hess.H

19×19 SparseArrays.SparseMatrixCSC{Float64, Int64} with 103 stored entries:
⎡⠑⣤⡂⠡⣤⡂⠡⠌⠂⠀⎤
⎢⠌⡈⠻⣦⡈⠻⣦⠠⠁⠀⎥
⎢⠠⠻⣦⡈⠻⣦⡈⠁⠄⠀⎥
⎢⡁⠆⠈⡛⠆⠈⡛⢌⠀⠀⎥
⎣⠈⠀⠁⠀⠀⠁⠀⠀⠑⠄⎦

and for the Jacobian

flp.jac

A AutoDiff Jacobian for ExaPF.ComposedExpressions{ExaPF.PolarBasis{Vector{Int64}, SparseArrays.SparseMatrixCSC{Float64, Int64}}, ExaPF.MultiExpressions}(PolarBasis (AbstractExpression), ExaPF.MultiExpressions(ExaPF.AutoDiff.AbstractExpression[PowerFlowBalance (AbstractExpression), PowerGenerationBounds (AbstractExpression), LineFlows (AbstractExpression)]))
Number of Jacobian colors: 8

generating the matrix

flp.jac.J

36×19 SparseArrays.SparseMatrixCSC{Float64, Int64} with 176 stored entries:
⎡⠑⣤⡂⠡⣤⡂⠡⠌⠊⠂⎤
⎢⠌⡈⠻⣦⡈⠻⣦⠠⠁⠀⎥
⎢⠠⠻⣦⡈⠻⣦⡈⠁⠄⠀⎥
⎢⠁⡆⠈⠛⡆⠈⠛⡌⠀⠀⎥
⎢⠑⣄⠂⠁⣄⠂⠁⠌⠂⠀⎥
⎢⠐⠈⢧⡀⠈⢧⡀⠀⠂⠀⎥
⎢⠁⡄⠀⠳⡄⠀⠳⡈⠀⠀⎥
⎢⠠⠙⣆⠀⠙⣆⠀⠀⠄⠀⎥
⎣⠂⡀⠈⢧⡀⠈⢧⠐⠀⠀⎦

Both AD backends use coloring to reduce the number of Hessian-vector and Jacobian-vector products required to evaluate the Hessian and Jacobian in sparse format.

To avoid dealing explicitly with the AD backends, Argos provides a function to query directly the Hessian and Jacobian in COO format:

# Query sparsity pattern:
j_I, j_J = Argos.jacobian_structure(flp)
nnzj = length(j_I)
j_V = zeros(nnzj)
Argos.jacobian_coo!(flp, j_V, x)
sparse(j_I, j_J, j_V) # build a SparseMatrixCSC

36×19 SparseArrays.SparseMatrixCSC{Float64, Int64} with 176 stored entries:
⎡⠑⣤⡂⠡⣤⡂⠡⠌⠊⠂⎤
⎢⠌⡈⠻⣦⡈⠻⣦⠠⠁⠀⎥
⎢⠠⠻⣦⡈⠻⣦⡈⠁⠄⠀⎥
⎢⠁⡆⠈⠛⡆⠈⠛⡌⠀⠀⎥
⎢⠑⣄⠂⠁⣄⠂⠁⠌⠂⠀⎥
⎢⠐⠈⢧⡀⠈⢧⡀⠀⠂⠀⎥
⎢⠁⡄⠀⠳⡄⠀⠳⡈⠀⠀⎥
⎢⠠⠙⣆⠀⠙⣆⠀⠀⠄⠀⎥
⎣⠂⡀⠈⢧⡀⠈⢧⠐⠀⠀⎦

and for the Hessian:

# Query sparsity pattern:
h_I, h_J = Argos.hessian_structure(flp)
nnzh = length(h_I)
h_V = zeros(nnzh)
y = rand(m)
Argos.hessian_lagrangian_coo!(flp, h_V, x, y, 1.0)
sparse(h_I, h_J, h_V) # build a SparseMatrixCSC

19×19 SparseArrays.SparseMatrixCSC{Float64, Int64} with 61 stored entries:
⎡⠑⣄⠀⠀⠀⠀⠀⠀⠀⠀⎤
⎢⠌⡈⠳⣄⠀⠀⠀⠀⠀⠀⎥
⎢⠠⠻⣦⡈⠳⣄⠀⠀⠀⠀⎥
⎢⡁⠆⠈⡛⠆⠈⡓⢄⠀⠀⎥
⎣⠈⠀⠁⠀⠀⠁⠀⠀⠑⠄⎦

Info

For the Hessian, only the lower-triangular are being returned.

Deport on CUDA GPU

Deporting all the operations on a CUDA GPU simply amounts to instantiating a FullSpaceEvaluator`](@ref) on the GPU, with

using CUDAKernels # suppose CUDAKernels has been downloaded
flp = Argos.FullSpaceEvaluator(datafile; device=CUDADevice())

Then, the API remains exactly the same as on the CPU.

When using device=CUDADevice(), the model is entirely instantiated on the device, without data left on the host (hence minimizing the communication costs). The computation of the derivatives is streamlined by propagating the tangents in parallel, leading to faster evaluations of the callbacks. As expected, the larger the model, the more significant the performance gain.