EAGO's Branch and Bound Routine

This component is meant to provide a flexible framework for implementing spatial branch-and-bound based optimization routines in Julia. All components of the branch-and-bound routine can be customized by the individual user: lower-bounding problem, upper-bounding problem.

Branch and Bound Node Storage

EAGO.NodeBBType
struct NodeBB

Store information associated with each node in the branch-and-bound tree.

  • lower_variable_bounds::Vector{Float64}: Lower bounds of variable box.

  • upper_variable_bounds::Vector{Float64}: Upper bounds of variable box.

  • is_integer::BitVector: Is dimension integer valued

  • continuous::Bool: Are all dimensions continuous (or fixed)

  • lower_bound::Float64: Lower bound of problem solution on nodeBB

  • upper_bound::Float64: Upper bound of problem solution on nodeBB

  • depth::Int64: Depth of node in B&B tree.

  • cont_depth::Int64: Depth of first parent in B&B tree that was continuously valued

  • id::Int64: Unique ID for each node.

  • branch_direction::EAGO.BranchDirection: Whether last branch was negative or positive in direction

  • last_branch::Int64: Dimension of last branch

  • branch_extent::Float64: Extent of last branch (using for psuedocost calculation)

source

The GlobalOptimizer structure holds all information relevant to branch-and-bound.

EAGO.GlobalOptimizerType
mutable struct GlobalOptimizer{Q, S, T<:ExtensionType} <: MathOptInterface.AbstractOptimizer

Optimizer internal to EAGO which holds information used to perform branch-and-bound in order to solve nonconvex MINLPs.

Descriptions of all fields available in extended help.

