Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

SDP relaxation with constraint decomposition #333

Merged
merged 33 commits into from
Sep 4, 2018

Conversation

kersulis
Copy link
Contributor

@kersulis kersulis commented Jul 18, 2018

This PR adds functionality for solving the SDP relaxation with decomposition of the semidefinite constraint. Replacing one large semidefinite constraint with a collection of small constraints (each associated with a maximal clique of a chordal extension of the system graph) can dramatically reduce solver times. See here for details. The application of this technique to the SDP OPF problem is due to Dan Molzahn, whose MATLAB implementation I have been comparing against.

Next steps:

  • Ongoing: quantify improvement by comparing MATLAB, PowerModels.jl without decomposition, and PowerModels.jl with decomposition. (Use Mosek for consistency with MATLAB code.) @ccoffrin how do you efficiently generate that big OPF results table? pandas has markdown table export, but I don't see that in DataFrames.jl.
  • Compare SCS and Mosek (results so far suggest SCS may not solve large instances in reasonable times even with decomposition).
  • Add clique merge functionality (application of this technique in SDP OPF setting also due to Dan Molzahn) to further reduce solver time. (currently debugging)
  • Fix infeasible cases
  • Add objective function scaling to improve numerics, Cost Scaling Parameter #368

@codecov
Copy link

codecov bot commented Jul 18, 2018

Codecov Report

Merging #333 into master will increase coverage by 0.04%.
The diff coverage is 97.44%.

Impacted file tree graph

@@            Coverage Diff             @@
##           master     #333      +/-   ##
==========================================
+ Coverage   96.96%   97.01%   +0.04%     
==========================================
  Files          51       51              
  Lines        6697     6891     +194     
==========================================
+ Hits         6494     6685     +191     
- Misses        203      206       +3
Impacted Files Coverage Δ
src/core/base.jl 99.47% <ø> (ø) ⬆️
test/opf.jl 100% <100%> (ø) ⬆️
src/form/wrm.jl 97.69% <97.05%> (-2.31%) ⬇️
src/util/obbt.jl 85.84% <0%> (+0.77%) ⬆️

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 9ecf261...098a568. Read the comment docs.

@ccoffrin
Copy link
Member

I don't see where you add the equality constraints for variables that are shared across the groups.

Also, feel free to put the implementation in wrm.jl

@kersulis
Copy link
Contributor Author

kersulis commented Jul 18, 2018

The WR and WI variables are generated the same as in SDPWRMForm, but in constraint_voltage I generate overlapping semidefinite constraints rather than one large constraint. Because I am just constraining groups of the same WR and WI variables, shouldn't the linking constraints be implicit? i.e. the solver should know there can be only one value for each variable, regardless of how many semidefinite constraints it shows up in.

Here is Mosek's description of the case30 OPF SDP relaxation problem without decomposition (and before presolve):

Constraints            : 1211            
Cones                  : 45              
Scalar variables       : 531             
Matrix variables       : 1               
Integer variables      : 0           

And with decomposition (again, before presolve):

Constraints            : 953             
Cones                  : 45              
Scalar variables       : 531             
Matrix variables       : 12              
Integer variables      : 0               

Same scalar variables, more matrix variables, fewer constraints. This seems reasonable in light of the smaller matrix constraints replacing the single big one.

There may be some benefit to making the linking constraints explicit, but the results I've observed so far are consistent with those of Dan's MATLAB implementation: SDPDecompForm reaches the same objective value as SDPWRMForm in roughly half the time (0.07s vs 0.19s for case30, for example).

When I use SCS with SDPDecompForm, I do get messages like this for some cases (case6 for example):

WARN: A->p (column pointers) not strictly increasing, column 1 empty
WARN: A->p (column pointers) not strictly increasing, column 2 empty
WARN: A->p (column pointers) not strictly increasing, column 3 empty
WARN: A->p (column pointers) not strictly increasing, column 5 empty
WARN: A->p (column pointers) not strictly increasing, column 8 empty
WARN: A->p (column pointers) not strictly increasing, column 9 empty
WARN: A->p (column pointers) not strictly increasing, column 13 empty
WARN: A->p (column pointers) not strictly increasing, column 17 empty
WARN: A->p (column pointers) not strictly increasing, column 19 empty

Despite these messages, SCS still reaches the same objective value, and faster than without decomposition (though the speedup is less for small networks like case6).

@kersulis
Copy link
Contributor Author

