update examples

Author: jalving

Manifest.toml

@@ -9,8 +9,12 @@ KaHyPar = "2a6221f6-aa48-11e9-3542-2d9e0ef01880"
 LightGraphs = "093fc24a-ae57-5d10-9952-331d41423f4d"
 MPI = "da04e1cc-30fd-572f-bb4f-1f8673147195"
 MPIClusterManagers = "e7922434-ae4b-11e9-05c5-9780451d2c66"
+Measures = "442fdcdd-2543-5da2-b0f3-8c86c306513e"
+PipsNLP = "0e2b999c-8a3a-11ea-3c4b-c36350f324a0"
 Plasmo = "d3f7391f-f14a-50cc-bbe4-76a32d1bad3c"
+PlasmoPlots = "a21664a3-f668-4e98-8a13-7892286ec16a"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 PowerModels = "c36e90e8-916a-50a6-bd94-075b64ef4655"
+PyPlot = "d330b81b-6aea-500a-939a-2ce795aea3ee"
 Requires = "ae029012-a4dd-5104-9daa-d747884805df"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"

Project.toml

@@ -1,12 +1,12 @@
 # PlasmoExamples
-Example scripts that showcase the core functionality of [Plasmo.jl](https://github.com/zavalab/Plasmo.jl).  This repository shows how to:
-* Use optigraphs to create and solve hierarchical optimization models
-* Solve a natural gas optimal control problem with [PipsSolver.jl](https://github.com/zavalab/PipsSolver.jl),
-* Solve a DC optimal power flow problem with [SchwarzSolver.jl](https://github.com/zavalab/SchwarzSolver.jl)
+Example scripts that showcase using [Plasmo.jl](https://github.com/zavalab/Plasmo.jl).  This repository shows how to:
+* Use `OptiGraphs` to model and solve optimization models
+* Solve a natural gas optimal control problem with [PipsNLP.jl](https://github.com/zavalab/PipsNLP.jl) using partitioning to obtain subproblems
+* Solve a DC optimal power flow problem with [SchwarzOpt.jl](https://github.com/zavalab/SchwarzOpt.jl) using overlapping schwarz decomposition
 
 
 ## Requirements
-* Julia 1.0.5
+* Julia 1.3 or later
 * The GR plotting backend to create plots for the example scripts.
 * PIPS-NLP with the ma57 linear solver. Use the [PIPS-NLP](https://github.com/Argonne-National-Laboratory/PIPS) installation instructions.
 * A working MPI installation for PIPS-NLP.  The scripts in this repository have been tested with [mpich](https://www.mpich.org/downloads/) version 3.3.2
@@ -20,7 +20,7 @@ The Julia packages needed to run the examples and case studies can be installed
 To setup the necessary Julia environment, first clone the PlasmoExamples directory with:
 
 ```
-git clone https://github.com/jalving/PlasmoExamples.git
+git clone https://github.com/zavalab/PlasmoExamples.git
 ```
 
 Once cloned, navigate to the `PlasmoExamples` folder and begin a Julia REPL session. The necessary Julia dependencies can then be downloaded as follows:
@@ -78,11 +78,9 @@ $ julia run_examples.jl
 This will run the 5 example scripts in the examples folder and create corresponding plots in `PlasmoExamples/examples/figures`. The example scripts perform the following:
 
 * example1.jl creates a simple optigraph, solves it with GLPK, and plots the optigraph structure.
-* example2.jl creates a hierarchical optigraph with a shared optiedge (linkconstraint), solves it with GLPK, and plots the optigraph structure.
-* example3.jl creates a hierarchical optigraph with a shared optinode, solves it with GLPK, and plots the optigraph structure.
-* example4.jl creates an optigraph for an optimal control model, partitions it with KaHyPar, and aggregates the resulting subgraphs into new optinodes.  This example
-also creates plots for the optigraph structure corresponding to the original optigraph, the partitioned optigraph, and the aggregated optigraph.
-* example5.jl builds on example4.jl and creates overlapping subgraphs.  It also plots the overlapping graph structures.
+* example2.jl creates a hierarchical optigraph with a shared optiedge and optinode (i.e. it contains linkconstraints), solves it with GLPK, and plots the optigraph structure.
+* example3.jl creates an optigraph for an optimal control model and produces the plot.
+* example4.jl builds on example3.jl and partitions the control problem with KaHyPar, aggregates the resulting subgraphs into new optinodes, and also shows how to create overlapping subgraphs. It also plots the resulting graph structures.
 
 
 ## Running Case Study 1
@@ -112,10 +110,10 @@ julia> include("run_casestudy1.jl")
 This will run the case study just like before.  However, it is now possible to modify partitioning parameters (`n_parts`, `n_processes`, and `max_imbalance`) in the `run_casestudy1.jl` script  
 (by modifying the script with any text editor) without restarting Julia. It is also possible to scale the model size by modifying `nx` (the number of space discretization points for each pipeline).
 After modifying the script, just run the `include("run_casestudy1.jl")` statement again.
-We recommend always keeping `n_parts` the same as `n_processes`.  We have observed some issues with `PipsSolver.jl` and non-uniform sub-problem allocations.
+We recommend always keeping `n_parts` the same as `n_processes`.  We have observed some issues with `PipsNLP.jl` and non-uniform sub-problem allocations.
 
 ## Running Case Study 2
-Case study 2 solves a DC optimal power flow problem with `SchwarzSolver.jl`. The solver can take advantage of multiple threads set with the
+Case study 2 solves a DC optimal power flow problem with `SchwarzOpt.jl`. The solver can take advantage of multiple threads set with the
 `JULIA_NUM_THREADS` environment variable.  Since the case study partitions the problem into 4 partitions, it can benefit from up to 4 threads which solve
 the sub-problems in parallel.  The following commands can be used to set the number of threads and execute the case study:
 
@@ -146,25 +144,4 @@ formulated (small partitions for a network can lead to feasibility issues).
 
 ## Other documentation
 Documentation describing Plasmo.jl and its underlying functions can be found at the [github pages](https://zavalab.github.io/Plasmo.jl/dev/) site.
-The source code for v0.3.0 can be found at: [https://github.com/zavalab/Plasmo.jl/tree/v0.3.0](https://github.com/zavalab/Plasmo.jl/tree/v0.3.0).
-
-### Overview of source code files
-The Plasmo.jl source code for version 0.3.0 is defined by modules located at [https://github.com/zavalab/Plasmo.jl/tree/v0.3.0/src](https://github.com/zavalab/Plasmo.jl/tree/v0.3.0/src).
-Here is a brief overview of the modules:
-
-* hypergraphs/hypergraph.jl - defines the `HyperGraph` object used to create a hypergraph representation of an optigraph. Notably, it extends a `LightGraphs.AbstractGraph`, so common lightgraphs
-operations work with the `HyperGraph` object (such as `neighbors` and `incidence_matrix`).
-* extras/kahypar.jl - defines a simple conveniece wrapper for the KaHyPar package if it is loaded. Makes it easier to directly partition an optigraph using KaHyPar.
-* extras/plots.jl - defines special plots for Plasmo.jl if the `Plots.jl` package is loaded. Extends the `Plots.plot` and `Plots.spy` functions for the `OptiGraph` object.
-* combine.jl - contains functions to aggregate an optigraph into an optinode, or to aggregate optigraph subgraphs into optinodes.
-* copy.jl - contains functions to copy objective and constraint expressions.  Used extensively for the aggregate functions defined in combine.jl.
-* graph_functions.jl - defines functions to perform graph operations on an optigraph such as querying the neighborhood around a set of nodes or performing subgraph expansion. Uses
-hypergraph.jl to perform said operations.
-* graph_interface.jl - currently only contains a function to generate a `HyperGraph` object from an `OptiGraph` object. Eventually other graph representations will be available.
-* macros.jl - contains the Julia macros used to model with optigraphs.  
-* optiedge.jl - defines the `OptiEdge` and `LinkConstraint` objects, as well as functions that work on these objects.
-* optigraph.jl - defines the `OptiGraph` objects, as well as functions that work on an optigraph. Notably, an `OptiGraph` extends a `JuMP.AbstractModel` so JuMP macros work with optigraphs.
-* optinode.jl - defines the `OptiNode` object, as well as functions that work on optinodes. The optinode also extends a `JuMP.AbstractModel`.
-* partition.jl - defines the `Partition` object which informs how to partition an optigraph.  The `Partition` object is used in the `make_subgraphs!` function to create subgraphs in an optigraph.
-* solve.jl - extends JuMP.jl functions to solve an optigraph using JuMP compatible solvers.
-* Plasmo.jl - imports and exports functions and defines abstract types.  
+The source code for v0.3.0 can be found at: [https://github.com/zavalab/Plasmo.jl/tree/v0.3.0](https://github.com/zavalab/Plasmo.jl/tree/v0.4.0).

README.md

@@ -6,63 +6,66 @@ using JLD2
 horizon = 24*3600   #the time horizon is in seconds
 nt = 24              #number of time points
 dt = horizon / (nt - 1) #time delta
-#nx = 3
+#nx = 10
 
 include((@__DIR__)*"/modelfunctions.jl")
 include((@__DIR__)*"/load_data.jl")
 
-gas_network = OptiGraph()
-for pipe in pipelines
-    add_subgraph!(gas_network,pipe)
-end
-for compressor in compressors
-    add_subgraph!(gas_network,compressor)
-end
-for junction in junctions
-    add_subgraph!(gas_network,junction)
-end
+function create_gas_network_problem()
+    gas_network = OptiGraph()
+    for pipe in pipelines
+        add_subgraph!(gas_network,pipe)
+    end
+    for compressor in compressors
+        add_subgraph!(gas_network,compressor)
+    end
+    for junction in junctions
+        add_subgraph!(gas_network,junction)
+    end
 
-for pipe in pipelines
-    junction_from,junction_to = pipe_map[pipe]
-    @linkconstraint(gas_network,[t = 1:nt],pipe[:pin][t] == junction_from[:time_nodes][t][:pressure])
-    @linkconstraint(gas_network,[t = 1:nt],pipe[:pout][t] == junction_to[:time_nodes][t][:pressure])
-end
+    for pipe in pipelines
+        junction_from,junction_to = pipe_map[pipe]
+        @linkconstraint(gas_network,[t = 1:nt],pipe[:pin][t] == junction_from[:time_nodes][t][:pressure])
+        @linkconstraint(gas_network,[t = 1:nt],pipe[:pout][t] == junction_to[:time_nodes][t][:pressure])
+    end
 
-for compressor in compressors
-    junction_from,junction_to = compressor_map[compressor]
-    nodes = compressor[:time_nodes]
-    @linkconstraint(gas_network,[t = 1:nt],nodes[t][:psuction] == junction_from[:time_nodes][t][:pressure])
-    @linkconstraint(gas_network,[t = 1:nt],nodes[t][:pdischarge]  == junction_to[:time_nodes][t][:pressure])
-end
+    for compressor in compressors
+        junction_from,junction_to = compressor_map[compressor]
+        nodes = compressor[:time_nodes]
+        @linkconstraint(gas_network,[t = 1:nt],nodes[t][:psuction] == junction_from[:time_nodes][t][:pressure])
+        @linkconstraint(gas_network,[t = 1:nt],nodes[t][:pdischarge]  == junction_to[:time_nodes][t][:pressure])
+    end
 
-for junction in junctions
-    devices_in = junction_map_in[junction]
-    devices_out = junction_map_out[junction]
+    for junction in junctions
+        devices_in = junction_map_in[junction]
+        devices_out = junction_map_out[junction]
 
-    if length(devices_in) > 0
-        flow_in = [sum(devices_in[i][:fout][t] for i = 1:length(devices_in)) for t = 1:nt]
-    else
-        flow_in = zeros(nt)
-    end
-    if length(devices_out) > 0
-        flow_out = [sum(devices_out[i][:fin][t] for i = 1:length(devices_out)) for t = 1:nt]
-    else
-        flow_out = zeros(nt)
-    end
-    total_supplied = [junction[:time_nodes][t][:total_supplied] for t = 1:nt]
-    total_delivered = [junction[:time_nodes][t][:total_delivered] for t = 1:nt]
+        if length(devices_in) > 0
+            flow_in = [sum(devices_in[i][:fout][t] for i = 1:length(devices_in)) for t = 1:nt]
+        else
+            flow_in = zeros(nt)
+        end
+        if length(devices_out) > 0
+            flow_out = [sum(devices_out[i][:fin][t] for i = 1:length(devices_out)) for t = 1:nt]
+        else
+            flow_out = zeros(nt)
+        end
+        total_supplied = [junction[:time_nodes][t][:total_supplied] for t = 1:nt]
+        total_delivered = [junction[:time_nodes][t][:total_delivered] for t = 1:nt]
 
-    @linkconstraint(gas_network,[t = 1:nt], flow_in[t] - flow_out[t] + total_supplied[t] - total_delivered[t] == 0)
-end
+        @linkconstraint(gas_network,[t = 1:nt], flow_in[t] - flow_out[t] + total_supplied[t] - total_delivered[t] == 0)
+    end
 
 
-#Fix Demands
-for (i,j_data) in junction_data
-    jmodel = jmap[i]
-    dvalues = j_data[:demand_values]
-    n_demands = length(dvalues)
-    nodes = jmodel[:time_nodes]
-    for (t,node) in enumerate(nodes)
-        @constraint(node,[d = 1:n_demands],node[:fdemand][d] == dvalues[d][t])
+    #Fix Demands
+    for (i,j_data) in junction_data
+        jmodel = jmap[i]
+        dvalues = j_data[:demand_values]
+        n_demands = length(dvalues)
+        nodes = jmodel[:time_nodes]
+        for (t,node) in enumerate(nodes)
+            @constraint(node,[d = 1:n_demands],node[:fdemand][d] == dvalues[d][t])
+        end
     end
+    return gas_network
 end

case_studies/gasnetwork/gasnetwork.jl

@@ -60,7 +60,7 @@ function create_pipeline_model(data,nt,nx)
         @variable(node, 1 <= px <= 100,start = 10)
         @variable(node, 0 <= fx <= 100, start = 10)
         @variable(node,slack >= 0, start = 1)
-        @NLnodeconstraint(node, slack*px - c3*fx*fx == 0)
+        @NLconstraint(node, slack*px - c3*fx*fx == 0)
     end
 
     #Setup boundary variables
@@ -95,7 +95,7 @@ function create_compressor_model(data,nt)
         @variable(node, flow >= 0, start = 10)
         @variable(node, 1 <= boost <= 30,start = 1)
         @constraint(node, pdischarge == psuction + boost)
-        @NLnodeconstraint(node, power == c4*flow*((pdischarge/psuction)^om-1) )
+        @NLconstraint(node, power == c4*flow*((pdischarge/psuction)^om-1) )
         @objective(node, Min, 0.1*power*(dt/3600.0))
     end
     @expression(graph,fin[t=1:nt],time_nodes[t][:flow])

case_studies/gasnetwork/modelfunctions.jl

@@ -1,11 +1,7 @@
-#create gasnetwork optigraph model
-include((@__DIR__)*"/gasnetwork.jl")
-
 #Obtain a hyergraph representation of the gas_network
-hypergraph,hyper_map = gethypergraph(gas_network)
+hypergraph,hyper_map = hyper_graph(gas_network)
 
 #Setup node and edge weights
-n_vertices = length(vertices(hypergraph))
 node_sizes = [num_variables(node) for node in all_nodes(gas_network)]
 edge_weights = [num_linkconstraints(edge) for edge in all_edges(gas_network)]
 
@@ -14,19 +10,18 @@ node_vector = KaHyPar.partition(hypergraph,n_parts,configuration = :edge_cut,
 imbalance = max_imbalance, node_sizes = node_sizes, edge_weights = edge_weights)
 
 #Create an optigraph partition
-partition = Partition(gas_network,node_vector,hyper_map)
+partition = Partition(node_vector,hyper_map)
 
 #Setup subgraphs based on partition
-make_subgraphs!(gas_network,partition)
+apply_partition!(gas_network,partition)
 
-#Combine the subgraphs into model-nodes
+#aggregate the subgraphs into optinodes
 combined_graph , combine_map  = aggregate(gas_network,0)
 
 #scale the objective function on each optinode
 for node in all_nodes(combined_graph)
     @objective(node,Min,1e-6*objective_function(node))
 end
-
 ##############################
 # Solve with PIPS-NLP
 ##############################
@@ -35,10 +30,11 @@ julia_workers = sort(collect(values(manager.mpi2j)))
 
 #Distribute the optigraph among the workers
 #Here, we create the variable `pipsgraph` on each worker
-remote_references = PipsSolver.distribute(combined_graph,julia_workers,remote_name = :pipsgraph)
+remote_references = PipsNLP.distribute_optigraph(combined_graph,julia_workers,remote_name = :pipsgraph)
 
+#NOTE: PIPS hits restoration phase; needs better initialization
 #Solve with PIPS-NLP on each mpi rank
 @mpi_do manager begin
     using MPI
-    PipsSolver.pipsnlp_solve(pipsgraph)
+    PipsNLP.pipsnlp_solve(pipsgraph)
 end

case_studies/gasnetwork/partition_and_solve.jl

@@ -0,0 +1,26 @@
+using Plasmo
+
+function create_dynamic_problem()
+    T = 100          #number of time points
+    d = sin.(1:T)    #disturbance vector
+
+    graph = OptiGraph()
+    @optinode(graph,state[1:T])
+    @optinode(graph,control[1:T-1])
+
+    for node in state
+        @variable(node,x)
+        @constraint(node, x >= 0)
+        @objective(node,Min,x^2)
+    end
+    for node in control
+        @variable(node,u)
+        @constraint(node, u >= -1000)
+        @objective(node,Min,u^2)
+    end
+
+    @linkconstraint(graph,[i = 1:T-1],state[i+1][:x] == state[i][:x] + control[i][:u] + d[i])
+    n1 = state[1]
+    @constraint(n1,n1[:x] == 0)
+    return graph
+end

examples/create_dynamic_problem.jl

@@ -0,0 +1,25 @@
+function simple_optigraph()
+    graph = OptiGraph()
+
+    @optinode(graph,n1)
+    @variable(n1, y >= 2)
+    @variable(n1,x >= 0)
+    @constraint(n1,x + y >= 3)
+    @objective(n1, Min, y)
+
+    @optinode(graph,n2)
+    @variable(n2, y)
+    @variable(n2,x >= 0)
+    @constraint(n2,x + y >= 3)
+    @objective(n2, Min, y)
+
+    @optinode(graph,n3)
+    @variable(n3, y )
+    @variable(n3,x >= 0)
+    @constraint(n3,x + y >= 3)
+    @objective(n3, Min, y)
+
+    #Create a link constraint linking the 3 models
+    @linkconstraint(graph, n1[:x] + n2[:x] + n3[:x] == 3)
+    return graph
+end

examples/create_optigraph.jl

@@ -1,40 +1,42 @@
 using Plasmo
-using Plots
-using PlasmoPlots
+using Plots, PlasmoPlots
 using GLPK
 
-graph1 = OptiGraph()
+graph = OptiGraph()
 
-@optinode(graph1,n1)
+@optinode(graph,n1)
 @variable(n1, y >= 2)
 @variable(n1,x >= 0)
 @constraint(n1,x + y >= 3)
 @objective(n1, Min, y)
 
-@optinode(graph1,n2)
+@optinode(graph,n2)
 @variable(n2, y)
 @variable(n2,x >= 0)
 @constraint(n2,x + y >= 3)
 @objective(n2, Min, y)
 
-@optinode(graph1,n3)
+@optinode(graph,n3)
 @variable(n3, y )
 @variable(n3,x >= 0)
 @constraint(n3,x + y >= 3)
 @objective(n3, Min, y)
 
 #Create a link constraint linking the 3 models
-@linkconstraint(graph1, n1[:x] + n2[:x] + n3[:x] == 3)
-set_optimizer(graph1,GLPK.Optimizer)
-optimize!(graph1)
+@linkconstraint(graph, n1[:x] + n2[:x] + n3[:x] == 3)
+set_optimizer(graph,GLPK.Optimizer)
+optimize!(graph)
 
 #Query Solution
 println("n1[:x] = ",value(n1[:x]))
 println("n2[:x] = ",value(n2[:x]))
 println("n3[:x] = ",value(n3[:x]))
-println("objective value = ",objective_value(graph1))
+println("objective value = ",objective_value(graph))
 
-
-plt_graph1 = Plots.plot(graph1,node_labels = true,markersize = 60,labelsize = 30,
+plt_graph = PlasmoPlots.layout_plot(graph,node_labels = true,markersize = 60,labelsize = 30,
 linewidth = 4,layout_options = Dict(:tol => 0.01,:iterations => 2));
-plt_matrix1 = Plots.spy(graph1,node_labels = true,markersize = 30);
+plt_matrix = PlasmoPlots.matrix_plot(graph,node_labels = true,markersize = 40,labelsize = 32);
+
+Plots.savefig(plt_graph,"figures/example1_layout.pdf")
+
+Plots.savefig(plt_matrix,"figures/example1_matrix.pdf")