Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 2 additions & 1 deletion CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,8 @@ PowerModels.jl Change Log
=========================

### Staged
- Improved support in `build_pf` for slack buses with multiple generators
- Added Jacobian calculation of basic network data (#762)
- Improved support in `build_pf` for slack buses with multiple generators (#785)
- Fixed bug in `compute_ac_pf` when reporting slack bus `va` values
- Fixed bug in `constraint_ohms_yt_from_on_off` for DCP models on branches with negative reactance values

Expand Down
2 changes: 2 additions & 0 deletions docs/src/basic-data-utilities.md
Original file line number Diff line number Diff line change
Expand Up @@ -48,6 +48,7 @@ calc_basic_incidence_matrix
calc_basic_admittance_matrix
calc_basic_susceptance_matrix
calc_basic_branch_susceptance_matrix
calc_basic_jacobian_matrix
calc_basic_ptdf_matrix
calc_basic_ptdf_row
```
Expand Down Expand Up @@ -104,6 +105,7 @@ basic network data using Julia's native linear equation solver,
compute_basic_dc_pf
```


!!! tip
By default PowerModels uses Julia's SparseArrays to ensure the best
performance of matrix operations on large power network datasets.
Expand Down
71 changes: 1 addition & 70 deletions src/core/data_basic.jl
Original file line number Diff line number Diff line change
Expand Up @@ -385,76 +385,6 @@ function calc_basic_ptdf_row(data::Dict{String,<:Any}, branch_index::Int)
return ptdf_column
end

"""
Computes the jacobian submatrices
H = dP/dδ
L = dQ/dV

The derivatives are ordered according to the bus number,
H[1, 1] is ∂P_1 / ∂δ_1
H[1, 2] is ∂P_1 / ∂δ_2

It does not exclude PV buses rows and columns
"""
function calc_basic_decoupled_jacobian_matrices(data::Dict{String,<:Any})
if !get(data, "basic_network", false)
Memento.warn(_LOGGER, "calc_basic_decoupled_jacobian_matrices requires basic network data and given data may be incompatible. make_basic_network can be used to transform data into the appropriate form.")
end
num_bus = length(data["bus"])
v = calc_basic_bus_voltage(data)
vm, va = abs.(v), angle.(v)
Y = calc_basic_admittance_matrix(data)
neighbors = [Set{Int}([i]) for i in 1:num_bus]
I, J, V = findnz(Y)
for nz in eachindex(V)
push!(neighbors[I[nz]], J[nz])
push!(neighbors[J[nz]], I[nz])
end
H0_I = Int64[]; H0_J = Int64[]; H0_V = Float64[];
L0_I = Int64[]; L0_J = Int64[]; L0_V = Float64[];
for i in 1:num_bus
for j in neighbors[i]
push!(H0_I, i); push!(H0_J, j); push!(H0_V, 0.0)
push!(L0_I, i); push!(L0_J, j); push!(L0_V, 0.0)
end
end
H = sparse(H0_I, H0_J, H0_V)
L = sparse(L0_I, L0_J, L0_V)
for i in 1:num_bus
for j in neighbors[i]
bus_type = data["bus"]["$(j)"]["bus_type"]
if bus_type == 1
if i == j
y_ii = Y[i,i]
H[i, j] = + vm[i] * sum( -real(Y[i,k]) * vm[k] * sin(va[i] - va[k]) + imag(Y[i,k]) * vm[k] * cos(va[i] - va[k]) for k in neighbors[i] if k != i )
L[i, j] = - 2*imag(y_ii)*vm[i] + sum( real(Y[i,k]) * vm[k] * sin(va[i] - va[k]) - imag(Y[i,k]) * vm[k] * cos(va[i] - va[k]) for k in neighbors[i] if k != i )
else
y_ij = Y[i,j]
H[i, j] = vm[i] * vm[j] * ( real(y_ij) * sin(va[i] - va[j]) - imag(y_ij) * cos(va[i] - va[j]) )
L[i, j] = vm[i] * ( real(y_ij) * sin(va[i] - va[j]) - imag(y_ij) * cos(va[i] - va[j]) )
end
elseif bus_type == 2
if i == j
y_ii = Y[i,i]
H[i, j] = vm[i] * sum( -real(Y[i,k]) * vm[k] * sin(va[i] - va[k]) + imag(Y[i,k]) * vm[k] * cos(va[i] - va[k]) for k in neighbors[i] if k != i )
L[i, j] = 1.0
else
y_ij = Y[i,j]
H[i, j] = vm[i] * vm[j] * ( real(y_ij) * sin(va[i] - va[j]) - imag(y_ij) * cos(va[i] - va[j]) )
L[i, j] = 0.0
end
elseif bus_type == 3
if i == j
H[i, j] = 1.0
L[i, j] = 1.0
end
else
@assert false
end
end
end
return H, L
end

"""
given a basic network data dict, returns a sparse real valued Jacobian matrix
Expand All @@ -465,6 +395,7 @@ function calc_basic_jacobian_matrix(data::Dict{String,<:Any})
if !get(data, "basic_network", false)
Memento.warn(_LOGGER, "calc_basic_jacobian_matrix requires basic network data and given data may be incompatible. make_basic_network can be used to transform data into the appropriate form.")
end

num_bus = length(data["bus"])
v = calc_basic_bus_voltage(data)
vm, va = abs.(v), angle.(v)
Expand Down
15 changes: 5 additions & 10 deletions test/common.jl
Original file line number Diff line number Diff line change
Expand Up @@ -72,7 +72,7 @@ end
"""
An AC Power Flow Solver from scratch.
"""
function compute_basic_ac_pf!(data::Dict{String, Any}; decoupled=false)
function compute_basic_ac_pf!(data::Dict{String, Any})
if !get(data, "basic_network", false)
Memento.warn(_LOGGER, "compute_basic_ac_pf requires basic network data and given data may be incompatible. make_basic_network can be used to transform data into the appropriate form.")
end
Expand Down Expand Up @@ -104,16 +104,11 @@ function compute_basic_ac_pf!(data::Dict{String, Any}; decoupled=false)
if LinearAlgebra.normInf([delta_P; delta_Q]) < tol
break
end

# STEP 2 and 3: Compute the jacobian and update step
if !decoupled
J = calc_basic_jacobian_matrix(data)
x = J \ [delta_P; delta_Q]
else
H, L = calc_basic_decoupled_jacobian_matrices(data)
va = H \ delta_P
vm = L \ delta_Q
x = [va; vm]
end
J = calc_basic_jacobian_matrix(data)
x = J \ [delta_P; delta_Q]

# STEP 4
# update voltage variables
for i in 1:bus_num
Expand Down
30 changes: 2 additions & 28 deletions test/data-basic.jl
Original file line number Diff line number Diff line change
Expand Up @@ -195,7 +195,7 @@ end
for (i, bus) in data["bus"] # All buses PQ
bus["bus_type"] = 1
end
J = PowerModels.calc_basic_jacobian_matrix(data);
J = calc_basic_jacobian_matrix(data)

num_bus = length(data["bus"])

Expand All @@ -212,7 +212,7 @@ end

@testset "24-bus-case" begin
data = make_basic_network(PowerModels.parse_file("../test/data/matpower/case24.m"))
J = PowerModels.calc_basic_jacobian_matrix(data)
J = calc_basic_jacobian_matrix(data)

num_bus = length(data["bus"])

Expand Down Expand Up @@ -244,15 +244,6 @@ end
@test isapprox(J[3, num_bus+9], -2.1624; atol=1e-4) # dP/dvm non-diagonal
@test isapprox(J[num_bus+3, num_bus+9], -7.9603; atol=1e-4) # dQ/dvm non-diagonal
end

@testset "30-bus-case-decoupled-matrices" begin
data = make_basic_network(PowerModels.parse_file("../test/data/matpower/case24.m"))
H, L = calc_basic_decoupled_jacobian_matrices(data)
J = calc_basic_jacobian_matrix(data)
n = length(data["bus"])
@test isapprox(J[1:n, 1:n], H; atol=1e-6)
@test isapprox(J[n+1:end, n+1:end], L; atol=1e-6)
end
end

@testset "basic ac power flow" begin
Expand All @@ -273,23 +264,6 @@ end
end
end

@testset "9-bus-case" begin
data = make_basic_network(PowerModels.parse_file("../test/data/matpower/case9.m"))
solution = compute_ac_pf(data)["solution"]

compute_basic_ac_pf!(data; decoupled=true)

for (i, bus) in data["bus"]
@test isapprox(solution["bus"][i]["va"], bus["va"]; atol=1e-4)
@test isapprox(solution["bus"][i]["vm"], bus["vm"]; atol=1e-4)
end

for (i, gen) in data["gen"]
@test isapprox(solution["gen"][i]["pg"], gen["pg"]; atol=1e-4)
@test isapprox(solution["gen"][i]["qg"], gen["qg"]; atol=1e-4)
end
end

@testset "14-bus-case" begin
data = make_basic_network(PowerModels.parse_file("../test/data/matpower/case14.m"))
solution = compute_ac_pf(data)["solution"]
Expand Down