This new implementation creates only those voltage product variables that are needed in the decomposed semidefinite constraints. It also makes the linking constraints explicit.

Here is a comparison of Mosek solver output for nesta_case30_ieee:

SDPDecompForm SDPWRMForm
Constraints 2261 3949
Cones 88 88
Scalar variables 834 1811
Matrix variables 26 1
Integer variables 0 0
Objective 204.96827 204.96835
Solver time 0.11 2.90

@kersulis
Copy link
Contributor Author

kersulis commented Aug 2, 2018

@molzahn, @bghaddar

@kersulis
Copy link
Contributor Author

After the cleanup, only SDPDecompForm remains. This implementation does not do clique merge. The only heuristic used is the sparsity-preserving Cholesky factorization that yields a chordal graph extension. I do not currently know of a better way to obtain a chordal extension, but I'll let you know if I find one :)

Copy link
Member

@ccoffrin ccoffrin left a comment

Choose a reason for hiding this comment

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

Suggesting some minor updates to function definitions, also I think there is a small bug in the RSOC case.


""
function variable_voltage(pm::GenericPowerModel{T}; nw::Int=pm.cnw, cnd::Int=pm.ccnd, bounded = true) where T <: AbstractWRMForm
function variable_voltage(pm::GenericPowerModel{SDPWRMForm}; nw::Int=pm.cnw, cnd::Int=pm.ccnd, bounded = true)
Copy link
Member

Choose a reason for hiding this comment

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

This function signature should be,

function variable_voltage(pm::GenericPowerModel{T}; nw::Int=pm.cnw, cnd::Int=pm.ccnd, bounded = true) where T <: SDPWRMForm

Otherwise it cannot be extended later on.

Copy link
Contributor Author

@kersulis kersulis Aug 14, 2018

Choose a reason for hiding this comment

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

Now that there are two methods for variable_voltage, which would you prefer to have the broader signature? I intentionally made this method specific to SDPWRMForm, but meant to use the broader signature at line 98.

end

function variable_voltage(pm::GenericPowerModel{SDPDecompForm}; nw::Int=pm.cnw, cnd::Int=pm.ccnd, bounded=true)
Copy link
Member

Choose a reason for hiding this comment

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

This function signature should be,

function variable_voltage(pm::GenericPowerModel{T}; nw::Int=pm.cnw, cnd::Int=pm.ccnd, bounded = true) where T <: SDPDecompForm

Otherwise it cannot be extended later on.


""
function constraint_voltage(pm::GenericPowerModel{T}, nw::Int, cnd::Int) where T <: AbstractWRMForm
function constraint_voltage(pm::GenericPowerModel{SDPWRMForm}, nw::Int, cnd::Int)
Copy link
Member

Choose a reason for hiding this comment

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

Same as above.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I made this method specific to SDPWRMForm because it creates the full voltage product matrix and adds a single SDconstraint; extensions would presumably incorporate constraint decomposition in some way, which would require the linking constraint method defined at 197.

src/form/wrm.jl Outdated
wi_ij = WI[1, 2]

# standard SOC form (Mosek doesn't like rotated form)
@constraint(pm.model, (wr_ii + wr_jj)/2 >= norm([(wr_ii - wr_jj)/2; wr_ij; wi_ij]))
Copy link
Member

Choose a reason for hiding this comment

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

I think this constraint is off by a factor of 2. Double check the math. In either case, this version should work,

@constraint(pm.model, (wr_ii + wr_jj) >= norm([(wr_ii - wr_jj); 2*wr_ij; 2*wi_ij]))

Copy link
Contributor Author

Choose a reason for hiding this comment

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

@molzahn: We checked this against the derivation yesterday, right?

@kersulis
Copy link
Contributor Author

Would it be more appropriate to rename SDPDecompForm to SDPWRMSparseForm?

@ccoffrin
Copy link
Member

It's up to you. Currently we don't have clear conventions on how to name these types, so I won't be too picky about the name at this point.

You might consider something like SDPWRMCFForm (CF is short for cholfact) or SDPWRMCFCEForm (CFCE is short for cholfact chordal extension). Try to think of a name that will be less likely to collide with other possible implementations.

Also, following in the new model documentation convention, consider adding a doc string with your preferred references.

Copy link
Member

@ccoffrin ccoffrin left a comment

Choose a reason for hiding this comment

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

Minor comments

src/form/wrm.jl Outdated

