Skip to content

Commit

Permalink
v0.5.2 (#40)
Browse files Browse the repository at this point in the history
* FIX: inconsistency between terminals / connections

* WIP: arbitrary phase pq inverter quick fix

* grid-following and grid-forming 1p&3p

* FIX: syntax in constraint_mc_pq_inverter

* ADD: JuMP v1 support

* v0.5.2 release prep

Co-authored-by: jtabarez <53569964+jtabarez@users.noreply.github.com>
  • Loading branch information
pseudocubic and jtabarez authored May 24, 2022
1 parent 2f65565 commit 884b55d
Show file tree
Hide file tree
Showing 10 changed files with 322 additions and 271 deletions.
6 changes: 6 additions & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,12 @@

- none

## v0.5.2

- Add JuMP v1 to compat
- Fixed for-loops over conductors in inverter constraints
- Added support for 1-phase solar devices in inverter constraints

## v0.5.1

- Removed duplicate `add_ct` function
Expand Down
4 changes: 2 additions & 2 deletions Project.toml
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "PowerModelsProtection"
uuid = "719c1aef-945b-435a-a240-4c2992e5e0df"
authors = ["Art Barnes", "Jose Tabarez"]
version = "0.5.1"
version = "0.5.2"

[deps]
InfrastructureModels = "2030c09a-7f63-5d83-885d-db604e0e9cc0"
Expand All @@ -15,7 +15,7 @@ Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"
[compat]
InfrastructureModels = "0.7.3"
Ipopt = "0.9, 1.0.2"
JuMP = "0.22, 0.23"
JuMP = "0.22, 0.23, 1"
Graphs = "1.4.1, 1.6"
PowerModels = "0.19.1"
PowerModelsDistribution = "~0.13.3, 0.14"
Expand Down
180 changes: 117 additions & 63 deletions src/core/constraint_inverter.jl
Original file line number Diff line number Diff line change
Expand Up @@ -249,7 +249,7 @@ end
Constraints for fault current contribution of multiconductor inverter in grid-following mode
"""
function constraint_mc_pq_inverter(pm::_PMD.AbstractUnbalancedIVRModel, nw::Int, i::Int, bus_id::Int, pg, qg, cmax)
function constraint_mc_pq_inverter(pm::_PMD.AbstractUnbalancedIVRModel, nw::Int, i::Int, bus_id::Int, pg, qg, qmax, qmin, imax, connections)
ar = -1/6
ai = sqrt(3)/6
a2r = -1/6
Expand All @@ -261,44 +261,98 @@ function constraint_mc_pq_inverter(pm::_PMD.AbstractUnbalancedIVRModel, nw::Int,
crg = _PMD.var(pm, nw, :crg, i)
cig = _PMD.var(pm, nw, :cig, i)

p_int = _PMD.var(pm, nw, :p_int, i)
q_int = _PMD.var(pm, nw, :q_int, i)
crg_pos= _PMD.var(pm, nw, :crg_pos, i)
cig_pos = _PMD.var(pm, nw, :cig_pos, i)
vrg_pos= _PMD.var(pm, nw, :vrg_pos, i)
vig_pos = _PMD.var(pm, nw, :vig_pos, i)
crg_pos_max = _PMD.var(pm, nw, :crg_pos_max, i)
cig_pos_max = _PMD.var(pm, nw, :cig_pos_max, i)
z = _PMD.var(pm, nw, :z_gfli, i)

cnds = _PMD.ref(pm, nw, :bus, bus_id)["terminals"]
ncnds = length(cnds)


# Zero-Sequence
JuMP.@constraint(pm.model, sum(crg[c] for c in cnds) == 0)
JuMP.@constraint(pm.model, sum(cig[c] for c in cnds) == 0)

# Negative-Sequence
JuMP.@constraint(pm.model, (1/3)*crg[1] + a2r*crg[2] - a2i*cig[2] + ar*crg[3] - ai*cig[3] == 0)
JuMP.@constraint(pm.model, (1/3)*cig[1] + a2r*cig[2] + a2i*crg[2] + ar*cig[3] + ai*crg[3] == 0)

# Positive-Sequence
JuMP.@constraint(pm.model, (1/3)*crg[1] + ar*crg[2] - ai*cig[2] + a2r*crg[3] - a2i*cig[3] == crg_pos)
JuMP.@constraint(pm.model, (1/3)*cig[1] + ar*cig[2] + ai*crg[2] + a2r*cig[3] + a2i*crg[3] == cig_pos)
JuMP.@constraint(pm.model, (1/3)*vr[1] + ar*vr[2] - ai*vi[2] + a2r*vr[3] - a2i*vi[3] == vrg_pos)
JuMP.@constraint(pm.model, (1/3)*vi[1] + ar*vi[2] + ai*vr[2] + a2r*vi[3] + a2i*vr[3] == vig_pos)

JuMP.@constraint(pm.model, 0.0 == crg_pos_max*cig_pos - cig_pos_max*crg_pos)
JuMP.@constraint(pm.model, crg_pos_max^2 + cig_pos_max^2 == cmax^2)
JuMP.@constraint(pm.model, crg_pos_max * crg_pos >= 0.0)
JuMP.@constraint(pm.model, cig_pos_max * cig_pos >= 0.0)
JuMP.@constraint(pm.model, crg_pos^2 + cig_pos^2 <= cmax^2)
JuMP.@NLconstraint(pm.model, (crg_pos^2 + cig_pos^2 - cmax^2)*z >= 0.0)
JuMP.@constraint(pm.model, p_int == vrg_pos*crg_pos + vig_pos*cig_pos)
JuMP.@constraint(pm.model, 0.0 == vig_pos*crg_pos - vrg_pos*cig_pos)
JuMP.@constraint(pm.model, p_int <= pg/3)
JuMP.@constraint(pm.model, p_int >= (1-z) * pg/3)
ncnds = length(connections)

m = 100

if ncnds == 3
crg_pos = _PMD.var(pm, nw)[:crg_pos] = JuMP.@variable(pm.model,
base_name = "$(nw)_crg_pos_$(i)",
start = 0.0,
)
cig_pos = _PMD.var(pm, nw)[:cig_pos] = JuMP.@variable(pm.model,
base_name = "$(nw)_cig_pos_$(i)",
start = 0.0,
)

vr_pos = _PMD.var(pm, nw)[:vr_pos] = JuMP.@variable(pm.model,
base_name = "$(nw)_vr_pos_$(i)",
start = 0.0,
)
vi_pos = _PMD.var(pm, nw)[:vi_pos] = JuMP.@variable(pm.model,
base_name = "$(nw)_vi_pos_$(i)",
start = 0.0,
)

z = _PMD.var(pm, nw)[:z_gfli] = JuMP.@variable(pm.model,
base_name = "$(nw)_z_gfli_$(i)",
lower_bound = 0.0,
upper_bound = 1.0,
start = 0.0
)

pg_int = _PMD.var(pm, nw)[:pg_int] = JuMP.@variable(pm.model,
base_name = "$(nw)_pg_int_$(i)",
lower_bound = 0.0,
upper_bound = pg[1],
start = 0.0,
)

qg_int = _PMD.var(pm, nw)[:qg_int] = JuMP.@variable(pm.model,
base_name = "$(nw)_qg_int_$(i)",
start = 0.0,
lower_bound = qmin[1],
upper_bound = qmax[1],
)

# Zero-Sequence
JuMP.@constraint(pm.model, sum(crg[c] for c in connections) == 0.0)
JuMP.@constraint(pm.model, sum(cig[c] for c in connections) == 0.0)

# Negative-Sequence
JuMP.@constraint(pm.model, (1/3)*crg[1] + a2r*crg[2] - a2i*cig[2] + ar*crg[3] - ai*cig[3] == 0.0)
JuMP.@constraint(pm.model, (1/3)*cig[1] + a2r*cig[2] + a2i*crg[2] + ar*cig[3] + ai*crg[3] == 0.0)

# Positive-Sequence
JuMP.@constraint(pm.model, (1/3)*crg[1] + ar*crg[2] - ai*cig[2] + a2r*crg[3] - a2i*cig[3] == crg_pos)
JuMP.@constraint(pm.model, (1/3)*cig[1] + ar*cig[2] + ai*crg[2] + a2r*cig[3] + a2i*crg[3] == cig_pos)
JuMP.@constraint(pm.model, (1/3)*vr[1] + ar*vr[2] - ai*vi[2] + a2r*vr[3] - a2i*vi[3] == vr_pos)
JuMP.@constraint(pm.model, (1/3)*vi[1] + ar*vi[2] + ai*vr[2] + a2r*vi[3] + a2i*vr[3] == vi_pos)

JuMP.@NLconstraint(pm.model, imax[1]^2 - crg_pos^2 - cig_pos^2 >= 0.0)
JuMP.@NLconstraint(pm.model, (imax[1]^2 - crg_pos^2 - cig_pos^2)*z <= 0.0)
JuMP.@NLconstraint(pm.model, vr_pos*crg_pos + vi_pos*cig_pos == pg[1] - pg_int*z)
JuMP.@NLconstraint(pm.model, vi_pos*crg_pos - vr_pos*cig_pos == 0.0)

else
z = _PMD.var(pm, nw)[:z_gfli] = JuMP.@variable(pm.model,
[c in connections], base_name = "$(nw)_z_gfli_$(i)",
lower_bound = 0.0,
upper_bound = 1.0,
start = 0.0
)
pg_int = _PMD.var(pm, nw)[:pg_int] = JuMP.@variable(pm.model,
[c in connections], base_name = "$(nw)_pg_int_$(i)",
lower_bound = 0.0,
upper_bound = pg[1],
start = 0.0,
)

qg_int = _PMD.var(pm, nw)[:qg_int] = JuMP.@variable(pm.model,
[c in connections], base_name = "$(nw)_qg_int_$(i)",
start = 0.0,
lower_bound = qmin[1],
upper_bound = qmax[1],
)

for (idx,c) in enumerate(connections)
JuMP.@NLconstraint(pm.model, imax[idx]^2 - crg[c]^2 - cig[c]^2 >= 0.0)
JuMP.@NLconstraint(pm.model, m*(imax[idx]^2 - crg[c]^2 - cig[c]^2)*z[c] <= 0.0)
JuMP.@NLconstraint(pm.model, vr[c]*crg[c] + vi[c]*cig[c] == pg[idx]-z[c]*pg_int[c])
JuMP.@NLconstraint(pm.model, vi[c]*crg[c] - vr[c]*cig[c] == qg[idx]-z[c]*qg_int[c])
end
end
end


Expand All @@ -321,17 +375,17 @@ function constraint_mc_grid_forming_inverter(pm::_PMD.AbstractUnbalancedIVRModel
cnds = _PMD.ref(pm, n, :bus, bus_id, "terminals")
ncnds = length(cnds)

vm = [vrstar[c]^2 + vistar[c]^2 for c in 1:ncnds]
vm = [vrstar[idx]^2 + vistar[idx]^2 for (idx,c) in enumerate(cnds)]

for c in 1:ncnds
for (idx,c) in cnds
# current limits
# z = 1 -> invert in in current limiting mode
JuMP.@constraint(pm.model, crg[c]^2 + cig[c]^2 <= cmax^2)
JuMP.@NLconstraint(pm.model, (crg[c]^2 + cig[c]^2 - cmax^2)*z[c] >= 0.0)

# terminal voltage mag
JuMP.@constraint(pm.model, vr[c]^2 + vi[c]^2 <= vm[c] * (1+z[c]))
JuMP.@constraint(pm.model, vr[c]^2 + vi[c]^2 >= vm[c] * (1-z[c]))
JuMP.@constraint(pm.model, vr[c]^2 + vi[c]^2 <= vm[idx] * (1+z[c]))
JuMP.@constraint(pm.model, vr[c]^2 + vi[c]^2 >= vm[idx] * (1-z[c]))

# terminal voltage phase
JuMP.@constraint(pm.model, 0.0 == vr[c]*vistar[c] - vi[c]*vrstar[c])
Expand Down Expand Up @@ -371,21 +425,21 @@ function constraint_mc_grid_formimg_inverter_impedance(pm::_PMD.AbstractUnbalanc
cnds = _PMD.ref(pm, nw, :bus, bus_id, "terminals")
ncnds = length(cnds)

vm = [vr0[c]^2 + vi0[c]^2 for c in 1:ncnds]
vm = [vr0[idx]^2 + vi0[idx]^2 for (idx,c) in enumerate(cnds)]

for c in 1:ncnds
for (idx,c) in enumerate(cnds)
# current limits
JuMP.@constraint(pm.model, crg[c]^2 + cig[c]^2 <= cmax^2)
JuMP.@NLconstraint(pm.model, (crg[c]^2 + cig[c]^2 - cmax^2)*z[c] >= 0.0)
JuMP.@NLconstraint(pm.model, crg[c]^2 + cig[c]^2 <= cmax[c]^2)
JuMP.@NLconstraint(pm.model, m*(cmax[c]^2 - crg[c]^2 - cig[c]^2)*z[c] <= 0.0)

# setpoint voltage mag
JuMP.@constraint(pm.model, vrstar[c]^2 + vistar[c]^2 <= vm[c] * (1+z[c]))
JuMP.@constraint(pm.model, vrstar[c]^2 + vistar[c]^2 >= vm[c] * (1-z[c]))
JuMP.@constraint(pm.model, vrstar[c]^2 + vistar[c]^2 <= vm[idx] * (1+z[c]))
JuMP.@constraint(pm.model, vrstar[c]^2 + vistar[c]^2 >= vm[idx] * (1-z[c]))

# setpoint voltage phase
JuMP.@constraint(pm.model, 0.0 == vrstar[c]*vi0[c] - vistar[c]*vr0[c])
JuMP.@constraint(pm.model, vrstar[c] * vr0[c] >= 0.0)
JuMP.@constraint(pm.model, vistar[c] * vi0[c] >= 0.0)
JuMP.@constraint(pm.model, 0.0 == vrstar[c]*vi0[idx] - vistar[c]*vr0[idx])
JuMP.@constraint(pm.model, vrstar[c] * vr0[idx] >= 0.0)
JuMP.@constraint(pm.model, vistar[c] * vi0[idx] >= 0.0)

# terminal voltage setpoint based virtual impedance
JuMP.@constraint(pm.model, vr[c] == vrstar[c] - r[c]*crg[c] + x[c]*cig[c])
Expand All @@ -405,7 +459,7 @@ end
Constraints for fault current contribution of multiconductor inverter in grid-forming mode with power matching
"""
function constraint_mc_grid_formimg_inverter_virtual_impedance(pm::_PMD.AbstractUnbalancedIVRModel, nw::Int, i::Int, bus_id::Int, vr0, vi0, pmax, cmax, smax, ang, terminals)
function constraint_mc_grid_formimg_inverter_virtual_impedance(pm::_PMD.AbstractUnbalancedIVRModel, nw::Int, i::Int, bus_id::Int, vr0, vi0, pmax, imax, smax, ang, connections)
vr = _PMD.var(pm, nw, :vr, bus_id)
vi = _PMD.var(pm, nw, :vi, bus_id)

Expand All @@ -424,20 +478,20 @@ function constraint_mc_grid_formimg_inverter_virtual_impedance(pm::_PMD.Abstract
p = _PMD.var(pm, nw, :p_solar, i)
q = _PMD.var(pm, nw, :q_solar, i)

vm = [vr0[c]^2 + vi0[c]^2 for c in terminals]
vm = [vr0[idx]^2 + vi0[idx]^2 for (idx,c) in enumerate(connections)]

for c in terminals
JuMP.@constraint(pm.model, crg[c]^2 + cig[c]^2 <= cmax^2)
JuMP.@NLconstraint(pm.model, (crg[c]^2 + cig[c]^2 - cmax^2)*z[c] >= 0.0)
for (idx,c) in enumerate(connections)
JuMP.@constraint(pm.model, crg[c]^2 + cig[c]^2 <= imax[c]^2)
JuMP.@NLconstraint(pm.model, (crg[c]^2 + cig[c]^2 - imax[c]^2)*z[c] >= 0.0)

JuMP.@constraint(pm.model, vrstar[c]^2 + vistar[c]^2 >= vm[c] * (1-z[c]))
JuMP.@constraint(pm.model, vrstar[c]^2 + vistar[c]^2 <= vm[c])
JuMP.@constraint(pm.model, vrstar[c]^2 + vistar[c]^2 >= vm[idx] * (1-z[c]))
JuMP.@constraint(pm.model, vrstar[c]^2 + vistar[c]^2 <= vm[idx])

if ang
# setpoint voltage phase
JuMP.@constraint(pm.model, 0.0 == vrstar[c]*vi0[c] - vistar[c]*vr0[c])
JuMP.@constraint(pm.model, vrstar[c] * vr0[c] >= 0.0)
JuMP.@constraint(pm.model, vistar[c] * vi0[c] >= 0.0)
JuMP.@constraint(pm.model, 0.0 == vrstar[c]*vi0[idx] - vistar[c]*vr0[idx])
JuMP.@constraint(pm.model, vrstar[c] * vr0[idx] >= 0.0)
JuMP.@constraint(pm.model, vistar[c] * vi0[idx] >= 0.0)
end

JuMP.@constraint(pm.model, rv[c] <= .10)
Expand All @@ -450,8 +504,8 @@ function constraint_mc_grid_formimg_inverter_virtual_impedance(pm::_PMD.Abstract
end

# DC-link power
JuMP.@constraint(pm.model, sum(vr[c]*crg[c] + vi[c]*cig[c] for c in terminals) == p)
JuMP.@constraint(pm.model, sum(vi[c]*crg[c] - vr[c]*cig[c] for c in terminals) == q)
JuMP.@constraint(pm.model, sum(vr[c]*crg[c] + vi[c]*cig[c] for c in connections) == p)
JuMP.@constraint(pm.model, sum(vi[c]*crg[c] - vr[c]*cig[c] for c in connections) == q)

JuMP.@constraint(pm.model, p^2 + q^2 <= smax^2)
JuMP.@constraint(pm.model, p <= pmax)
Expand Down
Loading

2 comments on commit 884b55d

@pseudocubic
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@JuliaRegistrator
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Registration pull request created: JuliaRegistries/General/60939

After the above pull request is merged, it is recommended that a tag is created on this repository for the registered package version.

This will be done automatically if the Julia TagBot GitHub Action is installed, or can be done manually through the github interface, or via:

git tag -a v0.5.2 -m "<description of version>" 884b55d39556f035d09235e5c315daf798c397ff
git push origin v0.5.2

Please sign in to comment.