Extended Help

  • _subsolvers::SubSolvers{Q, S} where {Q, S}: Storage for relaxed and upper optimizers to use, and any custom extensions

  • _parameters::EAGOParameters: Parameters that do not change during a global solve

  • _input_problem::InputProblem: Expressions and constraints added to the EAGO model (not directly used for relaxations)

  • _working_problem::ParsedProblem: Expressions and problem descriptions that EAGO uses to formulate relaxed problems

  • _auxiliary_variable_info::Union{Nothing, EAGO._AuxVarData}: Information on any auxiliary variables

  • obbt_variable_values::Vector{Bool}: Variables to perform OBBT on (default: all variables in nonlinear expressions)

  • enable_optimize_hook::Bool: Specifies that the optimize_hook! function should be called rather than throw the problem to the standard routine

  • ext::Any: (Deprecated, use _subsolvers instead) Storage for custom extension types

  • _end_state::EAGO.GlobalEndState: The completion status code for the branch-and-bound algorithm

  • _termination_status_code::MathOptInterface.TerminationStatusCode: The MathOptInterface-compliant completion status code

  • _result_status_code::MathOptInterface.ResultStatusCode: Value indicating the feasibility status of the result

  • _obj_mult::Float64: Multiplier used internally to convert objective sense from Max to Min. Only takes on values of {-1.0, 1.0}

  • _obj_var_slack_added::Bool: Flag to indicate if a slack variable was added for the objective function. This is done in some epigraph reformulations (see reform_epigraph_min!)

  • _stack::DataStructures.BinaryMinMaxHeap{NodeBB}: A heap of all nodes in the branch-and-bound tree

  • _current_node::NodeBB: The individual node being examined at any particular time. Nodes are removed from the stack and placed here, evaluated, and then sent back to the stack

  • _first_relax_point_set::Bool: (Unused) Flag for relaxation points

  • _current_xref::Vector{Float64}: (Unused) Variable values of a particular point

  • _candidate_xref::Vector{Float64}: (Unused) Variable values of a candidate point

  • _use_prior_objective_xref::Bool: (Unused) Flag to use variable values from previous evaluation on the current step

  • _current_objective_xref::Vector{Float64}: (Unused) Variable values for objective evaluation

  • _prior_objective_xref::Vector{Float64}: (Unused) Variable values for previous objective evaluation

  • _user_branch_variables::Bool: Flag for if the user has specified branch variables (see label_branch_variables!)

  • _fixed_variable::Vector{Bool}: Variables that are fixed in place

  • _branch_variable_count::Int64: Number of variables that can be branched on

  • _branch_to_sol_map::Vector{Int64}: Mapping from the branch variables to the full set of variables in the problem

  • _sol_to_branch_map::Vector{Int64}: Mapping from the full set of variables in the problem to the branch variables

  • _continuous_solution::Vector{Float64}: The final (or intermediate) variable values of the solution

  • _preprocess_feasibility::Bool: Flag to ensure preprocessing result is feasible

  • _preprocess_termination_status::MathOptInterface.TerminationStatusCode: Status codes for use in bounds tightening

  • _preprocess_primal_status::MathOptInterface.ResultStatusCode: Status codes for use in bounds tightening

  • _preprocess_dual_status::MathOptInterface.ResultStatusCode: Status codes for use in bounds tightening

  • _lower_primal_status::MathOptInterface.ResultStatusCode: Primal status of the lower problem

  • _lower_dual_status::MathOptInterface.ResultStatusCode: Dual status of the lower problem

  • _lower_termination_status::MathOptInterface.TerminationStatusCode: Termination status of the lower problem

  • _lower_feasibility::Bool: Flag for lower problem feasibility

  • _lower_objective_value::Float64: Objective value result from the lower problem

  • _lower_solution::Vector{Float64}: Variable values of the lower problem solution

  • _lower_lvd::Vector{Float64}: Lower variable duals for use in duality-based bound tightening

  • _lower_uvd::Vector{Float64}: Upper variable duals for use in duality-based bound tightening

  • _last_cut_objective::Float64: Objective value associated with the previous cut in the cutting planes algorithm

  • _upper_result_status::MathOptInterface.ResultStatusCode: Primal status of the upper problem

  • _upper_termination_status::MathOptInterface.TerminationStatusCode: Termination status of the upper problem

  • _upper_feasibility::Bool: Flag for upper problem feasibility

  • _upper_objective_value::Float64: Objective value result from the upper problem

  • _upper_variables::Vector{MathOptInterface.VariableIndex}:

  • _upper_solution::Vector{Float64}:

  • _postprocess_feasibility::Bool: (Unused) Flag to ensure postprocessing result is feasible

  • _time_left::Float64: Time remaining for the optimization algorithm. This is set in initial_parse! to the user-defined time limit and is decremented throughout global_solve!

  • _start_time::Float64: Storage for the time() when optimization began

  • _run_time::Float64: Current run time, incremented using time()-_start_time

  • _parse_time::Float64: A field to keep track of time spent on initial problem parsing

  • _presolve_time::Float64: Used in optimize_nonconvex.jl to track how long the presolve step takes

  • _last_preprocess_time::Float64: Updated each iteration to track the time of the preprocess step

  • _last_lower_problem_time::Float64: Updated each iteration to track the time of the lower problem step

  • _last_upper_problem_time::Float64: Updated each iteration to track the time of the upper problem step

  • _last_postprocessing_time::Float64: Updated each iteration to track the time of the postprocess step

  • _min_converged_value::Float64: A field to track convergence progress across iterations

  • _global_lower_bound::Float64: The best-known lower bound

  • _global_upper_bound::Float64: The best-known upper bound

  • _maximum_node_id::Int64: The total number of nodes that have been created

  • _iteration_count::Int64: The number of iterations the branch-and-bound algorithm has completed

  • _node_count::Int64: The number of nodes in the stack

  • _solution_value::Float64: (Unused) The best-known solution value

  • _feasible_solution_found::Bool: A flag for if a feasible solution was identified. Updated if preprocessing, lower problem, and upper problem all return feasible values

  • _solution_node::Int64: The node ID of the best-known feasible upper problem solution (default = -1, if no feasible solution is found)

  • _best_upper_value::Float64: The best-known upper bound

  • _obbt_working_lower_index::Vector{Bool}: Indices of variables to perform OBBT on

  • _obbt_working_upper_index::Vector{Bool}: Indices of variables to perform OBBT on

  • _lower_indx_diff::Vector{Bool}: Tracker for changes in obbtworkinglowerindex across iterations

  • _upper_indx_diff::Vector{Bool}: Tracker for changes in obbtworkingupperindex across iterations

  • _old_low_index::Vector{Bool}: Storage for indices prior to OBBT step

  • _old_upp_index::Vector{Bool}: Storage for indices prior to OBBT step

  • _new_low_index::Vector{Bool}: New indices following OBBT step; compared with _old_low_index

  • _new_upp_index::Vector{Bool}: New indices following OBBT step; compared with _old_upp_index

  • _obbt_variables::Vector{MathOptInterface.VariableIndex}: (Deprecated) Variables to perform OBBT on. Replaced by _obbt_working_lower_index and _obbt_working_upper_index

  • _obbt_variable_count::Int64: The number of variables to perform OBBT on

  • _obbt_performed_flag::Bool: (Unused) Flag to indicate whether OBBT has been performed

  • _lower_fbbt_buffer::Vector{Float64}: Buffer for FBBT lower bounds. Set in presolve, used in preprocess

  • _upper_fbbt_buffer::Vector{Float64}: Buffer for FBBT upper bounds. Set in presolve, used in preprocess

  • _cp_improvement::Float64: (Unused) Improvement in constraint propagation

  • _cp_evaluation_reverse::Bool: (Unused) Flag for if constraint propagation results need to be reversed

  • _cut_iterations::Int64: Iterations of the cutting planes algorithm completed

  • _cut_add_flag::Bool: (Unused) Flag to check if cuts should be added

  • _node_repetitions::Int64: Counter for number of times a node is evaluated. If the repeat_check function is overloaded to return true, a node will not be branched on, but will instead be added back into the stack using single_storage!. In this case, _node_repetitions is incremented

  • _log::Log: Storage for logging information during a branch-and-bound run

  • _affine_relax_ci::Vector{MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.LessThan{Float64}}}: Storage for affine constraints

  • _affine_objective_cut_ci::Union{Nothing, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarAffineFunction{Float64}, MathOptInterface.LessThan{Float64}}, MathOptInterface.ConstraintIndex{MathOptInterface.VariableIndex, MathOptInterface.LessThan{Float64}}}: Storage for a linear objective cut constraint

  • _relaxed_variable_number::Int64: (Unused) Number of relaxed variables

  • _relaxed_variable_index::Vector{MathOptInterface.VariableIndex}: Indices of relaxed variables

  • _relaxed_variable_et::Vector{Tuple{MathOptInterface.ConstraintIndex{MathOptInterface.VariableIndex, MathOptInterface.EqualTo{Float64}}, Int64}}: Stored EqualTo constraints

  • _relaxed_variable_lt::Vector{Tuple{MathOptInterface.ConstraintIndex{MathOptInterface.VariableIndex, MathOptInterface.LessThan{Float64}}, Int64}}: Stored LessThan constraints

  • _relaxed_variable_gt::Vector{Tuple{MathOptInterface.ConstraintIndex{MathOptInterface.VariableIndex, MathOptInterface.GreaterThan{Float64}}, Int64}}: Stored GreaterThan constraints

  • _relaxed_variable_integer::Vector{MathOptInterface.ConstraintIndex{MathOptInterface.VariableIndex, MathOptInterface.Integer}}: Stored Integer constraints

  • _branch_variables::Vector{Bool}: List of variables that can be branched on. If not user-specified, branch variables are identified in label_branch_variables!

  • _nonbranching_int::Bool: (Unused) Flag for non-branching integers

  • _new_eval_constraint::Bool: Flag indicating if an initial evaluation of the constraints has occurred

  • _new_eval_objective::Bool: Flag indicating if the objective expression was evaluated

  • _node_to_sv_leq_ci::Dict{Int64, MathOptInterface.ConstraintIndex{MathOptInterface.VariableIndex, MathOptInterface.LessThan{Float64}}}: Storage for carrying LessThan constraint information. Used in obbt! and update_relaxed_problem_box!

  • _node_to_sv_geq_ci::Dict{Int64, MathOptInterface.ConstraintIndex{MathOptInterface.VariableIndex, MathOptInterface.GreaterThan{Float64}}}: Storage for carrying GreaterThan constraint information. Used in obbt! and update_relaxed_problem_box!

  • _nonlinear_evaluator_created::Bool: Flag to check for nonlinear evaluators. Set to true in add_nonlinear_evaluator!

  • _branch_cost::EAGO.BranchCostStorage{Float64}: (FUTURE FEATURE, NOT CURRENTLY IMPLEMENTED) Storage for pseudocost branching

  • _branch_variable_sparsity::SparseArrays.SparseMatrixCSC{Bool, Int64}: (FUTURE FEATURE, NOT CURRENTLY IMPLEMENTED) Sparsity information of the branch variables

  • _constraint_infeasiblity::Vector{Float64}: (FUTURE FEATURE, NOT CURRENTLY IMPLEMENTED) Information on the infeasibility of each constraint

