Linearly Constrained Optimization
Mathwrist C++ API Series

Copyright ©Mathwrist LLC 2023
(January 20, 2023)
Abstract

This document details the C++ programming interface of linearly constrained (LC) optimization features in Mathwrist’s Numerical Programming Library (NPL). Please refer to our LC white paper for technical overview.

1 Introduction

Let ψ⁢(𝐱):ℝn→ℝ be a general smooth function with gradient 𝐠⁢(𝐱) and Hessian 𝐇⁢(𝐱). In other words, ψ⁢(𝐱) is at least twice continuously differentiable. Linearly constrained optimization solves the following minimization problem subject to a set of general linear constraints and simple bound constraints,

minx∈ℝn⁡ψ⁢(𝐱), s.t. ⁢𝐛l≤𝐀𝐱≤𝐛u⁢ and ⁢𝐥≤𝐱≤𝐮 (1)

In header file ”gs_active_set.h”, we provide two general LC solver classes

  • •

    GS_ActiveSetNewton

  • •

    GS_ActiveSetQuasiNewton

that embde line search methods within the active-set framework to solve problem (1).

A special case of linearly constrained optimization is when there are only simple bound constraints imposed on 𝐱, also known as box constrained optimization,

minx∈ℝn⁡ψ⁢(𝐱), s.t. ⁢𝐥≤𝐱≤𝐮 (2)

For efficiency consideration, this special case is separately treated in header file ”gs_simple_bound.h” by the following two solver classes.

  • •

    GS_SimpleBoundNewton

  • •

    GS_SimpleBoundQuasiNewton

2 Active-set Solvers

2.1 Inheritance Hierarchy

Class GS_ActiveSet in header file gs_active_set.h is the common base class for general LC solvers. It combines line search method with active-set framework. The concrete Newton and Quasi-Newton solvers derived from GS_ActiveSet compute descent directions using modified Newton or Quasi-Newton method respectively.

1 class ActiveSetWithLagrange : public ActiveSet
2 {
3 };
4 class GS_ActiveSet : public ActiveSetWithLagrange
5 {
6 };
7 class GS_ActiveSetNewton : public GS_ActiveSet
8 {
9 };
10 class GS_ActiveSetQuasiNewton : public GS_ActiveSet
11 {
12 };

Please refer to our linear programming API documentation for the details of ActiveSet base class.

2.2 Step Length Control

1 void set_max_step(double step_length);
2 double max_step() const;

Description:

  1. 1.

    Set the max step length in line search method. After computing a descent direction 𝐩, we use the line search method to determine an acceptable step length α^ along the direction 𝐩. Users can set a max step length α¯ upper bound based on the domain knowledge about ψ⁢(𝐱) so that α^≤α¯.

  2. 2.

    Return the max step length used in line search method.

2.3 Solver Operator

1 double operator() (Matrix& x);

Description:

Given an initial guess, compute the solution 𝐱* of the general LC problem (1) and return ψ⁢(𝐱*). Parameter x holds initial guess on call and solution 𝐱* on return.

2.4 Newton Method

1 GS_ActiveSetNewton(
2   const FunctionND& obj,
3   const LinearProg::Constraints& constr,
4   bool discrete);

Description:

This is the constructor to create a solver object using active-set and line search modified Newton methods to solve LC problem (1).

Parameter obj is a constant reference to the n-d objective function ψ⁢(𝐱). The solver class stores the reference as a data member. REQUIRED: The n-d function object referred by obj needs be alive in the lifecycle of the LC solver.

Parameter constr is a constant reference to the constraint specification object in LC problem formulation (1). The solver’s ActiveSet base class stores the reference as a data member. Please refer to our LP API documentation for the details of constraint specification. REQUIRED: The constraint object constr needs be alive in the lifecycle of the LC solver.

Parameter discrete controls whether the Hessian matrix 𝐇⁢(𝐱) is approximated by a finite difference method or explicitly computed by the objective function ψ⁢(𝐱). When discrete is true, FunctionND::hess_fwd_diff_grad() is used to approximate 𝐇⁢(𝐱). Otherwise, users need explicitly compute 𝐇⁢(𝐱) in the implementation of FunctionND::eval() function.

