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.NodeBB
— Typestruct 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 valuedcontinuous::Bool
: Are all dimensions continuous (or fixed)lower_bound::Float64
: Lower bound of problem solution on nodeBBupper_bound::Float64
: Upper bound of problem solution on nodeBBdepth::Int64
: Depth of node in B&B tree.cont_depth::Int64
: Depth of first parent in B&B tree that was continuously valuedid::Int64
: Unique ID for each node.branch_direction::EAGO.BranchDirection
: Whether last branch was negative or positive in directionlast_branch::Int64
: Dimension of last branchbranch_extent::Float64
: Extent of last branch (using for psuedocost calculation)
The GlobalOptimizer
structure holds all information relevant to branch-and-bound.
EAGO.GlobalOptimizer
— Typemutable 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 variablesobbt_variable_values::Vector{Bool}
: Variables to perform OBBT on (default: all variables in nonlinear expressions)enable_optimize_hook::Bool
: Specifies that theoptimize_hook!
function should be called rather than throw the problem to the standard routineext::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 (seereform_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 (seelabel_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 ininitial_parse!
to the user-defined time limit and is decremented throughoutglobal_solve!
_start_time::Float64
: Storage for thetime()
when optimization began_run_time::Float64
: Current run time, incremented usingtime()-_start_time
_parse_time::Float64
: A field to keep track of time spent on initial problem parsing_presolve_time::Float64
: Used inoptimize_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_last_printed_iteration::Int64
: Updated each time an iteration is printed_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 therepeat_check
function is overloaded to returntrue
, a node will not be branched on, but will instead be added back into the stack usingsingle_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 inlabel_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 inobbt!
andupdate_relaxed_problem_box!
_node_to_sv_geq_ci::Dict{Int64, MathOptInterface.ConstraintIndex{MathOptInterface.VariableIndex, MathOptInterface.GreaterThan{Float64}}}
: Storage for carrying GreaterThan constraint information. Used inobbt!
andupdate_relaxed_problem_box!
_nonlinear_evaluator_created::Bool
: Flag to check for nonlinear evaluators. Set totrue
inadd_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
Customizable Subroutines
Stack Management Subroutines
EAGO.branch_node!
— Methodbranch_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.
EAGO.select_branch_variable
— Methodselect_branch_variable(
t::ExtensionType,
m::GlobalOptimizer
) -> Any
Choose a variable to branch on. A maximum relative width branching rule is used by default.
EAGO.select_branch_point
— Methodselect_branch_point(
t::ExtensionType,
m::GlobalOptimizer,
i
) -> Float64
Select a point xb
within the domain of the i
th 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.
EAGO.node_selection!
— Methodnode_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).
EAGO.fathom!
— Methodfathom!(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.
EAGO.initialize_stack!
— Methodinitialize_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.
EAGO.single_storage!
— Methodsingle_storage!(t::ExtensionType, m::GlobalOptimizer)
Store the current node to the stack, without branching, after updating lower/upper bounds.
Internal Subproblem Status Codes and Subsolver Management
EAGO.RelaxResultStatus
— TypeRelaxResultStatus
Status code used internally to determine how to interpret the results from the solution of a relaxed problem.
EAGO.LocalResultStatus
— TypeLocalResultStatus
Status code used internally to determine how to interpret the results from the solution of a local problem solve.
EAGO.Incremental
— Typemutable 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.
EAGO.SubSolvers
— Typemutable 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 usingr = [...]
(<:MOI.AbstractOptimizer
) (default =Cbc.Optimizer()
)upper_optimizer::MathOptInterface.AbstractOptimizer
: Optimizer used to solve upper bounding problems. Set usingu = [...]
(<:MOI.AbstractOptimizer
) (default =Ipopt.Optimizer()
)ext::ExtensionType
: User-defined extension to use. Set usingt = [...]
(<:EAGO.ExtensionType
)
EAGO.set_default_config!
— Methodset_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
).
Main Subproblem and Termination Subroutines
EAGO.convergence_check
— Methodconvergence_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.
EAGO.cut_condition
— Methodcut_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.
EAGO.lower_problem!
— Methodlower_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.
EAGO.preprocess!
— Methodpreprocess!(
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.
EAGO.postprocess!
— Methodpostprocess!(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
.
EAGO.repeat_check
— Methodrepeat_check(t::ExtensionType, m::GlobalOptimizer) -> Bool
Check to see if current node should be reprocessed. Without any custom extension, return false
by default.
EAGO.termination_check
— Methodtermination_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.
EAGO.upper_problem!
— Methodupper_problem!(t::ExtensionType, m::GlobalOptimizer)
Default upper bounding problem which simply calls solve_local_nlp!
to solve the NLP locally.
EAGO.parse_global!
— MethodBasic parsing for global solutions (no extensive manipulation). By default, does nothing.
EAGO.optimize_hook!
— Methodoptimize_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.
Internal Subroutines
EAGO.is_integer_subproblem
— Methodis_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.
EAGO.is_integer_feasible_local
— Methodis_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
.
EAGO.is_integer_feasible_relaxed
— Methodis_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).
EAGO.interval_bound
— Functioninterval_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}
EAGO.lower_interval_bound
— Functionlower_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}
EAGO.same_box
— Methodsame_box(x::NodeBB, y::NodeBB, r::Float64) -> Bool
Check that node x
and y
have equal domains within an absolute tolerance of r
.
EAGO.solve_local_nlp!
— MethodConstructs and solves the problem locally on node y
updated the upper solution informaton in the optimizer.
EAGO.set_dual!
— Methodset_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.
EAGO.update_relaxed_problem_box!
— Functionupdate_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.
EAGO.reform_epigraph_min!
— Methodreform_epigraph_min!(m::GlobalOptimizer)
Perform an epigraph reformulation assuming the working_problem is a minimization problem.
EAGO.label_fixed_variables!
— Methodlabel_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.
EAGO.label_branch_variables!
— Methodlabel_branch_variables!(m::GlobalOptimizer)
Detect any variables participating in nonconvex terms and populate the _branch_variables
storage array.
EAGO.add_nonlinear!
— Methodadd_nonlinear!(m::GlobalOptimizer)
Add an Evaluator and nonlinear functions and populate each appropriately.
EAGO.parse_classify_problem!
— Methodparse_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 tooptimize_lp.jl
- 'MILP
: Mixed integer linear program; sent to
optimize_lp.jl` SOCP
: Second-order cone program; sent tooptimize_conic.jl
MINCVX
: Mixed-integer nonconvex; sent tooptimize_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
.
EAGO.local_problem_status
— Methodlocal_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.
Functions for Generating Console Output
EAGO.print_iteration!
— Functionprint_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.
EAGO.print_node!
— Functionprint_node!
Print information about the current node. Includes node ID, lower bound, upper bound, and interval box.
EAGO.print_results!
— Functionprint_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.
EAGO.print_solution!
— Functionprint_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.
Support for Log Output at Each Iteration
EAGO.Log
— Typemutable 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.
EAGO.log_iteration!
— Methodlog_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
.
Interval Representations of Expressions
EAGO.AbstractEAGOConstraint
— Typeabstract type AbstractEAGOConstraint
An abstract super-type used for representing constraints built by EAGO's backend.
EAGO.AffineFunctionEq
— Typemutable struct AffineFunctionEq <: EAGO.AbstractEAGOConstraint
Representation of an affine equality. Currently only used for bound tightening.
terms::Vector{Tuple{Float64, Int64}}
constant::Float64
len::Int64
EAGO.AffineFunctionIneq
— Typemutable struct AffineFunctionIneq <: EAGO.AbstractEAGOConstraint
Representation of an affine inequality. Currently only used for bound tightening.
terms::Vector{Tuple{Float64, Int64}}
constant::Float64
len::Int64
EAGO.BufferedQuadraticIneq
— Typemutable 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
EAGO.BufferedQuadraticEq
— Typemutable 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
EAGO.NonlinearExpression
— Typemutable struct NonlinearExpression{V, N, T<:RelaxTag} <: EAGO.AbstractEAGOConstraint
Stores a general quadratic function with a buffer.
EAGO.BufferedNonlinearFunction
— Typemutable 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.