source

Customizable Subroutines

Stack Management Subroutines

EAGO.branch_node!Method
branch_node!(t::ExtensionType, m::GlobalOptimizer)

Create two nodes from current_node and store them on the stack. Call select_branch_variable(t, m) and select_branch_point(t, m, k) to determine the variable that should be branched on and the point at which branching should occur, respectively.

source
EAGO.select_branch_variableMethod
select_branch_variable(
    t::ExtensionType,
    m::GlobalOptimizer
) -> Any

Choose a variable to branch on. A maximum relative width branching rule is used by default.

source
EAGO.select_branch_pointMethod
select_branch_point(
    t::ExtensionType,
    m::GlobalOptimizer,
    i
) -> Float64

Select a point xb within the domain of the ith branching variable. By default, this point is a convex combination of the solution to the relaxation and the midpoint of the node (branch_cvx_factor*xmid + (1-branch_cvx_factor)*xsol). If the solution lies within branch_offset of a bound, then the branch point is moved to a distance of branch_offset from that bound.

source
EAGO.node_selection!Method
node_selection!(t::ExtensionType, m::GlobalOptimizer)

Select the next node in the stack to evaluate. By default, perform best-first node selection (select the node with the lowest lower bound in the stack).

source
EAGO.fathom!Method
fathom!(t::ExtensionType, m::GlobalOptimizer)