2.5 Quasi-Newton Method

1 GS_ActiveSetQuasiNewton(
2   const FunctionND& obj,
3   const LinearProg::Constraints& constr);

Description:

This is the constructor to create a solver object using active-set and line search Quasi-Newton methods to solve LC problem (1). Parameters obj and constr are exactly same as the GS_ActiveSetNewton constructor.

3 Simple Bound Solvers

3.1 Inheritance Hierarchy

Class GS_SimpleBound in header file gs_simple_bound.h is the common base class for box constrained optimization solvers. GS_SimpleBound inherits from IterativeMethod. Please refer to our “1-d and n-d Functions” API documentation for iteration and convergence control of IterativeMethod.

1 class GS_SimpleBound : public IterativeMethod
2 {
3 };
4 class GS_SimpleBoundNewton : public GS_SimpleBound
5 {
6 };
7 class GS_SimpleBoundQuasiNewton : public GS_SimpleBound
8 {
9 };

The GS_SimpleBound base class implements a simplified version of active-set framework that just handles simple bound operations. It also uses the line search method to compute an acceptable step length α^ along a descent direction 𝐩. The concrete simple bound solvers GS_SimpleBoundNewton and GS_SimpleBoundQuasiNewton compute the direction 𝐩 using modified Newton and Quasi-Newton methods respectively.

3.2 Step Length Control

1 void set_max_step(double step_length);
2 double max_step() const;

Description:

  1. 1.

    Set the max step length in line search method. After computing a descent direction 𝐩, we use the line search method to determine an acceptable step length α^ along the direction 𝐩. Users can set a max step length α¯ upper bound based on the domain knowledge about ψ⁢(𝐱) so that α^≤α¯.

  2. 2.

    Return the max step length used in line search method.

3.3 Stationary Tolerance Control

1 void set_stationary_tolerance(double eps);
2 double stationary_tolerance() const;

Description:

  1. 1.

    Set the tolerance level ϵs to determine if stationary condition is satisfied. At an iteration point 𝐱k, if |𝐠~⁢(𝐱k)|∞<ϵs or 𝐠⁢(𝐱k)T⁢𝐩<ϵs, we consider 𝐱k is stationary wrt the current working set of constraints.

  2. 2.

    Return the stationary tolerance control ϵs.

3.4 Solver Operator

1 double operator() (Matrix& x);

Description:

Given an initial guess, compute the solution 𝐱* of the box constrained problem (2) and return ψ⁢(𝐱*). Parameter x holds initial guess on call and solution 𝐱* on return.

3.5 Newton Method

1 GS_SimpleBoundNewton(
2   const FunctionND& obj,
3   const ftl::Buffer<Bound>& simple_bounds,
4   bool discrete);

Description:

This is the constructor to create a solver object that uses active-set and line search modified Newton methods to solve box constrained problem (2).

Parameter obj is a constant reference to the n-d objective function ψ⁢(𝐱). The solver class stores the reference as a data member. REQUIRED: The n-d function object referred by obj needs be alive in the lifecycle of the solver.

Parameter simple_bounds is a buffer of simple bounds in problem formulation (2). All simple bounds are copied to the solver’s GS_SimpleBound base class as data member. Please refer to our LP API documentation for the details of Bound specification. REQUIRED: The size of the simple bound buffer should be same as the dimension size of objective function ψ⁢(𝐱).

Parameter discrete controls whether the Hessian matrix 𝐇⁢(𝐱) is approximated by a finite difference method or explicitly computed by the objective function ψ⁢(𝐱). When discrete is true, FunctionND::hess_fwd_diff_grad() is used to approximate 𝐇⁢(𝐱). Otherwise, users need explicitly compute 𝐇⁢(𝐱) in the implementation of FunctionND::eval() function.

3.6 Quasi-Newton Method

1 GS_SimpleBoundQuasiNewton(
2   const FunctionND& obj,
3   const ftl::Buffer<Bound>& simple_bounds);

Description:

This is the constructor to create a solver object that uses active-set and line search Quasi-Newton methods to solve box constrained problem (2). Parameters obj and constr are exactly same as the GS_SimpleBoundQuasiNewton constructor.