# standard SOC form (Mosek doesn't like rotated form)
@constraint(pm.model, (wr_ii + wr_jj)/2 >= norm([(wr_ii - wr_jj)/2; wr_ij; wi_ij]))
@constraint(pm.model, (wr_ii + wr_jj) >= norm([(wr_ii - wr_jj); 2*wr_ij; 2*wi_ij]))
@constraint(pm.model, wi_ij == -wi_ji)
Copy link
Member

Choose a reason for hiding this comment

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

FYI, in my past experience adding this symmetry constraint was causing some numerical issues in the solvers.

@ccoffrin ccoffrin merged commit 5262201 into lanl-ansi:master Sep 4, 2018
ccoffrin pushed a commit that referenced this pull request Oct 12, 2018
* prepare changelog for next release

* remove redefinition of InfrastructureModels.ismultinetwork

* Discrepancy BIM vs BFM (#286)

* found for discrepancy SOC BFM vs BIM
* Documentation fix
* Looser bounds on square current magnitude variable

* added note to changelog on bfm changes

* Tidy Matpower Test Cases (#288)

* update case names to be consistent with file names
* remove areas block from some test cases

* add note to change log

* ADD: Warning when fields are missing during PTI parse (#294)

Adds a warning enumerating the missing data fields when parsing a PTI
file.

* Update change log for minor release

* update to memento v0.7 (#299)

removes depreciation warnings.

* prepare change log for v0.7.2 release

* Core Multiphase Features (#297)

* adding multiphase flag to data parsers
* adding ismultiphase test to build_generic_model
* making network and phase kwargs in constraint templates
* making network and phase kwargs in variable constructors
* change usage of pm.ref[:nw] to nws()
* exporting multinetwork and multiphase helper functions
* minor fix to psse test
* multi-phase solution nominally working
* adding JSON output for multiphase data
* making phases optional to support array-based single phase models
* adding isapprox so compare_dict works with multi-phase
* adding ref file and making data functions work with multiphase
* adding MultiPhaseMatrix to support branch parameters
* add support for multiphase duals
* minor bugfix to apply_func
* adding area, zone, and base_kv to phaseless parameters
* adding basic multiphase tests
* adding multinetwork multiphase tests
* adding run_mn_pf formulation for testing
* removed multinetwork models defined in tests
* fixing dcline cost function unit conversion
* replace getstart with getval, closes #252
* ADD: operators for MultiPhaseValues
* REF: multiphasevalues operations
* ADD: Unit tests for MultiPhaseValue operations
* FIX: DistFlow formulation for multiphase (#289)
* Adds more `getmpv` variants
* ADD: Unit tests for Multphase and Multinetwork (#300)
* ADD: Unit tests for Multphase and Multinetwork
* REF: Updated Memento log tests
* ADD: Missing unit tests for check_connectivity (#301)
* ADD: Missing unit tests for check_connectivity
* ADD: Missing ismultiphase test
* UPD: updates network checks for multiphase networks (#302)
* UPD: updates network checks for multiphase networks
* update to changelog

* fix spelling of LICENSE

* adding basic mathmatical model to docs

* minor fix to math-model

* making objective functions addative over phases (#308)

* Matpower export (#309)

* Matpower Export Branch

Intermediate backup of this branch

* Iterative code backup

* intermediate backup

* Matpower export

Complete

* removing .project

* add .project to .gitignore

* Suggested modifications

* updated test due to change to case3.m

* add note to changelog [skip ci]

* adding Matpower export to acknowledgements

[ci skip]

* minor fix to acknowledgments [skip ci]

* File Parser Enhancements (#311)

Closes #310 

* move `check_network_data` from `parse_file`

into the format-specific parsers: `parse_matpower`, `parse_psse` and
`parse_json`

ref
#310 (comment)
87343

* bugfix

* parser support for iostream

implemented for `parse_json`, `parse_matpower`, `parse_psse`, and
`parse_pti`.

ref
#310 (comment)
82468

* remove filename kwarg

there’s a `.name` attribute in IOStream that I wasn’t previously aware
of

* add explicit return statements

* use in memory buffers instead of creating temp files

`::IOStream` -> `::IO`

* move matpower_export tests from test/data.jl to test/matpower.jl

* adding note to changelog

* Remove Series Variables from Branch Flow Formulation (#312)

* removing series variables
* removing dead code and adding pf hvdc test
* rename df.jl to bf.jl, consistency
* rename SOCDF to SOCBF, consistency
* update AbstractDFForm to AbstractBFForm, consistency

* Extended QC Model to Asymmetric Ling Charging (#315)

* adding test case for asymmetric line charging
* adding opf tests for asymmetric line charging
* adding ots tests for asymmetric line charging
* extending constraint_power_magnitude_link to asymetric line charge, closes #284

* Fix Parser for HVDC lines with PWL Costs (#317)

* do not add HVDC line costs when no HVDC lines exist
* fix non-increasing domain on default hvdc pwl costs, and added test case
* closes #316

* Update VAD BF Constraint (#318)

* remove use of getmpv in constraint_voltage_angle_difference
* apply complement to vad constraint

* bugfix calc_branch_y TPPM#45 (#319)

Uses the pseudo-inverse to avoid issues when certain phases are absent in the impedance matrix.

Added corresponding unit test

* fix MultiPhaseValue isapprox, closes #320

* FIX: Vector and Matrix Operations on MultiPhaseValues (#321)

Fixes problems with MultiPhaseVector and MultiPhaseMatrix operations
not being symmetrical, e.g. `Number + MultiPhaseValue` might not work
whereas `MultiPhaseValue + Number` might. This has been fixed.

Added `inv` function for matrices.

Added `real` and `imag` functions for MutliPhaseValues

Updated `calc_branch_y`, and one OTS test whose objective value
changed due to the difference in rounding errors between doing
the calculation in complex vs float types.

Closes #313

* fix subscripts in mathmatical model doc

* FIX: Line splitting in PTI parser (#323)

In String fields, like "name", there can be commas, which was causing
special cases (e.g. transformers) to split in the wrong places. Switched
all cases to use `get_line_elements` function, which relies on RegEx
and handles this case where a comma is in a String field surrounded by
single quotes.

* add note to changelog

* FIX: starbus index (#325)

Fixes bug in PSSE parser during the creation of a starbus for 3-winding
transformers where the starbus number would keep increasing in
magnitude due to `find_max_bus_id` running at every iteration and
finding a new, higher max bus id. Fixes this by excluding any busname
with `starbus` in it, which is added to the name during the creation of
any starbus.

* REF: Disables some warnings about missing fields in PTI (#326)

Specfically excludes warnings about missing fields in the sections
"IMPEDANCE CORRECTION, "MULTI-SECTION LINE", and "SWITCHED
SHUNT", which all have optional fields at the end, starting with "T#",
"DUM#", and "N#", respectively.

* adding video links to readme

* FIX: Single-quote in busnames (#329)

* FIX: Single-quote in busnames

Fixes problem where there could be a quotation mark inside of a bus
name by taking advantage of the fact that all busnames in PTI v33.3
must be exactly 12 characters long to be valid.

* ADD: Update PTI file to test change

Changes one of the busnames to `'1' DROP     '` to test change to RegEx

* Multi-conductor Renaming (#324)

* phase to conductor renaming
* update MultiPhase datatypes to MultiConductor
* add standard ids names to developer docs
* updating file names
* rename function args and tmp vars to reflect conductor short names

* update for info default logging level, closes #296 (#330)

* Some Minor Fixes (#331)

* update calc_branch_t to use function dispatch instead of map (close #271 )
* rename mp functions to mc functions
* rename mpv to mcv and reduce use of getmcv (close #295)
* adding doc string to post_tp_opf in mc tests

* adding link to network model video

* adding con help to exports

* update changelog for v0.8 release

* adding implicit single phase to support to buspairs

* make add_setpoint more flexible with different types

* adding atol to check_tnep_status

* update changelog for minor release

* fix typo in changelog

* bus fix TNEP voltage variable definitions

* updating tests for numerical stablity

* Tightening QCWRTriModel for opf (#332)

* update case5 to test non-contiguous bus ids
* tightening QC Tri
* update changelog

* update MINLP solvers used in testing

* OBBT functionality (#335)

* obbt initial commit
* add docs and stats

* add utilities to the docs menu

* adding kaarthik sundar to acknowledgments

* FIX: import_remaining array eltype check (#336)

import_remaining! did not check for eltype of arrays, which broke
behavior in ThreePhasePowerModels. Now checks to ensure that arrays
have Dict elements.

Closes ThreePhasePowerModels#51

* FIX: check_thermal_limits was returning NaN (#337)

During check thermal limits, `r` or `x` are zero, `rate_{a,b,c}` would get
set to NaN. Using complex formulation and `pinv` to solve this problem.

* Adding AbstractBFConicForm (#346)

* adding AbstractBFConicForm
* adding conic variant of BF model for testing

* ADD: Error handling for UnicodeError (#347)

Adds Unicode Error message that informs user only UTF-8 and ASCII files
are supported, and to re-encode file into a supported format.

Adds unit test

Closes #328

* Adding #328 note to changelog

* update changelog for v0.8.2

* Updated wr.jl

Removed a rogue + sign in line 811

* fixing bug in conic objtive (#348)

* Add Parallel Time Stat to OBBT (#349)

* improve update_data warning message
* adding sim_parallel_run_time to the obbt stats
* add note to changelog

* changlog typo fix

* Added Simplification for Piecewise Linear Cost Functions (#350)

* simplified multiconductor support in check_cost_functions
* adding pwl cost simplification code
* adding cost functions checks to matpower exporter
* add note to changelog

* expanding math model + link to constraint template (#351)

* Updated docs (#352)

* improving math model docs
* Adding references for formulations
* References for formulations
* Following recommendations in  #343
* Updated changelog

* Updating docs to address #8 (#356)

* improving math model docs
* Adding references for formulations
* References for formulations
* Following recommendations in  #343
* Updated changelog
* quick fixes per feedback
* Addressing #8 + tweaking docs
* adding Carpentier reference

* Add Support for Current Limits (#354)

* add current rating to the branch model
* making rates optional
* fixing tnep test
* add note on optional ratings to docs
* adding base current constraint implementation
* adding tests for current based opf
* add note to changelog
* adding test for specifying current limits in matpower data
* adding ka_base and rescale_ampere and to unit conversions

* Improve OBBT Solver Robustness  (#358)

* converted error to warnings

* cosmetic updated to bibtex entries (#367)

* OBBT Change more Errors to Warnings (#370)

* Conic form of SOC relaxation (#371)

* Add SOC conic form
* Tests for SOC conic form

* adding SOCWRConicPowerModel to changelog

* remove soc conic 6 bus test case due to numerics

* use InfrastructureModels.relaxation_complex_product_conic, closes #373 (#374)

* Basic Network Flow and Copper Plate Models (#365)

* adding network flow approximation formulation
* fixing bug in case6 so that there are two components
* adding optimal power balance problem

* Add version number to Project.TOML

* add v0.8.3 version number to change log

* Update uuid number according to the registry

* Add uuid from Registry Repo

* adding julia version upper bound

* updating publication reference

* Obbt exp (#379)

* revised obbt bound update criteria
* Message updates
* Add note to change log

* SDP relaxation with constraint decomposition (#333)

* working SDP decomposition
* add tests (apply existing SDP tests to SDPDecompForm)
* explicit generation of voltage vars and linking constraints
* move constraint decomp into struct to keep pm.ext clean
* no more infeasible cases! Fixes kersulis#1
* keep only the best decomp implementation, fixes kersulis#2
* bools and ints, not ints and floats
* broaden signatures for variable_voltage and constraint_voltage
* store graph ordering (suff. to reconstruct minimal chordal ext, cliques)
* do not store full adjacency matrix

* Minor Tidying of WRM Forms (#380)

* rename sdp decomp to sparse
* reorg SDP WRM functions
* adding sparse sdp to docs
* add note to change log
* add sdp thanks to readme

* adding indentation to ACR doc string

* Update Manifest.TOML to newer versions

* Citations for sparsity-exploiting SDP OPF (#385)

* fix indentation of sparse sdp docstring

* update change log for v0.8.4 release

* update jump constructor

* updating ipopt constructor

* update to MOI 0.6 and JuMP master

* omit newly broken tests

* add link to power flow forms docs to readme

* restor power flow tests

* adding explicit starting points to variables, restored pwl and qc tri models

* Typo fix [ci skip]

Fixed a spelling error

* WIP, SDP support

* minor fix to subscripts in math model docs

* update to latest scs version, v0.4.0

* update SCS settings

* remove comment

* removing extra square from polynomial cost function

* update base to use optimization factories

* small sdp models working

* Update Testing to load default packages

* Update testing

* update travis.yml

* update travis.yml

* Fix LinearAlgebra deprecations

* fix deprecation of assert to @Assert

* Enable testing of Matpower to find deprecations

* Update to MOI

* Update CI

* Update Manifest

* Fix typo

* add SparseArrays

* remove Base.var overloading

* Update overloading from Base. to LinearAlgebra.

* Remove more deprecated syntax

* Remove whitespace deprecations
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants