diff --git a/examples/cut_kKaHyPar_sea20.ini b/examples/cut_kKaHyPar_sea20.ini
new file mode 100644
index 0000000..9fa1d6b
--- /dev/null
+++ b/examples/cut_kKaHyPar_sea20.ini
@@ -0,0 +1,70 @@
+# general
+mode=direct
+objective=cut
+seed=-1
+cmaxnet=1000
+vcycles=0
+# main -> preprocessing -> min hash sparsifier
+p-use-sparsifier=true
+p-sparsifier-min-median-he-size=28
+p-sparsifier-max-hyperedge-size=1200
+p-sparsifier-max-cluster-size=10
+p-sparsifier-min-cluster-size=2
+p-sparsifier-num-hash-func=5
+p-sparsifier-combined-num-hash-func=100
+# main -> preprocessing -> community detection
+p-detect-communities=true
+p-detect-communities-in-ip=true
+p-reuse-communities=false
+p-max-louvain-pass-iterations=100
+p-min-eps-improvement=0.0001
+p-louvain-edge-weight=hybrid
+# main -> coarsening
+c-type=ml_style
+c-s=1
+c-t=160
+# main -> coarsening -> rating
+c-rating-score=heavy_edge 
+c-rating-use-communities=true
+c-rating-heavy_node_penalty=no_penalty
+c-rating-acceptance-criterion=best_prefer_unmatched
+c-fixed-vertex-acceptance-criterion=fixed_vertex_allowed
+# main -> initial partitioning
+i-mode=recursive
+i-technique=multi
+# initial partitioning -> coarsening
+i-c-type=ml_style
+i-c-s=1
+i-c-t=150
+# initial partitioning -> coarsening -> rating
+i-c-rating-score=heavy_edge 
+i-c-rating-use-communities=true
+i-c-rating-heavy_node_penalty=no_penalty
+i-c-rating-acceptance-criterion=best_prefer_unmatched
+i-c-fixed-vertex-acceptance-criterion=fixed_vertex_allowed
+# initial partitioning -> initial partitioning
+i-algo=pool
+i-runs=20
+# initial partitioning -> bin packing
+i-bp-algorithm=worst_fit
+i-bp-heuristic-prepacking=false
+i-bp-early-restart=true
+i-bp-late-restart=true
+# initial partitioning -> local search
+i-r-type=twoway_fm
+i-r-runs=-1
+i-r-fm-stop=simple
+i-r-fm-stop-i=50
+# main -> local search
+r-type=kway_fm_hyperflow_cutter
+r-runs=-1
+r-fm-stop=adaptive_opt
+r-fm-stop-alpha=1
+r-fm-stop-i=350
+# local_search -> flow scheduling and heuristics
+r-flow-execution-policy=exponential
+# local_search -> hyperflowcutter configuration
+r-hfc-size-constraint=mf-style
+r-hfc-scaling=16
+r-hfc-distance-based-piercing=true
+r-hfc-mbc=true
\ No newline at end of file
diff --git a/examples/optimal_control.jl b/examples/optimal_control.jl
index fe90cc2..6c842eb 100644
--- a/examples/optimal_control.jl
+++ b/examples/optimal_control.jl
@@ -5,9 +5,9 @@ using SchwarzOpt
 
 T = 100              #number of time points
 d = sin.(1:T)        #a disturbance vector
-imbalance = 0.1      #partition imbalance
-distance = 5         #expand distance
-n_parts = 6          #number of partitions
+imbalance = 0.05      #partition imbalance
+distance = 2      #expand distance
+n_parts = 4          #number of partitions
 
 #Create the model (an optigraph)
 graph = OptiGraph()
@@ -34,35 +34,55 @@ end
 #Partition the optigraph using recrusive bisection over a hypergraph
 hypergraph,hyper_map = hyper_graph(graph) #create hypergraph object based on graph
 
-partition_vector = KaHyPar.partition(hypergraph,n_parts,
-configuration = (@__DIR__)*"/cut_rKaHyPar_sea20.ini",
-imbalance = imbalance)
-
-partition = Partition(hypergraph,partition_vector,hyper_map)
-apply_partition!(graph,partition)
+partition_vector = KaHyPar.partition(
+    hypergraph,
+    n_parts,
+    configuration = (@__DIR__)*"/cut_kKaHyPar_sea20.ini",
+    imbalance = imbalance
+)
 
-#Inspect the graph structure. It should be RECURSIVE_GRAPH
-# println(Plasmo.graph_structure(graph))
+partition = Partition(hypergraph, partition_vector, hyper_map)
+apply_partition!(graph, partition)
 
 #calculate subproblems using expansion distance
 subgraphs = getsubgraphs(graph)
-expanded_subgraphs = Plasmo.expand.(graph,subgraphs,distance)
-sub_optimizer = optimizer_with_attributes(Ipopt.Optimizer,"print_level" => 0)
+expanded_subgraphs = Plasmo.expand.(graph, subgraphs, distance)
+sub_optimizer = optimizer_with_attributes(Ipopt.Optimizer, "print_level" => 0)
 
 
+# optimizer = SchwarzOpt.Optimizer(
+#     graph,
+#     expanded_subgraphs; 
+#     sub_optimizer=sub_optimizer, 
+#     max_iterations=50,
+#     mu=1.0
+# )
+# SchwarzOpt._initialize_optimizer!(optimizer)
 
-optimizer = SchwarzOpt.Optimizer(
-    graph,
-    expanded_subgraphs; 
-    sub_optimizer=sub_optimizer, 
-    max_iterations=50,
-    mu=1.0
-)
-SchwarzOpt._initialize_optimizer!(optimizer)
+# for subgraph in optimizer.subproblem_graphs
+#     JuMP.set_optimizer(subgraph, optimizer.sub_optimizer)
+# end
+
+# for subproblem_graph in optimizer.subproblem_graphs
+#     SchwarzOpt._update_subproblem!(optimizer, subproblem_graph)
+# end
+
+# # sg1 = optimizer.subproblem_graphs[1]
+# # xk,lk = SchwarzOpt._do_iteration(sg1)
+
+# for subproblem_graph in optimizer.subproblem_graphs
+#     xk,lk = SchwarzOpt._do_iteration(subproblem_graph)
+#     println(xk)
+#     println(lk)
+# end
+
+# inc_edges = SchwarzOpt._find_boundaries(graph, optimizer.subproblem_graphs)
 
 #optimize using schwarz overlapping decomposition
