Matlab Interface Tutorial

The Matlab interface was entirely redesigned to correct a few long-standing bugs and to eliminate the limitations of the original Fortran 77 interface.

To start using the CUTEr tools and access problems from within Matlab, make sure all your environment variables are set up properly, as explained in the CUTEr documentation. Examples are found in your $CUTER/common/doc directory. In the tutorial below, commands typed at the shell prompt are prefixed with shell# and commands typed at the Matlab prompt are prefixed with matlab>>. The optimization problem is written mathematically in the form

minimize f(x)
subject to cL ≤ c(x) ≤ cU
bL ≤ x ≤ bU ,

where x, a vector of Rn, are the variables of the problem, f, a function from Rn into R, is the objective function, and c, a vector-valued function from Rn into Rm, represents the constraint functions, which can be either equalities, inequalities, or a mixture. The final constraints are explicit bounds on the variables. The vector x represents the set of design variables of the problem, the objective f models the cost to be minimized or maximized while the constraints c represent the conditions that acceptable values of the design variables must satisfy.

Problem Setup

Before problem data may be accessed from Matlab, the problem must be decoded and linked with the main Matlab gateway interface. To do this, run the following command at the shell prompt. We use problem LUBRIFC as an example.

  shell# runcuter --package mx --decode LUBRIFC

Note that, by contrast with earlier versions of CUTEr, you need not specify whether the problem is constrained or unconstrained any longer. The above command created a MEX library in the current directory. Depending on your architecture and platform, this library will have different names. On a 32-bit Pentium machine running Linux, the name is typically mcuter.mexglx. You are now ready to start up Matlab.

At the Matlab command line, problem data is loaded using the command

  matlab>> prob = cuter_setup();

This creates a Matlab structure stored in the variable prob, which holds problem-related data. Asking Matlab to print this structure, we get

  matlab>> prob

  prob = 

            n: 751
            m: 500
         nnzh: 1749
         nnzj: 64494
            x: [751x1 double]
           bl: [751x1 double]
           bu: [751x1 double]
            v: [500x1 double]
           cl: [500x1 double]
           cu: [500x1 double]
       equatn: [500x1 logical]
       linear: [500x1 logical]
         name: 'LUBRIFC'

The first four fields of prob are, in order, the number of variables of the problem, the number of constraints, the number of nonzero elements in a triangle of the Hessian of the Lagrangian and the number of nonzero elements in the Jacobian matrix of the constraints. Fields in Matlab are accessed using the dot operator:

  matlab>> prob.nnzh

  ans =

          1749

The next three fields of prob are double precision arrays containing the starting point, the vector of lower bounds and the vector of upper bounds. These vectors have size prob.n x 1. The arrays prob.bl and prob.bu are the direct counterparts of the vectors bL and bU in the mathematical programming problem above.

Just in the same way, the three fields prob.v, prob.cl and prob.cu are double precision arrays containing the initial Lagrange multipliers and the lower and upper constraint bounds. Note that the multipliers prob.v only concern the general constraints, and not the bound constraints.

Finally, the field prob.equatn is a logical array such that the i-th constraint is an equality constraint if and only if prob.equatn(i) = true. Similarly, the i-th constraint is linear if and only if prob.linear(i) = true. The last field, prob.name, should be self-explanatory.

Evaluating Problem Data

The CUTEr gateway interface offers a number of convenient functions that allow you to evaluate problem-related data, such as the objective function, the constraints, only selected constraints, the Hessian of the Lagrangian, etc. These commands are quite flexible as the following examples will show.

Evaluating the objective function at a given point is done by means of the cuter_obj function. This function receives any appropriately sized double precision array. For instance, prob.x is such an array. The command

  matlab>> f = cuter_obj( prob.x )
  
  f =

       0

thus evaluates the objective function at the starting point (or at whichever vector represented by prob.x if you changed it since the creation of prob.) Of course, prob.x is a Matlab placeholder just like any other. You may use it as a variable and you are free to change its value. Typically, an algorithm would update prob.x using, e.g., prob.x = prob.x + step * direction. The cuter_obj function is flexible in the sense that it will also return the objective gradient at the given point if you request it:

  matlab>> [f,g] = cuter_obj( prob.x );