Remove nodes from the stack. By default, delete nodes from the stack if their lower bounds are greater than the current global upper bound.

source
EAGO.initialize_stack!Method
initialize_stack!(t::ExtensionType, m::GlobalOptimizer)

Prepare the stack for the branch-and-bound routine. By default, create an initial node with the variable bounds as box constraints and add it to the stack.

source
EAGO.single_storage!Method
single_storage!(t::ExtensionType, m::GlobalOptimizer)

Store the current node to the stack, without branching, after updating lower/upper bounds.

source

Internal Subproblem Status Codes and Subsolver Management

EAGO.RelaxResultStatusType
RelaxResultStatus

Status code used internally to determine how to interpret the results from the solution of a relaxed problem.

source
EAGO.LocalResultStatusType
LocalResultStatus

Status code used internally to determine how to interpret the results from the solution of a local problem solve.

source
EAGO.IncrementalType
mutable struct Incremental{S<:MathOptInterface.AbstractOptimizer} <: MathOptInterface.AbstractOptimizer

A type-stable cache used to wrapper for an optimizer that enables incremental modification of solvers that don't inherently suppport this. Explicitly checks support of MOI functionality used in EAGO.

(Deprecated) For Q = Val{true}, the subsolver supports incremental loading. For Q = Val{false}, the subsolver does not.

source
EAGO.SubSolversType
mutable struct SubSolvers{Q<:MathOptInterface.AbstractOptimizer, S<:MathOptInterface.AbstractOptimizer, T<:ExtensionType}

A structure containing the relaxed and upper optimizers to be used, as well as any user-defined extension.

  • relaxed_optimizer::MathOptInterface.AbstractOptimizer: Optimizer used to solve relaxed subproblems. Set using r = [...] (<: MOI.AbstractOptimizer) (default = Cbc.Optimizer())

  • upper_optimizer::MathOptInterface.AbstractOptimizer: Optimizer used to solve upper bounding problems. Set using u = [...] (<: MOI.AbstractOptimizer) (default = Ipopt.Optimizer())

  • ext::ExtensionType: User-defined extension to use. Set using t = [...](<: EAGO.ExtensionType)

source
EAGO.set_default_config!Method
set_default_config!

Configures subsolver tolerances based on tolerance parameters provided to EAGO (provided that a specialized subsolver configuration routine has been provided and m.user_solver_config = false).

source

Main Subproblem and Termination Subroutines

EAGO.convergence_checkMethod
convergence_check(
    t::ExtensionType,
    m::GlobalOptimizer
) -> Bool

Check for problem convergence.

By default, check if the lower and upper bounds have converged to within absolute and/or relative tolerances.

source
EAGO.cut_conditionMethod
cut_condition(t::ExtensionType, m::GlobalOptimizer) -> Bool