-# SchwarzOpt.optimize!(graph;
-# subgraphs = expanded_subgraphs,
-# sub_optimizer = sub_optimizer,
-# max_iterations = 50
-# )
+SchwarzOpt.optimize!(
+    graph;
+    subgraphs = expanded_subgraphs,
+    sub_optimizer = sub_optimizer,
+    max_iterations = 200,
+    mu=0.001
+)
diff --git a/src/optimizer.jl b/src/optimizer.jl
deleted file mode 100644
index 6cd0804..0000000
--- a/src/optimizer.jl
+++ /dev/null
@@ -1,537 +0,0 @@
-@kwdef mutable struct Timers
-    start_time::Float64 = 0.
-    initialize_time::Float64 = 0.
-    eval_objective_time::Float64 = 0.
-    eval_primal_feasibility_time::Float64 = 0.
-    eval_dual_feasibility_time::Float64 = 0.
-    update_subproblem_time::Float64 = 0.
-    solve_subproblem_time::Float64 = 0.
-    total_time::Float64 = 0.
-end
-
-mutable struct Optimizer <: MOI.AbstractOptimizer
-    graph::OptiGraph                              #Subgraphs should not have any overlap and should cover all the nodes
-    subproblem_graphs::Vector{OptiGraph}          #These are expanded subgraphs (i.e. subgraphs with overlap)
-
-    sub_optimizer::Any
-    tolerance::Float64
-    max_iterations::Int64
-
-    x_in_indices::Dict{OptiGraph,Vector{Int64}}   #primals into subproblem
-    l_in_indices::Dict{OptiGraph,Vector{Int64}}   #duals into subproblem
-    x_out_indices::Dict{OptiGraph,Vector{Int64}}  #primals out of subproblem
-    l_out_indices::Dict{OptiGraph,Vector{Int64}}  #duals out of subproblem
-
-    node_subgraph_map::Dict                       #map optinodes to "owning" subgraph
-    expanded_subgraph_map::Dict
-    incident_variable_map::Dict                   #map variable indices in optigraphs to ghost (copy) variables
-    primal_links::Vector                          #link constraints that do primal updates
-    dual_links::Vector                            #link constraints that do dual updates
-
-    x_out_vals::Dict{OptiGraph,Vector{Float64}}   #primal values out of each subproblem
-    l_out_vals::Dict{OptiGraph,Vector{Float64}}   #dual values out of each subproblem
-    x_vals::Vector{Float64}                       #current primal values that get communicated to subproblems
-    l_vals::Vector{Float64}                       #current dual values that get communicated to subproblems
-
-    status::MOI.TerminationStatusCode
-    err_pr::Float64
-    err_du::Float64
-    objective_value::Float64
-
-    iteration::Int64
-    primal_error_iters::Vector{Float64}
-    dual_error_iters::Vector{Float64}
-    objective_iters::Vector{Float64}
-    solve_time::Float64
-
-    #plasmo_optimizer_hook::Function
-    timers::Timers
-    loaded::Bool
-
-    function Optimizer()
-        optimizer = new()
-
-        optimizer.subproblem_graphs = Vector{OptiGraph}()
-        optimizer.x_in_indices = Dict{OptiGraph,Vector{Int64}}()
-        optimizer.l_in_indices = Dict{OptiGraph,Vector{Int64}}()
-        optimizer.x_out_indices = Dict{OptiGraph,Vector{Int64}}()
-        optimizer.l_out_indices = Dict{OptiGraph,Vector{Int64}}()
-
-        optimizer.node_subgraph_map = Dict()
-        optimizer.expanded_subgraph_map = Dict()
-        optimizer.primal_links = LinkConstraintRef[]
-        optimizer.dual_links = LinkConstraintRef[]
-        optimizer.incident_variable_map = Dict()
-
-        optimizer.x_vals = Vector{Float64}()  #Primal values for communication
-        optimizer.l_vals = Vector{Float64}()  #Dual values for communication
-        optimizer.x_out_vals = Dict{OptiGraph,Vector{Float64}}()
-        optimizer.l_out_vals = Dict{OptiGraph,Vector{Float64}}()
-
-        optimizer.err_pr = Inf
-        optimizer.err_du = Inf
-        optimizer.objective_value = 0
-        optimizer.status = MOI.OPTIMIZE_NOT_CALLED
-        optimizer.iteration = 0
-        optimizer.primal_error_iters = Float64[]
-        optimizer.dual_error_iters = Float64[]
-        optimizer.objective_iters = Float64[]
-
-        #optimizer.plasmo_optimizer_hook = SchwarzOpt.optimize!
-        optimizer.timers = Timers()
-        optimizer.loaded = false
-
-        return optimizer
-    end
-end
-
-function Optimizer(graph::OptiGraph,subgraphs::Vector{OptiGraph};
-    sub_optimizer = nothing,
-    primal_links::Vector{LinkConstraintRef} = LinkConstraintRef[],
-    dual_links::Vector{LinkConstraintRef} = LinkConstraintRef[],
-    tolerance = 1e-6,
-    max_iterations = 100)
-
-    optimizer = Optimizer()
-    for subgraph in subgraphs
-        sub = Plasmo.optigraph_reference(subgraph)
-        push!(optimizer.subproblem_graphs,sub)
-        optimizer.x_out_vals[sub] = Float64[]
-        optimizer.l_out_vals[sub] = Float64[]
-        optimizer.x_in_indices[sub] = Int64[]
-        optimizer.l_in_indices[sub] = Int64[]
-        optimizer.x_out_indices[sub] = Int64[]
-        optimizer.l_out_indices[sub] = Int64[]
-    end
-
-    optimizer.graph = graph
-    # optimizer.subproblem_graphs = subgraphs
-    optimizer.sub_optimizer = sub_optimizer
-    optimizer.primal_links = primal_links
-    optimizer.dual_links = dual_links
-    optimizer.tolerance = tolerance
-    optimizer.max_iterations = max_iterations
-
-    return optimizer
-end
-
-function _check_valid_subgraphs(graph::OptiGraph,subgraphs::Vector{OptiGraph})
-    subgraph_nodes = union(all_nodes.(subgraphs)...)
-    length(subgraph_nodes) == num_all_nodes(graph) || error("Invalid subgraphs given to optimizer.  The number of nodes in subgraphs does not match the number of nodes
-    in the optigraph.")
-    all((node) -> node in all_nodes(graph),subgraph_nodes) || error("Invalid subgraphs given to optimizer.  At least one provided subgraph
-    constrains an optinode that is not in the optigraph.")
-
-    #TODO: check for overlap of at least 1.
-
-    #TODO: check for non-contiguous partitions and raise warning if they are not
-    return true
-end
-
-function _initialize_optimizer!(optimizer::Optimizer)
-    if optimizer.sub_optimizer == nothing
-        error("No optimizer set for the subproblems.  Please provide an optimizer constructor to use to solve subproblem optigraphs")
-    end
-
-    optimizer.timers.initialize_time = @elapsed begin
-        graph = optimizer.graph
-        overlap_subgraphs = optimizer.subproblem_graphs
-        primal_links = optimizer.primal_links
-        dual_links = optimizer.dual_links
-
-        _check_valid_subgraphs(graph,overlap_subgraphs)
-        n_subproblems = length(optimizer.subproblem_graphs)
-        #MAP OPTINODES TO ORIGNAL SUBGRAPHS
-        original_subgraphs = getsubgraphs(graph)
-        for sub in original_subgraphs
-            for node in all_nodes(sub)
-                optimizer.node_subgraph_map[node] = sub
-            end
-        end
-
-        #Map subgraphs to their expanded versions
-        for i = 1:length(optimizer.subproblem_graphs)
-            original_subgraph = original_subgraphs[i]
-            expanded_subgraph = optimizer.subproblem_graphs[i]
-            @assert intersect(all_nodes(original_subgraph),all_nodes(expanded_subgraph)) == all_nodes(original_subgraph)
-            optimizer.expanded_subgraph_map[original_subgraph] = expanded_subgraph #TODO: make sure these match up
-            expanded_subgraph.ext[:restricted_subgraph] = original_subgraph
-            expanded_subgraph.ext[:restricted_objective] = sum(objective_function(node) for node in all_nodes(original_subgraph))
-        end
-
-        #FIND SUBGRAPH BOUNDARIES AND ASSIGN LINKS AS EITHER PRIMAL OR DUAL
-        subgraph_boundary_edges = _find_boundaries(graph,optimizer.subproblem_graphs)
-        primal_links,dual_links = _assign_links(optimizer.subproblem_graphs,subgraph_boundary_edges,primal_links,dual_links)
-        optimizer.primal_links = primal_links
-        optimizer.dual_links = dual_links
-        @assert length(primal_links) == n_subproblems
-        @assert length(dual_links) == n_subproblems
-        ########################################################
-        #INITIALIZE SUBPROBLEM DATA
-        ########################################################
-        for subgraph in optimizer.subproblem_graphs
-            #Primal data
-            subgraph.ext[:x_in_map] = Dict{Int64,VariableRef}()                              #variables into subproblem (these are local copies)
-            subgraph.ext[:x_out_map] = Dict{Int64,VariableRef}()                             #variables out of subproblem
-            subgraph.ext[:incident_variable_map] = Dict{VariableRef,VariableRef}()           #map incident variables to local variables
-
-            #Dual data
-            subgraph.ext[:l_in_map] = Dict{Int64,GenericAffExpr{Float64,VariableRef}}()    #duals for penalty objective term
-            subgraph.ext[:l_out_map] = Dict{Int64,LinkConstraintRef}()                     #duals from linkconstraint
-
-            #Original objective function
-            JuMP.set_objective(subgraph,MOI.MIN_SENSE,sum(objective_function(node) for node in all_nodes(subgraph)))
-            obj = objective_function(subgraph)
-            subgraph.ext[:original_objective] = obj
-        end
-
-        #TODO:pre-allocate vectors based on dual and primal links
-        #inspect to figure out dual links, incident variables, and source and target variables
-        #######################################################
-        #INITIALIZE SUBPROBLEMS
-        #######################################################
-        for i = 1:length(optimizer.subproblem_graphs)
-            subproblem_graph = optimizer.subproblem_graphs[i]
-            subproblem_dual_links = dual_links[i]
-            subproblem_primal_links = primal_links[i]
-
-            #DUAL LINKS
-            #check in/out for each link
-            for link_reference in subproblem_dual_links
-                #INPUTS to subproblem
-                dual_start = 0
-                push!(optimizer.l_vals,dual_start)
-                idx = length(optimizer.l_vals)
-                push!(optimizer.l_in_indices[subproblem_graph],idx)                                #add index to subproblem l inputs
-                _add_subproblem_dual_penalty!(subproblem_graph,link_reference,idx)   #add penalty to subproblem objective
-
-                #push!(subproblem_graph.ext[:l_in],link_reference)
-                #subproblem_graph.ext[:l_in][idx] = link_reference
-
-                #OUTPUTS from subproblem
-                link = constraint_object(link_reference)
-                vars = collect(keys(link.func.terms))                                     #variables in linkconsstraint
-                incident_variables = setdiff(vars,all_variables(subproblem_graph))
-                incident_node = Plasmo.attached_node(link)                                         #this is the owning node for the link
-                @assert !(incident_node in all_nodes(subproblem_graph))
-                original_subgraph = optimizer.node_subgraph_map[incident_node]                        #get the restricted subgraph
-                target_subgraph = optimizer.expanded_subgraph_map[original_subgraph]                     #the subproblem that "owns" this link_constraint
-                push!(optimizer.l_out_indices[target_subgraph],idx)                            #add index to target subproblem outputs
-                target_subgraph.ext[:l_out_map][idx] = link_reference
-            end
-
-            #PRIMAL LINKS
-            #check in/out for each link
-            for link_reference in subproblem_primal_links
-                link = constraint_object(link_reference)
-                vars = collect(keys(link.func.terms))
-
-                local_variables = intersect(vars,all_variables(subproblem_graph))
-                for var in local_variables
-                    subproblem_graph.ext[:incident_variable_map][var] = var
-                end
-
-                incident_variables = setdiff(vars,all_variables(subproblem_graph))
-                for incident_variable in incident_variables
-                    if !(incident_variable in keys(optimizer.incident_variable_map))                  #if incident variable hasn't been counted yet, create a new one
-                        JuMP.start_value(incident_variable) == nothing ? start = 0 : start = JuMP.start_value(incident_variable)
-                        push!(optimizer.x_vals,start)                                                 #increment x_vals
-
-                        idx = length(optimizer.x_vals)                                                #get index
-                        optimizer.incident_variable_map[incident_variable] = idx                      #map index for external variable
-
-                        incident_node = getnode(incident_variable)                                    #get the node for this external variable
-                        original_subgraph = optimizer.node_subgraph_map[incident_node]                #get the restricted subgraph
-                        source_subgraph = optimizer.expanded_subgraph_map[original_subgraph]                    #get the subproblem that owns this external variable
-
-                        #OUTPUTS
-                        push!(optimizer.x_out_indices[source_subgraph],idx)                           #add index to source subproblem outputs
-                        #push!(source_subgraph.ext[:x_out],incident_variable)                          #map incident variable to source problem primal outputs
-                        source_subgraph.ext[:x_out_map][idx] = incident_variable
-                    else
-                        idx = optimizer.incident_variable_map[incident_variable]
-                    end
-                    #If this subproblem needs to make a local copy of the incident variable
-                    if !(incident_variable in keys(subproblem_graph.ext[:incident_variable_map]))
-                        _add_subproblem_variable!(subproblem_graph,incident_variable,idx)
-                        push!(optimizer.x_in_indices[subproblem_graph],idx)
-                    end
-                end
-                _add_subproblem_constraint!(subproblem_graph,link_reference)                                    #Add link constraint to the subproblem. This problem "owns" the constraint
-            end
-        end
-        optimizer.loaded = true
-    end
-end
-
-function _add_subproblem_variable!(subproblem_graph::OptiGraph,incident_variable::VariableRef,idx::Int64)
-    JuMP.start_value(incident_variable) == nothing ? start = 0 : start = JuMP.start_value(incident_variable)
-    copy_node = @optinode(subproblem_graph)
-    copy_variable = @variable(copy_node,start = start)
-    JuMP.set_name(copy_variable,name(incident_variable)*"_copy")
-    if JuMP.has_lower_bound(incident_variable)
-        JuMP.set_lower_bound(copy_variable,lower_bound(incident_variable))
-    end
-    if JuMP.has_upper_bound(incident_variable)
-        JuMP.set_upper_bound(copy_variable,upper_bound(incident_variable))
-    end
-
-    #JuMP.fix(copy_variable,start)
-    subproblem_graph.ext[:incident_variable_map][incident_variable] = copy_variable
-    subproblem_graph.ext[:x_in_map][idx] = copy_variable
-    return nothing
-end
-
-#Add link constraint to optigraph connecting to copy variable
-function _add_subproblem_constraint!(subproblem_graph::OptiGraph,link_reference::LinkConstraintRef)
-    varmap = subproblem_graph.ext[:incident_variable_map]
-
-    agg_map = Plasmo.AggregateMap(varmap, Dict{Plasmo.ConstraintRef,Plasmo.ConstraintRef}(), Dict{Plasmo.LinkConstraint, Plasmo.ConstraintRef}())
-
-    link = constraint_object(link_reference)
-    copy_link = Plasmo._copy_constraint(link, agg_map)
-    copy_linkref = JuMP.add_constraint(subproblem_graph, copy_link) #this is a linkconstraint
-    #push!(graph.ext[:incident_constraints], copy_linkref) #this isn't used anywhere?
-    return nothing
-end
-
-#Add penalty to subproblem objective
-function _add_subproblem_dual_penalty!(subproblem_graph::OptiGraph,link_reference::LinkConstraintRef,idx::Int64)
-    link = constraint_object(link_reference)
-    link_variables = collect(keys(link.func.terms))
-    local_link_variables = intersect(link_variables, all_variables(subproblem_graph))
-    penalty = sum(local_link_variables)
-    set_objective_function(subproblem_graph,objective_function(subproblem_graph) - penalty)
-    subproblem_graph.ext[:l_in_map][idx] = penalty
-    return nothing
-end
-
-function _do_iteration(subproblem_graph::OptiGraph)
-    Plasmo.optimize!(subproblem_graph)
-    term_status = termination_status(subproblem_graph)
-    !(term_status in [MOI.TerminationStatusCode(4),MOI.TerminationStatusCode(1),MOI.TerminationStatusCode(10)]) && @warn("Suboptimal solution detected for subproblem with status $term_status")
-    #has_values(subproblem_graph) || error("Could not obtain values for problem $subproblem_graph with status $term_status")
-
-    x_out = subproblem_graph.ext[:x_out_map]           #primal variables to communicate
-    l_out = subproblem_graph.ext[:l_out_map]           #dual variables (on linkconstraints) to communicate
-
-    xk = Dict(key => value(subproblem_graph,val) for (key,val) in x_out)
-    lk = Dict(key => dual(subproblem_graph,val) for (key,val) in l_out)
-
-
-    return xk, lk
-end
-
-function _update_subproblem!(subproblem_graph::OptiGraph,x_in_vals::Vector{Float64},l_in_vals::Vector{Float64},x_in_inds::Vector{Int64},l_in_inds::Vector{Int64})
-    #Fix primal inputs
-    @assert length(x_in_vals) == length(x_in_inds)
-    @assert length(l_in_vals) == length(l_in_inds)
-    for (i,idx) in enumerate(x_in_inds)
-        variable = subproblem_graph.ext[:x_in_map][idx]
-        JuMP.fix(variable, x_in_vals[i]; force = true) #fix variable for this subproblem.  variable should be the copy made for this subgraph
-    end
-    #This is slow
-    for (i,idx) in enumerate(l_in_inds)
-        penalty = subproblem_graph.ext[:l_in_map][idx]
-        #TODO: just update objective function coefficient in the backend.
-        JuMP.set_objective_function(subproblem_graph,subproblem_graph.ext[:original_objective] - l_in_vals[i]*penalty)
-    end
-    return nothing
-end
-
-function _calculate_objective_value(optimizer)
-    return sum(value(optimizer.subproblem_graphs[i].ext[:restricted_objective]) for i = 1:length(optimizer.subproblem_graphs))
-end
-
-#TODO: figure out these mappings ahead of time.  This is bottlenecking
-function _calculate_primal_feasibility(optimizer)
-    linkrefs = getlinkconstraints(optimizer.graph)
-    prf = []
-    for linkref in linkrefs
-        val = 0
-        linkcon = constraint_object(linkref)
-        terms = linkcon.func.terms
-        for (term,coeff) in terms
-            node = getnode(term)
-            graph = optimizer.node_subgraph_map[node]
-            subproblem_graph = optimizer.expanded_subgraph_map[graph]
-            val += coeff*value(subproblem_graph, term)
-        end
-        push!(prf, val - linkcon.set.value)
-    end
-    return prf
-end
-
-#NOTE: Need at least an overlap of one
-function _calculate_dual_feasibility(optimizer)
-    linkrefs = getlinkconstraints(optimizer.graph)
-    duf = []
-    for linkref in linkrefs
-        lambdas = []
-        linkcon = constraint_object(linkref)
-        graphs = []
-        for node in getnodes(linkcon)
-            graph = optimizer.node_subgraph_map[node]
-            subproblem_graph = optimizer.expanded_subgraph_map[graph]
-            push!(graphs,subproblem_graph)
-        end
-        graphs = unique(graphs)
-        for subproblem_graph in graphs
-            l_val = dual(subproblem_graph,linkref)   #check each subproblem's dual value for this linkconstraint
-            push!(lambdas,l_val)
-        end
-        #@assert length(lambdas) == 2
-        dual_res = maximum(diff(lambdas))#lambdas[1] - lambdas[2]
-        push!(duf,dual_res)
-    end
-    return duf
-end
-
-function optimize!(optimizer::Optimizer)
-    if !optimizer.loaded
-        println("Initializing SchwarzOpt...")
-        _initialize_optimizer!(optimizer)
-    end
-
-    println("###########################################################")
-    println("Optimizing with SchwarzOpt v0.1.0 using $(Threads.nthreads()) threads")
-    println("###########################################################")
-    println()
-    println("Number of variables: $(num_all_variables(optimizer.graph))")
-    println("Number of constraints: $(num_all_constraints(optimizer.graph) + num_all_linkconstraints(optimizer.graph))")
-    println("Number of subproblems: $(length(optimizer.subproblem_graphs))")
-    println("Overlap: ")
-    println("Subproblem variables:   $(num_all_variables.(optimizer.subproblem_graphs))")
-    println("Subproblem constraints: $(num_all_constraints.(optimizer.subproblem_graphs))")
-    println()
-
-    optimizer.timers = Timers()
-    optimizer.timers.start_time = time()
-
-    for subgraph in optimizer.subproblem_graphs
-        JuMP.set_optimizer(subgraph,optimizer.sub_optimizer)
-    end
-
-    optimizer.err_pr = Inf
-    optimizer.err_du = Inf
-    optimizer.iteration = 0
-    while optimizer.err_pr > optimizer.tolerance || optimizer.err_du > optimizer.tolerance
-        optimizer.iteration += 1
-        if optimizer.iteration > optimizer.max_iterations
-            optimizer.status = MOI.ITERATION_LIMIT
-            break
-        end
-
-        #Do iteration for each subproblem
-        optimizer.timers.update_subproblem_time += @elapsed begin
-            if optimizer.iteration > 1 #don't fix variables in first iteration
-                for subproblem_graph in optimizer.subproblem_graphs
-                    x_in_inds = optimizer.x_in_indices[subproblem_graph]
-                    l_in_inds = optimizer.l_in_indices[subproblem_graph]
-                    x_in_vals = optimizer.x_vals[x_in_inds]
-                    l_in_vals = optimizer.l_vals[l_in_inds]
-                    _update_subproblem!(subproblem_graph,x_in_vals,l_in_vals,x_in_inds,l_in_inds)
-                end
-            end
-        end
-
-        optimizer.timers.solve_subproblem_time += @elapsed begin
-            Threads.@threads for subproblem_graph in optimizer.subproblem_graphs
-                #Returns primal and dual information we need to communicate to other subproblems
-                xk,lk = _do_iteration(subproblem_graph)
-
-                #Updates primal and dual information for other subproblems.
-                for (idx,val) in xk
-                    optimizer.x_vals[idx] = val
-                end
-                for (idx,val) in lk
-                    optimizer.l_vals[idx] = val
-                end
-            end
-        end
-
-        #Evaluate residuals
-        optimizer.timers.eval_primal_feasibility_time += @elapsed prf = _calculate_primal_feasibility(optimizer)
-        optimizer.timers.eval_dual_feasibility_time += @elapsed duf = _calculate_dual_feasibility(optimizer)
-        optimizer.err_pr = norm(prf[:],Inf)
-        optimizer.err_du = norm(duf[:],Inf)
-
-
-        #NOTE: This is the wrong way to calculate this
-        #TODO: Calculate objective value correctly
-
-        optimizer.timers.eval_objective_time += @elapsed optimizer.objective_value = _calculate_objective_value(optimizer)
-        push!(optimizer.primal_error_iters,optimizer.err_pr)
-        push!(optimizer.dual_error_iters,optimizer.err_du)
-        push!(optimizer.objective_iters,optimizer.objective_value)
-
-        #Print iteration
-        if optimizer.iteration % 20 == 0 || optimizer.iteration == 1
-            @printf "%4s | %8s | %8s | %8s" "Iter" "Obj" "Prf" "Duf\n"
-        end
-        @printf("%4i | %7.2e | %7.2e | %7.2e\n",optimizer.iteration,optimizer.objective_value,optimizer.err_pr,optimizer.err_du)
-
-        optimizer.timers.update_subproblem_time += @elapsed begin
-            for subproblem in optimizer.subproblem_graphs
-                JuMP.set_start_value.(Ref(subproblem),all_variables(subproblem),value.(Ref(subproblem),all_variables(subproblem)))
-            end
-        end
-    end
-
-    #Point variables to restricted solutions
-    for (node,subgraph) in optimizer.node_subgraph_map
-        subproblem_graph = optimizer.expanded_subgraph_map[subgraph]
-        backend(node).last_solution_id = subproblem_graph.id
-    end
-
-    optimizer.timers.total_time =  time() - optimizer.timers.start_time
-    if optimizer.status != MOI.ITERATION_LIMIT
-        optimizer.status = termination_status(optimizer.subproblem_graphs[1])
-    end
-    # optimizer.graph.optimizer = optimizer
-    println()
-    println("Number of Iterations: ",length(optimizer.iteration))
-
-    println("Time spent in subproblems: ",optimizer.timers.solve_subproblem_time)
-    println("Solution Time: ",optimizer.timers.total_time)
-    println("EXIT: SchwarzOpt Finished with status: ",optimizer.status)
-end
-
-#call with SchwarzOpt.optimize!(graph)
-"""
-    SchwarzOpt.optimize!(graph::OptiGraph;kwargs...)
-
-Optimize an optigraph with overlapping schwarz decomposition.
-"""
-function optimize!(graph::OptiGraph;
-    subgraphs = Plasmo.OptiGraph[],
-    sub_optimizer = Ipopt.Optimizer,
-    overlap = 1,
-    max_iterations = 50,
-    primal_links = LinkConstraintRef[],
-    dual_links = LinkConstraintRef[])
-
-    #TODO: heuristic to calculate overlap
-    # Plasmo.graph_structure(graph) == Plasmo.RECURSIVE_GRAPH || error("Invalid optigraph.  SchwarzOpt requires a RECURSIVE_GRAPH_STRUCTURE.  Consider using partition_to_subgraphs to obtain
-    # a valid optigraph.")
-    #check for subgraphs, or apply overlap if no subgraphs defined
-    if length(subgraphs) == 0
-        println("Optimizing with overlap of $overlap")
-        subgraphs = expand(graph,subgraphs,overlap)
-    else
-        println("Optimizing with user provided overlap")
-    end
-
-    optimizer = Optimizer(
-        graph,
-        subgraphs; 
-        primal_links = primal_links, 
-        dual_links = dual_links,
-        sub_optimizer = sub_optimizer, 
-        max_iterations = max_iterations
-    )
-
-    #set the optimizer for querying results on the optigraph
-    # graph.optimizer = optimizer
-    optimize!(optimizer)
-end
diff --git a/src/optimizer2.jl b/src/optimizer2.jl
index 90d9d17..8a4783b 100644
--- a/src/optimizer2.jl
+++ b/src/optimizer2.jl
@@ -20,6 +20,7 @@ mutable struct Optimizer <: MOI.AbstractOptimizer
     node_subgraph_map::Dict                       # map optinodes to "owning" subgraph
     expanded_subgraph_map::Dict
     incident_variable_map::Dict                   # map variable indices in optigraphs to ghost (copy) variables
+    incident_constraint_map::Dict
     subproblem_links::Vector                      # link constraints that do dual updates
 
     x_vals::Vector{Float64}                       # current primal values that get communicated to subproblems
@@ -47,7 +48,9 @@ mutable struct Optimizer <: MOI.AbstractOptimizer
         optimizer.node_subgraph_map = Dict()
         optimizer.expanded_subgraph_map = Dict()
         optimizer.subproblem_links = Vector{Vector{LinkConstraintRef}}()
-        optimizer.incident_variable_map = Dict()
+
+        optimizer.incident_variable_map = OrderedDict()
+        optimizer.incident_constraint_map = OrderedDict()
 
         optimizer.x_vals = Vector{Float64}()  #Primal values for communication
         optimizer.l_vals = Vector{Float64}()  #Dual values for communication
@@ -137,7 +140,7 @@ function _initialize_optimizer!(optimizer::Optimizer)
         end
 
         # gather the boundary links for each subproblem
-        subgraph_boundary_edges = _find_boundaries(graph,optimizer.subproblem_graphs)
+        subgraph_boundary_edges = _find_boundaries(graph, optimizer.subproblem_graphs)
         subproblem_links = _gather_links(optimizer.subproblem_graphs, subgraph_boundary_edges)
         optimizer.subproblem_links = subproblem_links
         @assert length(subproblem_links) == n_subproblems
@@ -145,12 +148,13 @@ function _initialize_optimizer!(optimizer::Optimizer)
         ########################################################
         #INITIALIZE SUBPROBLEM DATA
         ########################################################
-        for subgraph in optimizer.subproblem_graphs
+        for (i,subgraph) in enumerate(optimizer.subproblem_graphs)
             # initialize primal-dual maps
             subgraph.ext[:x_map] = Dict{Int64,VariableRef}()
             subgraph.ext[:x_idx_map] = Dict{VariableRef,Int64}()
             subgraph.ext[:l_map] = Dict{Int64,LinkConstraintRef}()
             subgraph.ext[:l_idx_map] = Dict{LinkConstraintRef,Int64}()
+            subgraph.ext[:incident_links] = optimizer.subproblem_links[i]
 
             # set original subproblem objective function
             JuMP.set_objective(subgraph, MOI.MIN_SENSE, sum(objective_function(node) for node in all_nodes(subgraph)))
@@ -163,76 +167,85 @@ function _initialize_optimizer!(optimizer::Optimizer)
         #######################################################
         for i = 1:length(optimizer.subproblem_graphs)
             subproblem_graph = optimizer.subproblem_graphs[i]
-            sub_links = subproblem_links[i]
+            sub_links = subproblem_graph.ext[:incident_links]
 
-            #map all of the boundary variables for each subproblem
+            #map each subproblem to the primals and duals it is responsible for updating
             for link_reference in sub_links
                 # build primal information map
                 link = constraint_object(link_reference)
                 vars = collect(keys(link.func.terms))
-
                 incident_variables = setdiff(vars, all_variables(subproblem_graph))
+
                 for incident_variable in incident_variables
-                    if !(incident_variable in keys(optimizer.incident_variable_map))        # if incident variable hasn't been counted yet, create a new one
+                    # if incident variable hasn't been counted yet, add to array
+                    if !(incident_variable in keys(optimizer.incident_variable_map))        
                         JuMP.start_value(incident_variable) == nothing ? start = 0 : start = JuMP.start_value(incident_variable)
-                        push!(optimizer.x_vals, start)                                                #increment x_vals
-                        idx = length(optimizer.x_vals) 
-                        optimizer.incident_variable_map[incident_variable] = idx                      #map index for external variable
+                        push!(optimizer.x_vals, start)               
+                        idx = length(optimizer.x_vals)
+                        optimizer.incident_variable_map[incident_variable] = idx
 
-                        incident_node = getnode(incident_variable)                                    #get the node for this external variable
-                        original_subgraph = optimizer.node_subgraph_map[incident_node]                #get the restricted subgraph
-                        exp_subgraph = optimizer.expanded_subgraph_map[original_subgraph]             #get the subproblem that owns this external variable
+                        incident_node = getnode(incident_variable)                                    # get the node for this external variable
+                        original_subgraph = optimizer.node_subgraph_map[incident_node]                # get the restricted subgraph
+                        exp_subgraph = optimizer.expanded_subgraph_map[original_subgraph]             # get the subproblem that owns this external variable
                         exp_subgraph.ext[:x_idx_map][incident_variable] = idx
-                        exp_subgraph.ext[:x_map][idx] = incident_variable                                  # this graph is used to calculate this variable
+                        exp_subgraph.ext[:x_map][idx] = incident_variable                             # this graph is used to calculate this variable
                     end
                 end
 
                 # build dual information map
-                dual_start = 0
-                push!(optimizer.l_vals, dual_start)
-                idx = length(optimizer.l_vals)
-                
-                link = constraint_object(link_reference)
-                vars = collect(keys(link.func.terms))                                     
-
-                # we use the attached node to determine which subproblem updates the dual
-                incident_node = Plasmo.attached_node(link)                                         #this is the owning node for the link
-                
-                # NOTE: I don't think it matters if the attached node is in the subproblem, we just grab the 'owning' subgraph regardless
-                # @assert !(incident_node in all_nodes(subproblem_graph))
-                original_subgraph = optimizer.node_subgraph_map[incident_node]                     #get the restricted subgraph
-                exp_subgraph = optimizer.expanded_subgraph_map[original_subgraph]                  #the subproblem that "owns" this link_constraint
-                exp_subgraph.ext[:l_idx_map][link_reference] = idx
-                exp_subgraph.ext[:l_map][idx] = link_reference                                     # this graph is used to calculate this dual
+                if !(link_reference in keys(optimizer.incident_constraint_map))
+                    dual_start = 0
+                    push!(optimizer.l_vals, dual_start)
+                    idx = length(optimizer.l_vals)
+                    optimizer.incident_constraint_map[link_reference] = idx
+                    
+                    link = constraint_object(link_reference)
+                    vars = collect(keys(link.func.terms))                                     
+
+                    # we use the attached node to determine which subproblem updates the dual
+                    incident_node = Plasmo.attached_node(link)                                         #this is the owning node for the link
+                    
+                    # NOTE: I don't think it matters if the attached node is in the subproblem, we just grab the 'owning' subgraph regardless
+                    # @assert !(incident_node in all_nodes(subproblem_graph))
+                    original_subgraph = optimizer.node_subgraph_map[incident_node]                     #get the restricted subgraph
+                    exp_subgraph = optimizer.expanded_subgraph_map[original_subgraph]                  #the subproblem that "owns" this link_constraint
+                    exp_subgraph.ext[:l_idx_map][link_reference] = idx
+                    exp_subgraph.ext[:l_map][idx] = link_reference                                     # this graph is used to calculate this dual
+                end
             end
         end
         optimizer.loaded = true
     end
 end
 
-#Add penalty to subproblem objective
+function _update_subproblem!(optimizer, subproblem_graph::OptiGraph)
+    for link in subproblem_graph.ext[:incident_links]
+        _formulate_penalty!(optimizer, subproblem_graph, link)
+    end
+    return nothing
+end
+
 function _formulate_penalty!(optimizer, subproblem_graph::OptiGraph, link_reference::LinkConstraintRef)
     link = constraint_object(link_reference)
     link_variables = collect(keys(link.func.terms))
-    
     local_link_variables = intersect(link_variables, all_variables(subproblem_graph))
-    
+
     external_variables = setdiff(link_variables, local_link_variables)
-    external_var_inds = [subproblem_graph.ext[:x_idx_map][var] for var in external_variables]
+    external_var_inds = [optimizer.incident_variable_map[var] for var in external_variables]
     external_values = [optimizer.x_vals[x_idx] for x_idx in external_var_inds]
 
-    local_link_coeffs = [link.func.terms[var] for var in link_variables]
+    local_link_coeffs = [link.func.terms[var] for var in local_link_variables]
     external_coeffs = [link.func.terms[var] for var in external_variables]
 
-    l_idx = subproblem_graph.ext[:l_idx_map][link_reference]
+    l_idx = optimizer.incident_constraint_map[link_reference]
     l_link = optimizer.l_vals[l_idx]
 
-    rhs = link.func.set.value
+    rhs = link.set.value
     penalty = local_link_coeffs'*local_link_variables + external_coeffs'*external_values - rhs
     augmented_penalty = penalty^2
 
     # TODO: add to objective instead of reset
-    set_objective_function(subproblem_graph, objective_function(subproblem_graph) - l_link*penalty + 0.5*optimizer.mu*augmented_penalty)
+    set_objective_function(subproblem_graph, subproblem_graph.ext[:original_objective] - l_link*penalty + 0.5*optimizer.mu*augmented_penalty)
     return nothing
 end
 
@@ -250,13 +263,6 @@ function _do_iteration(subproblem_graph::OptiGraph)
     return xk, lk
 end
 
-function _update_subproblem!(optimizer, subproblem_graph::OptiGraph)
-    for link in subproblem_graph.ext[:incident_links]
-        _formulate_penalty!(optimizer, subproblem_graph, link)
-    end
-    return nothing
-end
-
 function _calculate_objective_value(optimizer)
     return sum(value(optimizer.subproblem_graphs[i].ext[:restricted_objective]) for i = 1:length(optimizer.subproblem_graphs))
 end
@@ -303,7 +309,7 @@ function _calculate_dual_feasibility(optimizer)
         end
 
         # dual residual between subproblems
-        dual_res = maximum(diff(lambdas))#lambdas[1] - lambdas[2]
+        dual_res = maximum(diff(duals))#lambdas[1] - lambdas[2]
         push!(duf, dual_res)
     end
     return duf
@@ -346,14 +352,8 @@ function optimize!(optimizer::Optimizer)
 
         #Do iteration for each subproblem
         optimizer.timers.update_subproblem_time += @elapsed begin
-            if optimizer.iteration > 1 #don't fix variables in first iteration
-                for subproblem_graph in optimizer.subproblem_graphs
-                    x_in_inds = optimizer.x_in_indices[subproblem_graph]
-                    l_in_inds = optimizer.l_in_indices[subproblem_graph]
-                    x_in_vals = optimizer.x_vals[x_in_inds]
-                    l_in_vals = optimizer.l_vals[l_in_inds]
-                    _update_subproblem!(subproblem_graph, x_in_vals, l_in_vals, x_in_inds, l_in_inds)
-                end
+            for subproblem_graph in optimizer.subproblem_graphs
+                _update_subproblem!(optimizer, subproblem_graph)
             end
         end
 
@@ -392,7 +392,7 @@ function optimize!(optimizer::Optimizer)
 
         optimizer.timers.update_subproblem_time += @elapsed begin
             for subproblem in optimizer.subproblem_graphs
-                JuMP.set_start_value.(Ref(subproblem),all_variables(subproblem),value.(Ref(subproblem),all_variables(subproblem)))
+                JuMP.set_start_value.(Ref(subproblem), all_variables(subproblem), value.(Ref(subproblem),all_variables(subproblem)))
             end
         end
     end
diff --git a/src/utils.jl b/src/utils.jl
deleted file mode 100644
index ae04425..0000000
--- a/src/utils.jl
+++ /dev/null
@@ -1,78 +0,0 @@
-#TODO: Update internal expand
-function _expand_subgraphs(graph::OptiGraph,overlap::Int64)
-    subproblem_graphs = []
-    boundary_linkedges_list = []
-    hypergraph,hyper_map = gethypergraph(graph)
-
-    #NOTE: This could all be calculated simultaneously
-    for subgraph in getsubgraphs(mg)
-
-        println("Performing neighborhood expansion...")
-        subnodes = all_nodes(subgraph)
-        hypernodes = [hyper_map[node] for node in subnodes]
-
-        #Return new hypernodes and hyperedges covered through the expansion.
-        overlap_hnodes = Plasmo.neighborhood(hypergraph,hypernodes,overlap)
-        overlap_hedges = Plasmo.induced_edges(hypergraph,overlap_hnodes)
-        boundary_hedges = Plasmo.incident_edges(hypergraph,overlap_hnodes)
-
-        overlap_nodes = [hyper_map[node] for node in overlap_hnodes]
-        overlap_edges = [hyper_map[edge] for edge in overlap_hedges]
-        boundary_edges = [hyper_map[edge] for edge in boundary_hedges]
-
-        #Setup subproblem graphs
-        subproblem_graph = OptiGraph()
-        subproblem_graph.modelnodes = overlap_nodes
-        subproblem_graph.linkedges = overlap_edges
-
-        push!(subproblem_graphs,subproblem_graph)
-        push!(boundary_linkedges_list,boundary_edges)
-    end
-
-    return subproblem_graphs,boundary_linkedges_list
-end
-
-#Links can be formulated as either constraints or penalties
-function _assign_links(subgraphs,subgraph_boundary_edges,input_primal_links,input_dual_links)
-    subgraph_primal_links = []
-    subgraph_dual_links = []
-
-    for (i,edge_set) in enumerate(subgraph_boundary_edges)
-        primal_links = LinkConstraintRef[]
-        dual_links = LinkConstraintRef[]
-
-        for edge in edge_set
-            linkrefs = edge.linkrefs
-            for linkref in linkrefs
-                if linkref in input_primal_links
-                    push!(primal_links,linkref)
-                elseif linkref in input_dual_links
-                    push!(dual_links,linkref)
-                else
-                    #use attached node to assign link
-                    target_node = constraint_object(linkref).attached_node
-                    #target_node = edge.nodes[end]
-                    if !(target_node in all_nodes(subgraphs[i]))
-                        push!(dual_links,linkref)  #send primal info to target node, receive dual info
-                    else
-                        push!(primal_links,linkref)   #receive primal info at target node, send back dual info
-                    end
-                end
-            end
-        end
-        push!(subgraph_primal_links,primal_links)
-        push!(subgraph_dual_links,dual_links)
-    end
-    return subgraph_primal_links,subgraph_dual_links
-end
-
-#Find boundary edges of expanded subgraphs
-function _find_boundaries(optigraph::OptiGraph,subgraphs::Vector{OptiGraph})
-    boundary_linkedges_list = []
-    for subgraph in subgraphs
-        subnodes = all_nodes(subgraph)
-        boundary_edges = Plasmo.incident_edges(optigraph,subnodes)
-        push!(boundary_linkedges_list,boundary_edges)
-    end
-    return boundary_linkedges_list
-end