Support P-Split Reformulation for (affine and quadratic constraints only)

Project.toml

@@ -12,11 +12,15 @@ Aqua = "0.8"
 JuMP = "1.18"
 Reexport = "1"
 julia = "1.6"
+Juniper = "0.9.3"
+Ipopt = "1.9.0"
 
 [extras]
 Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595"
 HiGHS = "87dc4568-4c63-4d18-b0c0-bb2238e4078b"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
+Ipopt = "b6b21f68-93f8-5de0-b562-5493be1d77c9"
+Juniper = "2ddba703-00a4-53a7-87a5-e8b9971dde84"
 
 [targets]
-test = ["Aqua", "HiGHS", "Test"]
+test = ["Aqua", "HiGHS", "Test", "Juniper", "Ipopt"]

README.md

@@ -177,6 +177,14 @@ The following reformulation methods are currently supported:
     - `optimizer`: Optimizer to use when solving subproblems to determine M values. This is a required value.
     - `default_M`: Default big-M value to use if no big-M is specified for a logical variable (1e9).
 
+5. [P-Split](https://arxiv.org/abs/2202.05198): This method reformulates each disjunct constraint into P constraints, each with a partitioned group defined by the user. This method requires that terms in the constraint be convex additively seperable with respect to each variable. The `PSplit` struct is created with the following required arguments:
+
+    - `partition`: Partition of the variables to be split. All variables must be in exactly one partition. (e.g., The variables `x[1:4]` can be partitioned into two groups ` partition = [[x[1], x[2]], [x[3], x[4]]]`)
+    - `PSplit(n_parts, model)`: Automatically partition all variables in the model into `n_parts` groups
+
+    All variables must be included in exactly one partition. For manual partitioning, ensure each variable appears in exactly one group. For automatic partitioning, variables are divided as evenly as possible among the specified number of partitions.
+
+
 ## Release Notes
 
 Prior to `v0.4.0`, the package did not leverage the JuMP extension capabilities and was not as robust. For these earlier releases, refer to [Perez, Joshi, and Grossmann, 2023](https://arxiv.org/abs/2304.10492v1) and the following [JuliaCon 2022 Talk](https://www.youtube.com/watch?v=AMIrgTTfUkI).
@@ -192,7 +200,7 @@ using HiGHS
 m = GDPModel(HiGHS.Optimizer)
 @variable(m, 0 ≤ x[1:2] ≤ 20)
 @variable(m, Y[1:2], Logical)
-@constraint(m, [i = 1:2], [2,5][i] ≤ x[i] ≤ [6,9][i], Disjunct(Y[1]))
+@constraint(m, [i = 1:2], [2,5][i] ≤ x[i] ≤ [6,9][i], Disjunct(Y[1])) 
 @constraint(m, [i = 1:2], [8,10][i] ≤ x[i] ≤ [11,15][i], Disjunct(Y[2]))
 @disjunction(m, Y)
 @objective(m, Max, sum(x))

src/DisjunctiveProgramming.jl

@@ -25,6 +25,7 @@ include("mbm.jl")
 include("indicator.jl")
 include("print.jl")
 include("utilities.jl")
+include("psplit.jl")
 
 # Define additional stuff that should not be exported
 const _EXCLUDE_SYMBOLS = [Symbol(@__MODULE__), :eval, :include]

src/datatypes.jl

@@ -431,6 +431,72 @@ mutable struct _Hull{V <: JuMP.AbstractVariableRef, T} <: AbstractReformulationM
     end
 end
 
+"""
+    PSplit <: AbstractReformulationMethod
+
+A type for using the P-split reformulation approach for disjunctive constraints.
+This method partitions variables into groups and handles each group separately.
+
+# Constructors
+- `PSplit(partition::Vector{Vector{V}})`: Create a PSplit with the given 
+partition of variables
+- `PSplit(n_parts::Int, model::JuMP.AbstractModel)`: Automatically partition 
+model variables into `n_parts` groups
+
+# Fields
+- `partition::Vector{Vector{V}}`: The partition of variables, where each inner 
+vector represents a group of variables that will be handled together
+"""
+struct PSplit{V <: JuMP.AbstractVariableRef} <: AbstractReformulationMethod
+    partition::Vector{Vector{V}}
+
+    function PSplit(partition::Vector{Vector{V}}) where 
+        {V <: JuMP.AbstractVariableRef}
+        new{V}(partition)
+    end
+
+    function PSplit(n_parts::Int, model::JuMP.AbstractModel)
+        n_parts > 0 || error("Number of partitions must be 
+        positive, got $n_parts")
+        variables = collect(JuMP.all_variables(model))
+        n_vars = length(variables)
+
+        n_parts = min(n_parts, n_vars)
+        n_parts > 0 || error("No variables found in the model")
+
+        base_size = n_vars ÷ n_parts
+        remaining = n_vars % n_parts
+
+        partition = Vector{Vector{eltype(variables)}}()
+        start_idx = 1
+
+        for i in 1:n_parts
+            part_size = i <= remaining ? base_size + 1 : base_size
+            end_idx = start_idx + part_size - 1
+            push!(partition, variables[start_idx:end_idx])
+            start_idx = end_idx + 1
+        end
+
+        return PSplit(partition)
+    end
+end
+
+# temp struct to store variable disaggregations (reset for each disjunction)
+mutable struct _PSplit{V <: JuMP.AbstractVariableRef, M <: JuMP.AbstractModel, T} <: AbstractReformulationMethod
+    partition::Vector{Vector{V}}
+    sum_constraints::Dict{LogicalVariableRef{M}, Vector{<:AbstractConstraint}}
+    hull::_Hull{V, T}
+    function _PSplit(method::PSplit{V}, model::M) where 
+        {V <: JuMP.AbstractVariableRef, M <: JuMP.AbstractModel}
+        T = JuMP.value_type(M)
+        new{V, M, T}(
+            method.partition, 
+            Dict{LogicalVariableRef{M}, Vector{<:AbstractConstraint}}(), 
+            _Hull(Hull(), Set{V}())
+        )
+    end
+end
+
 """
     Indicator <: AbstractReformulationMethod
 
@@ -539,4 +605,4 @@ function VariableProperties(vref::JuMP.GenericVariableRef{T}) where T
     name = JuMP.name(vref)
     set = JuMP.is_variable_in_set(vref) ? JuMP.moi_set(JuMP.constraint_object(JuMP.VariableInSetRef(vref))) : nothing
     return VariableProperties(info, name, set, nothing)
-end
+end

src/hull.jl

@@ -22,14 +22,14 @@ function _disaggregate_variables(
     end
 end
 function _disaggregate_variable(
-    model::JuMP.AbstractModel, 
+    model::M, 
     lvref::LogicalVariableRef, 
     vref::JuMP.AbstractVariableRef, 
     method::_Hull
-    )
+    ) where {M <: JuMP.AbstractModel}
     #create disaggregated vref
     lb, ub = variable_bound_info(vref)
-    T = JuMP.value_type(typeof(model))
+    T = JuMP.value_type(M)
     info = JuMP.VariableInfo(
         true,      # has_lb = true
         lb,        # lower_bound = lb