Returns true if a cut should be added and computes a new reference point to add the cut at. By default, checks that cut_max_iterations are not exceeded and that the improvement in the objective value associated with the previous cut is greater than both an absolute tolerance cut_ϵ_abs and a relative tolerance cut_ϵ_rel. Returns false otherwise.

source
EAGO.lower_problem!Method
lower_problem!(
    t::ExtensionType,
    m::GlobalOptimizer{R, S, Q<:ExtensionType}
)

Constructs a relaxation of the MINLP on node y and solves it using the default EAGO relaxation scheme. By default, EAGO applies Kelley's algorithm (from Kelley Jr., J.E.: The cutting-plane method for solving convex programs. J. Soc. Ind. Appl. Math. 8(4), 703 to 712 (1960)) while cut_condition(m) returns true then activates the integrality constraints of the relaxed problems and solves the resulting MILP relaxation. results are stored to the _lower_solution, _lower_termination_status, _lower_primal_status, _lower_dual_status, _lower_objective_value, and _lower_feasibility. Further, lower and upper variable duals are stored _lower_lvd and _lower_uvd, respectively, for use in duality based bound tightening. If relaxation-based bounds are weaker or cutting-planes are numerically poorly ill-posed, then interval bounds are used instead. If the problem is dual feasible but the primal status is ambiguous the dual objective value is used for the lower bound to avoid numerical issues.

source
EAGO.preprocess!Method
preprocess!(
    t::ExtensionType,
    m::GlobalOptimizer{R, S, Q<:ExtensionType}
)

Runs contractor methods prior to solving lower bounding problem. By default linear and quadratic contractor methods followed by interval constraint propagation then optimization-based bound tightening for a specified number of iterations while the subproblem at current node n has not been proven infeasible.

source
EAGO.postprocess!Method
postprocess!(t::ExtensionType, m::GlobalOptimizer)

Default postprocess perfoms duality-based bound tightening (Tawarmalani, M., Sahinidis, N.V.: Global optimization of mixed-integer nonlinear programs: a theoretical and computational study. Math. Progr. 99, 563–591 (2004).) up to an iteration limit set by dbbt_depth.

source
EAGO.repeat_checkMethod
repeat_check(t::ExtensionType, m::GlobalOptimizer) -> Bool

Check to see if current node should be reprocessed. Without any custom extension, return false by default.

source
EAGO.termination_checkMethod
termination_check(m::GlobalOptimizer)
termination_check(t::ExtensionType, m::GlobalOptimizer) -> Bool

Check for termination of the branch-and-bound algorithm.

If only the GlobalOptimizer is given as an argument, termination_check dispatches to the other form using the ExtensionType given in the SubSolvers. If there is no user-defined extension, then by default, this will check for satisfaction of absolute or relative tolerances, solution infeasibility, and other specified limits. Returns true if any conditions are met and branch-and-bound should end, and false otherwise.

source
EAGO.upper_problem!Method
upper_problem!(t::ExtensionType, m::GlobalOptimizer)

Default upper bounding problem which simply calls solve_local_nlp! to solve the NLP locally.

source
EAGO.parse_global!Method

Basic parsing for global solutions (no extensive manipulation). By default, does nothing.

source
EAGO.optimize_hook!Method
optimize_hook!(t::ExtensionType, m::Optimizer)

Provide a hook for extensions to EAGO.

The user-defined extension of optimize_hook! is used in EAGO's overloading of MOI.optimize! (see EAGO.jl/src/eago_optimizer/optimize/optimize.jl). Without the optimize_hook! specified, EAGO will run initial_parse!, parse_classify_problem!, and then optimize! using the parsed problem type. The user-specified optimize_hook! should thus take the new extension and Optimizer as inputs and will execute when the user writes optimize!(model).

Example

Here, optimize_hook! is used to bypass EAGO's problem parsing and treat every problem using its branch-and-bound routine. This is done in this example by telling EAGO to treat the problem as a mixed integer nonconvex problem, which normally dispatches to branch-and-bound.

struct MyNewExtension <: EAGO.ExtensionType end
import EAGO: optimize_hook!
function EAGO.optimize_hook!(t::MyNewExtension, m::Optimizer)
    initial_parse!(m)
    optimize!(EAGO.MINCVX(), m)
end

The same functionality could be accomplished by setting the EAGOParameter field force_global_solve to be true.

source

Internal Subroutines

EAGO.is_integer_subproblemMethod
is_integer_subproblem(m)

