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,

minxnψ(𝐱), 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,

minxnψ(𝐱), 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.