The vector g now contains the gradient of f at prob.x`. We did not print it since in this example it is a vector of size 751 but we can view a few elements:

  matlab >> g(503:507)

  ans =

      0.0004
      0.0008
      0.0012
      0.0016
      0.0020

Evaluating constraint values is done with the cuter_cons function. This function can be used to evaluate all constraints at once or a single constraint by using one of the following forms:

  matlab>> c = cuter_cons( prob.x );         % Evaluate all constraints
  matlab>> c12 = cuter_cons( prob.x, 12 );   % Evaluate the 12-th constraint
  matlab>> c(12)

  ans =

     -0.5950

  matlab>> c12

  c12 =

     -0.5950

Again, this function is flexible and returns the constraints Jacobian or constraint gradient if you request it:

  matlab>> [c, J] = cuter_cons( prob.x );           % Obtain c(x) and Jacobian
  matlab>> size(J)

  ans =

     500   751

  matlab>> [c12, gc12] = cuter_cons( prob.x, 12 );  % Obtain c12(x) and gradient
  matlab>> size(gc12)

  ans =

     751     1

  matlab>> sum( gc12 ~= J(12,:)' )

  ans =

       0

Most other functions available follow the same pattern as cuter_obj and cuter_cons. They all share the prefix cuter_ to distinguish them from other functions which might be available in your environment and avoid name resolution conflicts. They are summarized in the following table.

Matlab Function CUTEr Library Functions Purpose
cuter_dims udimen / cdimen Discover problem dimensions (may be called before cuter_setup())
cuter_setup usetup / csetup Setup problem data structure
cuter_varnames varnames Discover variable names
cuter_connames connames Discover constraint names
cuter_obj uofg / cofg Evaluate objective function value and its gradient if requested
cuter_objcons cfn Evaluate objective and constraints
cuter_cons ccfg / ccifg Evaluate constraint bodies and their gradients if requested. Evaluate a single constraint value and its gradient if requested.
cuter_scons ccfsg / ccifsg Evaluate constraint bodies and Jacobian in sparse format. Evaluate a single constraint value and its gradient as a sparse vector.
cuter_lagjac cgr Evaluate Jacobian and gradient of either objective or Lagrangian.
cuter_slagjac csgr Evaluate Jacobian in sparse format and gradient of either objective or Lagrangian as a sparse vector.
cuter_jprod cjprod Evaluate the matrix-vector product between the Jacobian and a vector.
cuter_jtprod cjprod Evaluate the matrix-vector product between the transpose Jacobian and a vector.
cuter_hess udh / cdh Evaluate the Hessian matrix of the Lagrangian, or of the objective if the problem is unconstrained. The matrix is stored as a dense array.
cuter_ihess udh / cdh Evaluate the Hessian matrix of the i-th problem function (i=0 is the objective function), or of the objective if problem is unconstrained.
cuter_hprod uprod / cprod Evaluate the matrix-vector product between the Hessian of the Lagrangian (or the objective if unconstrained) and a vector.
cuter_gradhess ugrdh / cgrdh Evaluate the gradient of either the objective or the Lagrangian, the Jacobian (or its transpose) and the Hessian of the Lagrangian in dense format.
cuter_sphess ush / csh Evaluate the Hessian matrix of the Lagrangian, or of the objective if the problem is unconstrained, in sparse format.
cuter_spihess ush / csh Evaluate the Hessian matrix of the i-th problem function (i=0 is the objective function), or of the objective if problem is unconstrained, in sparse format.

Getting Help

At all times, (a bit of) help is available on the commands related to the CUTEr-Matlab interface using the Matlab help command. For example:

  matlab>> help cuter_cons
    Return constraint bodies and Jacobian if requested.
    or return a single constraint value and its gradient if requested
    Usage:  c = cuter_cons(x)    or  [c,J]   = cuter_cons(x)
           ci = cuter_cons(x,i)  or  [ci,gi] = cuter_cons(x,i)

Interfacing with Matlab Solvers

Existing solvers typically require the user to provide a function to evaluate the objective function, a function to evaluate the gradient, and so forth. Using the CUTEr/Matlab interface, it is easy to write a simple wrapper that may be passed to an existing solver. Here is an example, which returns the objective value, or the objective value and gradient, or the objective value, gradient and Hessian.

  function [varargout] = obj( x );
  % Evaluate objective function
  % Usage:       f = obj(x)      evaluates function value only
  %          [f,g] = obj(x)  evaluates function value and gradient
  %        [f,g,H] = obj(x)  evaluates function value, gradient and Hessian

  if nargout > 3
      error( 'obj: too many output arguments' );
  end

  if nargout == 1
      % Compute objective function value
      varargout{1} = cuter_obj(x);
  end

  if nargout > 1 
      % Gradient is requested
      [varargout{1}, varargout{2}] = cuter_obj(x);
      if nargout > 2
          % Hessian matrix is requested
          varargout{3} = cuter_sphess(x);
      end
  end

Tips

In order to avoid exiting and restarting Matlab every time you wish to work on a different problem, the bang operator in Matlab lets you execute a shell command without exiting. You can use this to decode a problem directly from the Matlab prompt:

  matlab>> !runcuter --package mx --decode BRAINPC1

Troubleshooting

  1. Gfortran is not officially supported to create MEX files and often results in a hard MATLAB crash. The solution is to use g95 instead. Binaries and source are available from http://www.g95.org/.
  1. At least on Mac Intel machines, the options to mex need to be changed to avoid messages of the form
   dyld: lazy symbol binding failed: Symbol not found: __gfortran_st_open
     Referenced from: /Users/dpo/tmp/test-cuter/mcuter.mexmaci
     Expected in: flat namespace

followed by a hard Matlab crash.

In the same directory as the mex executable, there is a shell script named mexopts.sh. Copy this file to ~/.matlab/R2008A/ (or whichever your current revision is) and edit it. The relevant section for Mac Intel machines is maci. Change the line

      LD="$CC"

to

      LD="$FC"

You will then need to rebuild your MEX file.