Returns true that the subproblem at the current node n has participating integer variables that have not been fixed to constant valued as the branch-and-bound algorithm progresses. Returns false otherwise.

source
EAGO.is_integer_feasible_localMethod
is_integer_feasible_local(m::GlobalOptimizer, d) -> Bool

Checks that the solution of a local solve is integer feasible to within the tolerances specified by integer_abs_tol and integer_rel_tol.

source
EAGO.is_integer_feasible_relaxedMethod
is_integer_feasible_relaxed(m::GlobalOptimizer) -> Bool

Check that the solution of the lower (relaxed problem) is integer feasible to within tolerances specified by the parameters: integer_abs_tol (absolute tolerance) and integer_rel_tol (relative tolerance).

source
EAGO.interval_boundFunction
interval_bound(::GlobalOptimizer, ::T)

Compute a tuple representing the lower and upper interval bounds for an AbstractEAGOConstraint representing an equality constraint.

Options for T (all are subtypes of AbstractEAGOConstraint):

  • AffineFunctionEq
  • AffineFunctionIneq
  • BufferedQuadraticEq
  • BufferedQuadraticIneq
  • BufferedNonlinearFunction{V,N,T} where {V,N,T}
source
EAGO.lower_interval_boundFunction
lower_interval_bound(::GlobalOptimizer, ::T)

Compute the lower interval bound for an AbstractEAGOConstraint representing an inequality constraint.

Options for T (all are subtypes of AbstractEAGOConstraint):

  • AffineFunctionIneq
  • BufferedQuadraticIneq
  • BufferedSOC
  • BufferedNonlinearFunction{V,N,T} where {V,N,T}
source
EAGO.same_boxMethod
same_box(x::NodeBB, y::NodeBB, r::Float64) -> Bool

Check that node x and y have equal domains within an absolute tolerance of r.

source
EAGO.solve_local_nlp!Method

Constructs and solves the problem locally on node y updated the upper solution informaton in the optimizer.

source
EAGO.set_dual!Method
set_dual!(m)

Retrieves the lower and upper duals for variable bounds from the relaxed_optimizer and sets the appropriate values in the _lower_lvd and _lower_uvd storage fields.

source
EAGO.update_relaxed_problem_box!Function
update_relaxed_problem_box!(m)

Update the relaxed constraint by setting the constraint set of v == x* , xL_i <= x_i, and x_i <= xU_i for each such constraint added to the relaxed optimizer. Resets integral valued constraints to either EqualTo or Interval constraints.

source
EAGO.reform_epigraph_min!Method
reform_epigraph_min!(m::GlobalOptimizer)

Perform an epigraph reformulation assuming the working_problem is a minimization problem.

source
EAGO.label_fixed_variables!Method
label_fixed_variables!(m::GlobalOptimizer) -> Vector{Bool}

Detect any variables set to a fixed value by equality or inequality constraints and populate the _fixed_variable storage array.

source
EAGO.label_branch_variables!Method
label_branch_variables!(m::GlobalOptimizer)

Detect any variables participating in nonconvex terms and populate the _branch_variables storage array.

source
EAGO.add_nonlinear!Method
add_nonlinear!(m::GlobalOptimizer)

Add an Evaluator and nonlinear functions and populate each appropriately.

source
EAGO.parse_classify_problem!Method
parse_classify_problem!(m::GlobalOptimizer)

Interprets the type/number of constraints and the type of objective function to infer a problem type. Current possible types include:

  • LP: Linear program; sent to optimize_lp.jl
  • 'MILP: Mixed integer linear program; sent tooptimize_lp.jl`
  • SOCP: Second-order cone program; sent to optimize_conic.jl
  • MINCVX: Mixed-integer nonconvex; sent to optimize_nonconvex.jl

If the force_global_solve parameter is set to true, parse_classify_problem! will set the problem type to MINCVX to pass the problem to optimize_nonconvex.jl.

source
EAGO.local_problem_statusMethod
local_problem_status(t, r)

Takes an MOI.TerminationStatusCode and a MOI.ResultStatusCode and returns true if this corresponds to a solution that is proven to be feasible. Returns false otherwise.

source

Functions for Generating Console Output

EAGO.print_iteration!Function

print_iteration!

Print status information based on iteration count. The header print frequency is based on the header_iterations setting, and the data print frequency is based on the output_iterations setting.

