Interface to ConicBundle for the Language C++ using Matrix Classes

Solve $min_{y\in\mathbf{R}^m} f_0(y) + f_1(y) + ... + f_k(y)$ for convex functions f_i, the y-variables may be bounded or box constrained. The most important steps are explained here. More...

Classes

class  ConicBundle::PrimalMatrix
 If in Lagrangean relaxation primal solutions are in the form of a real vector or, more generally a matrix, then an approximate primal solution can be generated by supplying primal information of this form for each epsilon subgradient within ConicBundle::MatrixFunctionOracle::evaluate(). More...
class  ConicBundle::MatrixFunctionOracle
 Oracle interface (abstract class). For each of your functions, provide an instance of a derived class. More...
class  ConicBundle::MatrixCBSolver
 Bundle method solver. More...
class  ConicBundle::MatrixBSolver
 This is the common abstract Matrix Class interface to the two real solvers ConicBundle::MatrixFCBSolver and ConicBundle::MatrixNBSolver. More...
class  ConicBundle::MatrixFCBSolver
 The Full Conic Bundle method solver invoked by ConicBundle::MatrixCBSolver(), it uses a separate cutting model for each function. More...
class  ConicBundle::MatrixNBSolver
 The minimal bundle method solver invoked by ConicBundle::MatrixCBSolver(true), it uses no bundle but only one aggregate and one new subgradient. More...

Detailed Description

Solve $min_{y\in\mathbf{R}^m} f_0(y) + f_1(y) + ... + f_k(y)$ for convex functions f_i, the y-variables may be bounded or box constrained. The most important steps are explained here.

Before starting, note that this interface relies on several classes defined in CH_Matrix_Classes, in particular on CH_Matrix_Classes::Matrix and CH_Matrix_Classes::Indexmatrix. The full functionality of ConicBundle is only available with these classes. If you prefer an interface without having to use these, please use the Interface to ConicBundle for the Language C++. We now give a short overview of the most important steps in using the ConicBundle solver.

Setting up the Problem, the Functions, and the Main Loop

First create a new problem/solver ConicBundle::MatrixCBSolver, let us call it solver for brevity. In invoking the constructor a boolean flag may be set to enforce the use of a minimal bundle model that employs only one common aggregate and one common subgradient for all functions, so basically "no bundle", which may be favorable if fast iterations and/or little memory consumption are essential.

Next, set the dimension of the design variables/argument as well as possible box constraints on these by ConicBundle::MatrixCBSolver::init_problem().

Now set up each of your functions f_i as a ConicBundle::MatrixFunctionOracle. Via the routine ConicBundle::MatrixFunctionOracle::evaluate() you will supply, for a given argument, the function value and a subgradient (=the gradient if the function is differentiable) to the solver. The function evaluate() is the only function that you definitely have to provide, see the miniature example below.

The function oracles have to be added to the solver using the routine ConicBundle::MatrixCBSolver::add_function().

Once all functions are added, the optimization process can be started. If you know a good starting point then set it with ConicBundle::MatrixCBSolver::set_new_center_point() now, otherwise the method will pick the zero vector or, in the case of box constraints, the point closest to zero as starting point.

Finally, set up a loop that calls ConicBundle::MatrixCBSolver::do_descent_step() until ConicBundle::MatrixCBSolver::termination_code() is nonzero.

Retrieving Some Solution Information

After the first call to ConicBundle::MatrixCBSolver::do_descent_step() you can retrieve, at any time, the current objective value by ConicBundle::MatrixCBSolver::get_objval() and the argument leading to this value by ConicBundle::MatrixCBSolver::get_center(). For some screen output, use ConicBundle::MatrixCBSolver::set_out().

Lagrangean Relaxation, Primal Approximations, and Cutting Planes

If you are optimizing the Lagrange multipliers of a Lagrangean relaxation, you might be interested in getting an approximation to your primal optimal solution. This can be done by specifying in each function for each (epsilon) subgradient the corresponding primal vectors that generate it, see the parameter primal_solutions in ConicBundle::MatrixFunctionOracle::evaluate() as a start. Then for each of your functions, you can retrieve the current primal approximation using ConicBundle::MatrixCBSolver::get_approximate_primal().

If, in addition, you plan to improve your primal relaxation via cutting planes, that are strongly violated by the current primal approximation, you should have a look at ConicBundle::MatrixCBSolver::append_variables(), ConicBundle::MatrixFunctionOracle::subgradient_extension() and ConicBundle::MatrixCBSolver::reinit_function_model(). If you want to get rid of primal constraints/dual variables, use ConicBundle::MatrixCBSolver::get_approximate_slacks() and ConicBundle::MatrixCBSolver::delete_variables().

//******************************************************************************
//*       Miniature Example in C++ for Convex Quadratic in Two Variables       * 
//******************************************************************************
#include "MatCBSolver.hxx"

using namespace std;
using namespace ConicBundle;
using namespace CH_Matrix_Classes;

// f(x)=x^TAx+b^Tx+c with A=[5 1;1 4], b=[-12;-10], c=3 
class QFunction: public MatrixFunctionOracle
{
  Matrix A,b;
  double c;

public:
  QFunction()
  {
    A.init(2,2,0.);  A(0,0) = 5;  A(1,1) = 4;  A(0,1) = A(1,0) = 1.;
    b.init(2,1,0.);  b(0) = -12;  b(1) = -10;
    c=3;
  }

  int evaluate(const Matrix& x, double relprec, 
               double& objective_value,
               Matrix&  cut_vals,Matrix&  subgradients,
               vector<PrimalData*>&  primal_solutions,
               PrimalExtender*&)
  {
    /* compute objective as inner product of x with A*x+b plus c */
    objective_value= ip(x,A*x+b)+c;
    /* compute and store one subgradient(=gradient) with two coordinates */
    cut_vals.init(1,1,objective_value);
    subgradients=2.*A*x+b;
    return 0;
  }
};

int main(void)
{
  QFunction fun;

  MatrixCBSolver solver; 
  solver.init_problem(2);        // 2 variables, no bounds
  solver.add_function(fun);      
  solver.set_out(&cout,1);
 
  do {
    solver.do_descent_step();
  } while (!solver.termination_code());

  solver.print_termination_code(cout);

  cout<<" objval="<<solver.get_objval()<<"\n optimizer: ";
  Matrix x;
  solver.get_center(x);
  x.display(cout);

  return 0;
}

Generated on Mon Nov 8 19:36:39 2010 for ConicBundle by  doxygen 1.5.6