Skip to content

Commit

Permalink
Adding SOC extended DistFlow (BFM) formulation (#229)
Browse files Browse the repository at this point in the history
* Adding SOC BIM / extended DistFlow formulation
* Updated current variable name + dev documentation for variable naming
* defining opf_bf and pf_bf problem type
* add arc_from to constraint_voltage_angle_difference
* change signature of constraint_power_losses
* updated docstrings
* Derivation of distflow formulation with asymmetric shunts
* moved/removed functions + derive bound for series current
* Complete DF derivation with asymmetric shunts, including conductance
* compatible with asymmetric pi-section
  • Loading branch information
frederikgeth authored and ccoffrin committed Apr 6, 2018
1 parent b33773e commit 57c3299
Show file tree
Hide file tree
Showing 19 changed files with 677 additions and 39 deletions.
131 changes: 130 additions & 1 deletion docs/src/developer.md
Original file line number Diff line number Diff line change
@@ -1,3 +1,132 @@
# Developer Documentation

Nothing yet.
## Variable and parameter naming scheme

### Suffixes

- `_fr`: from-side ('i'-node)
- `_to`: to-side ('j'-node)

### Power

Defining power $s = p + j \cdot q$ and $sm = |s|$
- `s`: complex power (VA)
- `sm`: apparent power (VA)
- `p`: active power (W)
- `q`: reactive power (var)

### Voltage

Defining voltage $v = vm \angle va = vr + j \cdot vi$:
- `vm`: magnitude of (complex) voltage (V)
- `va`: angle of complex voltage (rad)
- `vr`: real part of (complex) voltage (V)
- `vi`: imaginary part of complex voltage (V)

### Current

Defining current $c = cm \angle ca = cr + j \cdot ci$:
- `cm`: magnitude of (complex) current (A)
- `ca`: angle of complex current (rad)
- `cr`: real part of (complex) current (A)
- `ci`: imaginary part of complex current (A)

### Voltage products

Defining voltage product $w = v_i \cdot v_j$ then
$w = wm \angle wa = wr + j\cdot wi$:
- `wm` (short for vvm): magnitude of (complex) voltage products (V$^2$)
- `wa` (short for vva): angle of complex voltage products (rad)
- `wr` (short for vvr): real part of (complex) voltage products (V$^2$)
- `wi` (short for vvi): imaginary part of complex voltage products (V$^2$)

### Current products

Defining current product $cc = c_i \cdot c_j$ then
$cc = ccm \angle cca = ccr + j\cdot cci$:
- `ccm`: magnitude of (complex) current products (A$^2$)
- `cca`: angle of complex current products (rad)
- `ccr`: real part of (complex) current products (A$^2$)
- `cci`: imaginary part of complex current products (A$^2$)

### Transformer ratio

Defining complex transformer ratio
$t = tm \angle ta = tr + j\cdot ti$:
- `tm`: magnitude of (complex) transformer ratio (-)
- `ta`: angle of complex transformer ratio (rad)
- `tr`: real part of (complex) transformer ratio (-)
- `ti`: imaginary part of complex transformer ratio (-)

### Impedance

Defining impedance
$z = r + j\cdot x$:
- `r`: resistance ($\Omega$)
- `x`: reactance ($\Omega$)

### Admittance

Defining admittance
$y = g + j\cdot b$:
- `g`: conductance ($S$)
- `b`: susceptance ($S$)


## DistFlow derivation

### For an asymmetric pi section
Following notation of [^1], but recognizing it derives the SOC BFM without shunts. In a pi-section, part of the total current $ I_{lij}$ at the from side flows through the series impedance, $I ^{s}_{lij}$, part of it flows through the from side shunt admittance $ I^{sh}_{lij}$. Vice versa for the to-side. Indicated by superscripts 's' (series) and 'sh' (shunt).
- Ohm's law: $U^{mag}_{j} \angle \theta_{j} = U^{mag}_{i}\angle \theta_{i} - z^{s}_{lij} \cdot I^{s}_{lij}$ $\forall lij$
- KCL at shunts: $ I_{lij} = I^{s}_{lij} + I^{sh}_{lij}$, $ I_{lji} = I^{s}_{lji} + I^{sh}_{lji} $
- Observing: $I^{s}_{lij} = - I^{s}_{lji}$, $ \vert I^{s}_{lij} \vert = \vert I^{s}_{lji} \vert $
- Ohm's law times its own complex conjugate: $(U^{mag}_{j})^2 = (U^{mag}_{i}\angle \theta_{i} - z^{s}_{lij} \cdot I^{s}_{lij})\cdot (U^{mag}_{i}\angle \theta_{i} - z^{s}_{lij} \cdot I^{s}_{lij})^*$
- Defining $S^{s}_{lij} = P^{s}_{lij} + j\cdot Q^{s}_{lij} = (U^{mag}_{i}\angle \theta_{i}) \cdot (I^{s}_{lij})^*$
- Working it out $(U^{mag}_{j})^2 = (U^{mag}_{i})^2 - 2 \cdot(r^{s}_{lij} \cdot P^{s}_{lij} + x^{s}_{lij} \cdot Q^{s}_{lij}) $ + $((r^{s}_{lij})^2 + (x^{s}_{lij})^2)\vert I^{s}_{lij} \vert^2$

Power flow balance w.r.t. branch *total* losses
- Active power flow: $P_{lij}$ + $ P_{lji} $ = $ g^{sh}_{lij} \cdot (U^{mag}_{i})^2 + r^{s}_{l} \cdot \vert I^{s}_{lij} \vert^2 + g^{sh}_{lji} \cdot (U^{mag}_{j})^2 $
- Reactive power flow: $Q_{lij}$ + $ Q_{lji} $ = $ -b^{sh}_{lij} \cdot (U^{mag}_{i})^2 + x^{s}_{l} \cdot \vert I^{s}_{lij} \vert^2 - b^{sh}_{lji} \cdot (U^{mag}_{j})^2 $
- Current definition: $ \vert S^{s}_{lij} \vert^2 $ $=(U^{mag}_{i})^2 \cdot \vert I^{s}_{lij} \vert^2 $

Substitution:
- Voltage from: $(U^{mag}_{i})^2 \rightarrow w_{i}$
- Voltage to: $(U^{mag}_{j})^2 \rightarrow w_{j}$
- Series current : $\vert I^{s}_{lij} \vert^2 \rightarrow l^{s}_{l}$
Note that $l^{s}_{l}$ represents squared magnitude of the *series* current, i.e. the current flow through the series impedance in the pi-model.

Power flow balance w.r.t. branch *total* losses
- Active power flow: $P_{lij}$ + $ P_{lji} $ = $ g^{sh}_{lij} \cdot w_{i} + r^{s}_{l} \cdot l^{s}_{l} + g^{sh}_{lji} \cdot w_{j} $
- Reactive power flow: $Q_{lij}$ + $ Q_{lji} $ = $ -b^{sh}_{lij} \cdot w_{i} + x^{s}_{l} \cdot l^{s}_{l} - b^{sh}_{lji} \cdot w_{j} $

Power flow balance w.r.t. branch *series* losses:
- Series active power flow : $P^{s}_{lij} + P^{s}_{lji}$ $ = r^{s}_{l} \cdot l^{s}_{l} $
- Series reactive power flow: $Q^{s}_{lij} + Q^{s}_{lji}$ $ = x^{s}_{l} \cdot l^{s}_{l} $

Valid equality to link $w_{i}, l_{lij}, P^{s}_{lij}, Q^{s}_{lij}$:
- Nonconvex current definition: $(P^{s}_{lij})^2$ + $(Q^{s}_{lij})^2$ $=w_{i} \cdot l_{lij} $
- SOC current definition: $(P^{s}_{lij})^2$ + $(Q^{s}_{lij})^2$ $\leq$ $ w_{i} \cdot l_{lij} $


### Adding an ideal transformer
Adding an ideal transformer at the from side implicitly creates an internal branch voltage, between the transformer and the pi-section.
- new voltage: $w^{'}_{l}$
- ideal voltage magnitude transformer: $w^{'}_{l} = \frac{w_{i}}{(t^{mag})^2}$

W.r.t to the pi-section only formulation, we effectively perform the following substitution in all the equations above:
- $ w_{i} \rightarrow \frac{w_{i}}{(t^{mag})^2}$

The branch's power balance isn't otherwise impacted by adding the ideal transformer, as such transformer is lossless.

### Adding total current limits
- Total current from: $ \vert I_{lij} \vert \leq I^{rated}_{l}$
- Total current to: $ \vert I_{lji} \vert \leq I^{rated}_{l}$

In squared voltage magnitude variables:
- Total current from: $ (P_{lij})^2$ + $(Q_{lij})^2 \leq (I^{rated}_{l})^2 \cdot w_{i}$
- Total current to: $ (P_{lji})^2$ + $(Q_{lji})^2 \leq (I^{rated}_{l})^2 \cdot w_{j}$




[^1] Gan, L., Li, N., Topcu, U., & Low, S. (2012). Branch flow model for radial networks: convex relaxation. 51st IEEE Conference on Decision and Control, 1–8. Retrieved from http://smart.caltech.edu/papers/ExactRelaxation.pdf
5 changes: 5 additions & 0 deletions docs/src/formulations.md
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@ We begin with the top of the hierarchy, where we can distinguish between AC and
AbstractACPForm <: AbstractPowerFormulation
AbstractDCPForm <: AbstractPowerFormulation
AbstractWRForm <: AbstractPowerFormulation
AbstractWForm <: AbstractPowerFormulation
```

From there, different forms for ACP and DCP are possible:
Expand All @@ -17,6 +18,8 @@ StandardDCPForm <: AbstractDCPForm

SOCWRForm <: AbstractWRForm
QCWRForm <: AbstractWRForm

SOCDFForm <: AbstractWForm
```

## Power Models
Expand All @@ -29,6 +32,8 @@ DCPPowerModel = GenericPowerModel{StandardDCPForm}

SOCWRPowerModel = GenericPowerModel{SOCWRForm}
QCWRPowerModel = GenericPowerModel{QCWRForm}

SOCDFPowerModel = GenericPowerModel{SOCDFForm}
```

For details on `GenericPowerModel`, see the section on [Power Model](@ref).
Expand Down
86 changes: 86 additions & 0 deletions docs/src/specifications.md
Original file line number Diff line number Diff line change
Expand Up @@ -37,6 +37,45 @@ for (i,dcline) in pm.ref[:dcline]
end
```

## Optimal Power Flow (OPF) using the Branch Flow Model

### Objective
```julia
objective_min_fuel_cost(pm)
```

### Variables
```julia
variable_voltage(pm)
variable_active_generation(pm)
variable_reactive_generation(pm)
variable_branch_flow(pm)
variable_branch_current(pm)
variable_dcline_flow(pm)
```

### Constraints
```julia
constraint_theta_ref(pm)
constraint_voltage(pm)
for (i,bus) in pm.ref[:bus]
constraint_kcl_shunt(pm, bus)
end
for (i,branch) in pm.ref[:branch]
constraint_flow_losses(pm, branch)
constraint_voltage_magnitude_difference(pm, branch)
constraint_branch_current(pm, branch)

constraint_voltage_angle_difference(pm, branch)

constraint_thermal_limit_from(pm, branch)
constraint_thermal_limit_to(pm, branch)
end
for (i,dcline) in pm.ref[:dcline]
constraint_dcline(pm, dcline)
end
```

## Optimal Transmission Switching (OTS)

### General Assumptions
Expand Down Expand Up @@ -127,6 +166,53 @@ for (i,dcline) in pm.ref[:dcline]
end
```


## Power Flow (PF) using the Branch Flow Model

### Assumptions

### Variables
```julia
variable_voltage(pm, bounded = false)
variable_active_generation(pm, bounded = false)
variable_reactive_generation(pm, bounded = false)
variable_branch_flow(pm, bounded = false)
constraint_branch_current(pm, bounded = false)
variable_branch_current(pm, bounded = false)
```

### Constraints
```julia
constraint_theta_ref(pm)
constraint_voltage_magnitude_setpoint(pm, pm.ref[:bus][pm.ref[:ref_bus]])
constraint_voltage(pm)


for (i,bus) in pm.ref[:bus]
constraint_kcl_shunt(pm, bus)

# PV Bus Constraints
if length(pm.ref[:bus_gens][i]) > 0 && i != pm.ref[:ref_bus]
# this assumes inactive generators are filtered out of bus_gens
@assert bus["bus_type"] == 2

constraint_voltage_magnitude_setpoint(pm, bus)
for j in pm.ref[:bus_gens][i]
constraint_active_gen_setpoint(pm, pm.ref[:gen][j])
end
end
end

for (i,branch) in pm.ref[:branch]
constraint_flow_losses(pm, branch)
constraint_voltage_magnitude_difference(pm, branch)
constraint_branch_current(pm, branch)
end
for (i,dcline) in pm.ref[:dcline]
constraint_active_dcline_setpoint(pm, dcline)
end
```

## Transmission Network Expansion Planning (TNEP)

### Objective
Expand Down
3 changes: 3 additions & 0 deletions src/PowerModels.jl
Original file line number Diff line number Diff line change
Expand Up @@ -36,12 +36,15 @@ include("form/acp.jl")
include("form/acr.jl")
include("form/act.jl")
include("form/dcp.jl")
include("form/df.jl")
include("form/wr.jl")
include("form/wrm.jl")
include("form/shared.jl")

include("prob/pf.jl")
include("prob/pf_bf.jl")
include("prob/opf.jl")
include("prob/opf_bf.jl")
include("prob/ots.jl")
include("prob/tnep.jl")

Expand Down
4 changes: 4 additions & 0 deletions src/core/constraint.jl
Original file line number Diff line number Diff line change
Expand Up @@ -101,3 +101,7 @@ function constraint_active_dcline_setpoint(pm::GenericPowerModel, n::Int, f_idx,
@constraint(pm.model, p_fr == pf)
@constraint(pm.model, p_to == pt)
end

"do nothing, this model does not have complex voltage constraints"
function constraint_voltage(pm::GenericPowerModel, n::Int)
end
53 changes: 52 additions & 1 deletion src/core/constraint_template.jl
Original file line number Diff line number Diff line change
Expand Up @@ -471,11 +471,12 @@ function constraint_voltage_angle_difference(pm::GenericPowerModel, n::Int, i::I
branch = ref(pm, n, :branch, i)
f_bus = branch["f_bus"]
t_bus = branch["t_bus"]
arc_from = (i, f_bus, t_bus)
pair = (f_bus, t_bus)
buspair = ref(pm, n, :buspairs, pair)

if buspair["branch"] == i
constraint_voltage_angle_difference(pm, n, f_bus, t_bus, buspair["angmin"], buspair["angmax"])
constraint_voltage_angle_difference(pm, n, arc_from, f_bus, t_bus, buspair["angmin"], buspair["angmax"])
end
end
constraint_voltage_angle_difference(pm::GenericPowerModel, i::Int) = constraint_voltage_angle_difference(pm, pm.cnw, i)
Expand Down Expand Up @@ -530,3 +531,53 @@ function constraint_loss_lb(pm::GenericPowerModel, n::Int, i::Int)
constraint_loss_lb(pm, n, f_bus, t_bus, f_idx, t_idx, g_fr, b_fr, g_to, b_to, tr)
end
constraint_loss_lb(pm::GenericPowerModel, i::Int) = constraint_loss_lb(pm, pm.cnw, i)

function constraint_flow_losses(pm::GenericPowerModel, n::Int, i)
branch = ref(pm, n, :branch, i)
f_bus = branch["f_bus"]
t_bus = branch["t_bus"]
f_idx = (i, f_bus, t_bus)
t_idx = (i, t_bus, f_bus)
r = branch["br_r"]
x = branch["br_x"]
tm = branch["tap"]

g_sh_fr = branch["g_fr"]
g_sh_to = branch["g_to"]
b_sh_fr = branch["b_fr"]
b_sh_to = branch["b_to"]
constraint_flow_losses(pm::GenericPowerModel, n::Int, i, f_bus, t_bus, f_idx, t_idx, r, x, g_sh_fr, g_sh_to, b_sh_fr, b_sh_to, tm)
end
constraint_flow_losses(pm::GenericPowerModel, i::Int) = constraint_flow_losses(pm, pm.cnw, i)

function constraint_voltage_magnitude_difference(pm::GenericPowerModel, n::Int, i)
branch = ref(pm, n, :branch, i)
f_bus = branch["f_bus"]
t_bus = branch["t_bus"]
f_idx = (i, f_bus, t_bus)
t_idx = (i, t_bus, f_bus)

r = branch["br_r"]
x = branch["br_x"]
g_sh_fr = branch["g_fr"]
b_sh_fr = branch["b_fr"]
tm = branch["tap"]

constraint_voltage_magnitude_difference(pm::GenericPowerModel, n::Int, i, f_bus, t_bus, f_idx, t_idx, r, x, g_sh_fr, b_sh_fr, tm)

end
constraint_voltage_magnitude_difference(pm::GenericPowerModel, i::Int) = constraint_voltage_magnitude_difference(pm, pm.cnw, i)


function constraint_branch_current(pm::GenericPowerModel, n::Int, i)
branch = ref(pm, n, :branch, i)
f_bus = branch["f_bus"]
t_bus = branch["t_bus"]
f_idx = (i, f_bus, t_bus)
tm = branch["tap"]
g_sh_fr = branch["g_fr"]
b_sh_fr = branch["b_fr"]

constraint_branch_current(pm::GenericPowerModel, n::Int, i, f_bus, f_idx, g_sh_fr, b_sh_fr, tm)
end
constraint_branch_current(pm::GenericPowerModel, i::Int) = constraint_branch_current(pm, pm.cnw, i)
Loading

0 comments on commit 57c3299

Please sign in to comment.