source
EAGO.print_node!Function

print_node!

Print information about the current node. Includes node ID, lower bound, upper bound, and interval box.

source
EAGO.print_results!Function
print_results!(m::GlobalOptimizer, lower_flag::Bool)

Print the results of a single (lower or upper) bounding problem. lower_flag=true prints information for the lower problem, lower_flag=false prints information for the upper problem.

source
EAGO.print_solution!Function

print_solution!

Print solution information for the B&B problem. Display node with the best solution, solution value, solution, and time spent solving subproblems. This print occurs following termination of the B&B algorithm.

source

Support for Log Output at Each Iteration

EAGO.LogType
mutable struct Log

A structure used to store information on the history of the solution procedure for generating convergence plots and other analyses.

  • current_lower_bound::Vector{Float64}: Storage for lower bound calculated for current node.

  • current_upper_bound::Vector{Float64}: Storage for upper bound calculated for current node.

  • preprocessing_time::Vector{Float64}: Storage for preprocessing time of each iteration.

  • lower_problem_time::Vector{Float64}: Storage for lower bounding time of each iteration.

  • upper_problem_time::Vector{Float64}: Storage for upper bounding time of each iteration.

  • postprocessing_time::Vector{Float64}: Storage for postprocessing time of each iteration.

  • preprocessing_feas::Vector{Bool}: Storage for preprocessing feasibility of each iteration.

  • lower_problem_feas::Vector{Bool}: Storage for lower bounding feasibility of each iteration.

  • upper_problem_feas::Vector{Bool}: Storage for upper bounding feasibility of each iteration.

  • postprocessing_feas::Vector{Bool}: Storage for postprocessing feasibility of each iteration.

  • global_lower_bound::Vector{Float64}: Storage for best (global) lower bound at each iteration.

  • global_upper_bound::Vector{Float64}: Storage for best (global) upper bound at each iteration.

  • node_count::Vector{Int64}: Number of nodes at each iteration.

  • run_time::Vector{Float64}: Run time at each iteration.

source
EAGO.log_iteration!Method
log_iteration!(m::GlobalOptimizer)

If log_on is true, the global_lower_bound, global_upper_bound, run_time, and node_count are stored every log_interval. If log_subproblem_info then the lower bound, feasibility and run times of the subproblems are logged every log_interval.

source

Interval Representations of Expressions

EAGO.AffineFunctionEqType
mutable struct AffineFunctionEq <: EAGO.AbstractEAGOConstraint

Representation of an affine equality. Currently only used for bound tightening.

  • terms::Vector{Tuple{Float64, Int64}}

  • constant::Float64

  • len::Int64

source
EAGO.AffineFunctionIneqType
mutable struct AffineFunctionIneq <: EAGO.AbstractEAGOConstraint

Representation of an affine inequality. Currently only used for bound tightening.

  • terms::Vector{Tuple{Float64, Int64}}

  • constant::Float64

  • len::Int64

source
EAGO.BufferedQuadraticIneqType
mutable struct BufferedQuadraticIneq <: EAGO.AbstractEAGOConstraint

Representation of a general quadratic inequality constraint with a buffer.

  • func::MathOptInterface.ScalarQuadraticFunction{Float64}

  • buffer::Dict{Int64, Float64}

  • saf::MathOptInterface.ScalarAffineFunction{Float64}

  • len::Int64

source
EAGO.BufferedQuadraticEqType
mutable struct BufferedQuadraticEq <: EAGO.AbstractEAGOConstraint

Representation of a general quadratic equality constraint with a buffer.

  • func::MathOptInterface.ScalarQuadraticFunction{Float64}

  • minus_func::MathOptInterface.ScalarQuadraticFunction{Float64}

  • buffer::Dict{Int64, Float64}

  • saf::MathOptInterface.ScalarAffineFunction{Float64}

  • len::Int64

source
EAGO.NonlinearExpressionType
mutable struct NonlinearExpression{V, N, T<:RelaxTag} <: EAGO.AbstractEAGOConstraint

Stores a general quadratic function with a buffer.

source
EAGO.BufferedNonlinearFunctionType
mutable struct BufferedNonlinearFunction{V, N, T<:RelaxTag} <: EAGO.AbstractEAGOConstraint

Stores a general nonlinear function with a buffer represented by the sum of a tape and a scalar affine function.

source