The Xpress NonLinear API Functions

Topics covered in this chapter:

Instead of writing an extended MPS file and reading in the model from the file, it is possible to embed Xpress NonLinear directly into your application, and to create the problem, solve it and analyze the solution entirely by using the Xpress NonLinear API functions. This example uses the C header files and API calls. We shall assume you have some familiarity with the Xpress Optimizer API functions.

The structure of the model and the naming system will follow that used in the previous section, so you should read the chapter first.

Header files

The header file containing the Xpress NonLinear definitions is xslp.h. This must be included together with the Xpress Optimizer header xprs.h, where xprs.h must come first.

#include "xprs.h"
#include "xslp.h"

Initialization

Xpress NonLinear and Xpress Optimizer both need to be initialized, and an empty problem created. All Xpress NonLinear functions return a code indicating whether the function completed successfully. A non-zero value indicates an error. For ease of reading, we have for the most part omitted the tests on the return codes, but a well-written program should always test the values.

XPRSprob mprob;
XSLPprob sprob;

if (ReturnValue=XPRSinit(NULL)) goto ErrorReturn;
if (ReturnValue=XSLPinit()) goto ErrorReturn;
if (ReturnValue=XPRScreateprob(&mprob)) goto ErrorReturn;
if (ReturnValue=XSLPcreateprob(&sprob, &mprob)) goto ErrorReturn;

Callbacks

It is good practice to set up at least a message callback, so that any messages produced by the system appear on the screen or in a file. The XSLPsetcbmessage function sets both the Xpress NonLinear and Xpress Optimizer callbacks, so that all messages appear in the same place.

XSLPsetcbmessage(sprob, XSLPMessage, NULL);

void XPRS_CC XSLPMessage(XSLPprob my_prob, void *my_object, char *msg, int len,
       int msg_type)
{
  switch (msg_type) {
  case 4: /* error */
  case 3: /* warning */
  case 2: /* dialogue */
  case 1: /* information */
    printf("%s\n", msg);
    break;
  default: /* exiting */
    fflush(stdout);
    break;
  }
}

This is a simple callback routine, which prints any message to standard output.

Creating the linear part of the problem

The linear part of the problem, and the definitions of the rows and columns of the problem are carried out using the normal Xpress Optimizer functions.

#define MAXROW 20
#define MAXCOL 20
#define MAXELT 50
  int nRow, nCol, nSide, nRowName, nColName;
  int Sin, Cos;
  char RowType[MAXROW];
  double RHS[MAXROW], OBJ[MAXCOL], Element[MAXELT];
  double Lower[MAXCOL], Upper[MAXCOL];
  int ColStart[MAXCOL+1], RowIndex[MAXELT];
  char RowNames[500], ColNames[500];

In this example, we have set the dimensions by using #define statements, rather than working out the actual sizes required from the number of sides and then allocating the space dynamically.

  nSide = 5;
  nRowName = 0;
  nColName = 0;

By making the number of sides a variable (nSide) we can create other polygons by changing its value.

It is useful – at least while building a model – to be able to see what has been created. We will therefore create meaningful names for the rows and columns. nRowName and nColName count along the character buffers RowNames and ColNames.

  nRow = nSide-2 + (nSide-1)*(nSide-2)/2 + 1;
  nCol = (nSide-1)*2 + 2;
  for (i=0; i<nRow; i++) RHS[i] = 0;

The number of constraints is:

nSide-2	for the relationships between adjacent thetas.
(nSide-1)*(nSide-2)/2	for the distances between pairs of vertices.
1	for the OBJEQ non-linear "objective function".

The number of columns is:

nSide-1	for the thetas.
nSide-1	for the rhos.
1	for the OBJX objective function column.

We are using "C"-style numbering for rows and columns, so the counting starts from zero.

  nRow = 0;
  RowType[nRow++] = 'E'; /* OBJEQ */
  nRowName = nRowName + 1 + sprintf(&RowNames[nRowName], "OBJEQ");
  for (i=1; i<nSide-1; i++) {
    RowType[nRow++] = 'G'; /* T2T1 .. T4T3 */
    RHS[i] = 0.001;
    nRowName = nRowName + 1 + sprintf(&RowNames[nRowName], "T%dT%d", i+1, i);
  }

This sets the row type indicator for OBJEQ and the theta relationships, with a right hand side of 0.001. We also create row names in the RowNames buffer. Each name is terminated by a NULL character (automatically placed there by the sprintf function). sprintf returns the length of the string written, excluding the terminating NULL character.

  for (i=1; i<nSide-1; i++) {
    for (j=i+1; j<nSide; j++) {
      RowType[nRow] = 'L';
      RHS[nRow++] = 1.0;
      nRowName = nRowName + 1 + sprintf(&RowNames[nRowName], "V%dV%d", i, j);
    }
  }

This defines the L-type rows which constrain the distances between pairs of vertices. The right hand side is 1.0 (the maximum value) and the names are of the form ViVj.

  for (i=0; i<nCol; i++) {
    OBJ[i] = 0;        /* objective function */
    Lower[i] = 0;      /* lower bound normally zero */
    Upper[i] = XPRS_PLUSINFINITY; /* upper bound = infinity */
  }

This sets up the standard column data, with objective function entries of zero, and default bounds of zero to plus infinity. We shall change these for the individual items as required.

  nCol = 0;
  nElement = 0;
  ColStart[nCol] = nElement;
  OBJ[nCol] = 1.0;
  Lower[nCol++] = XPRS_MINUSINFINITY; /* free column */
  Element[nElement] = -1.0;
  RowIndex[nElement++] = 0;
  nColName = nColName + 1 + sprintf(&ColNames[nColName], "OBJX");

This starts the construction of the matrix elements. nElement counts through the Element and RowIndex arrays, nCol counts through the ColStart, OBJ, Lower and Upper arrays. The first column, OBJX, has the objective function value of +1 and a value of -1 in the OBJEQ row. It is also defined to be "free", by making its lower bound equal to minus infinity.

  iRow = 0
  for (i=1; i<nSide; i++) {
    nColName = nColName + 1 + sprintf(&ColNames[nColName], "THETA%d", i);
    ColStart[nCol++] = nElement;
    if (i < nSide-1) {
      Element[nElement] = -1;
      RowIndex[nElement++] = iRow+1;
    }
    if (i > 1) {
      Element[nElement] = 1;
      RowIndex[nElement++] = iRow;
    }
    iRow++;
  }

This creates the relationships between adjacent thetas. The tests on i are to deal with the first and last thetas which do not have relationships with both their predecessor and successor.

  Upper[nCol-1] = 3.1415926;

This sets the bound on the final theta to be π. The column index is nCol-1 because nCol has already been incremented.

  for (i=1; i<nSide; i++) {
    Lower[nCol] = 0.01;         /* lower bound */
    Upper[nCol] = 1;
    ColStart[nCol++] = nElement;
    nColName = nColName + 1 + sprintf(&ColNames[nColName], "RHO%d", i);
  }
  ColStart[nCol] = nElement;

The remaining columns – the rho variables – have only non-linear formulas and so they do not appear in the linear section except as empty columns. They are bounded between 0.01 and 1.0 but have no entries. The final entry in ColStart is one after the end of the last column.

  XPRSsetintcontrol(mprob, XPRS_MPSNAMELENGTH, 16);

If you are creating your own names – as we are here – then you need to make sure that Xpress Optimizer can handle both the names you have created and the names that will be created by Xpress NonLinear. Typically, Xpress NonLinear will create names which are three characters longer than the names you have used. If the longest name would be more than 8 characters, you should set the Xpress Optimizer name length to be larger – it comes in multiples of 8, so we have used 16 here. If you do not make the name length sufficiently large, then the XPRSaddnames function will return an error either here or during the Xpress NonLinear "construct" phase.

  XPRSloadlp(mprob, "Polygon", nCol, nRow, RowType, RHS, NULL,
    OBJ, ColStart, NULL, RowIndex, Element, Lower, Upper);

This actually loads the model into Xpress Optimizer. We are not using ranges or column element counts, which is why the two arguments are NULL.

  XPRSaddnames(mprob, 1, RowNames, 0, nRow-1);
  XPRSaddnames(mprob, 2, ColNames, 0, nCol-1);

The row and column names can now be added.

Adding the non-linear part of the problem

Be warned – this section is complicated, but it is the most efficient way – from SLP's point of view – to input formulae. See the next section for a much easier (but less efficient) way of inputting the formulae directly.

#define MAXTOKEN 200
#define MAXFORM 20
...
  int Sin, Cos;
  FormulaStart[MAXFORM];
  Type[MAXTOKEN];
  double Value[MAXTOKEN];

The arrays for the non-linear part can often be re-used from the linear part. The new arrays are FormulaStart for the nonlinear formulas and Type and Value to hold the internal forms of the formulae.

  XSLPgetindex(sprob, XSLP_INTERNALFUNCNAMES, "SIN", &Sin);
  XSLPgetindex(sprob, XSLP_INTERNALFUNCNAMES, "COS", &Cos);

We will be using the Xpress NonLinear internal functions SIN and COS. The XSLPgetindex function finds the index of an Xpress NonLinear entity (character variable, internal or user function).

  nToken = 0;
  nForm = 0;
  RowIndex[nForm] = 0;
  FormulaStart[nForm++] = nToken;

For each nonlinear formula, the following information is required:

RowIndex	the index of the row.
FormulaStart	the beginning of the internal formula array for the nonlinear formula.

  for (i=1; i<nSide-1; i++) {
    Type[nToken] = XSLP_COL;
    Value[nToken++] = nSide+i+1;
    Type[nToken] = XSLP_COL;
    Value[nToken++] = nSide+i;
    Type[nToken] = XSLP_OP;
    Value[nToken++] = XSLP_MULTIPLY;
    Type[nToken] = XSLP_RB;
    Value[nToken++] = 0;
    Type[nToken] = XSLP_COL;
    Value[nToken++] = i+1;
    Type[nToken] = XSLP_COL;
    Value[nToken++] = i;
    Type[nToken] = XSLP_OP;
    Value[nToken++] = XSLP_MINUS;
    Type[nToken] = XSLP_IFUN;
    Value[nToken++] = Sin;
    Type[nToken] = XSLP_OP
    Value[nToken++] = XSLP_MULTIPLY;
    if (i>1) {
      Type[nToken] = XSLP_OP;
      Value[nToken++] = XSLP_PLUS;
    }
  }

This looks very complicated, but it is really just rather large. We are using the "reverse Polish" or "parsed" form of the formula for area. The original formula, written in the normal way, would look like this:
RHO2 * RHO1 * SIN ( THETA2 - THETA1 ) + .......
In reverse Polish notation, tokens are pushed onto the stack or popped from it. Typically, this means that a binary operation A x B is written as A B x (push A, push B, pop A and B and push the result). The first term of our area formula then becomes:
RHO2 RHO1 * ) THETA2 THETA1 - SIN *
Notice that the right hand bracket appears as an explicit token. This allows the SIN function to identify where its argument list starts – and incidentally allows functions to have varying numbers of arguments.

Each token of the formula is written as two items – Type and Value.
Type is an integer and is one of the defined types of token, as given in the xslp.h header file. XSLP_CON, for example, is a constant; XSLP_COL is a column.
Value is a double precision value, and its meaning depends on the corresponding Type. For a Type of XSLP_CON, Value is the constant value; for XSLP_COL, Value is the column number; for XSLP_OP (arithmetic operation), Value is the operand number as defined in xslp.h; for a function (type XSLP_IFUN for internal functions, XSLP_FUN for user functions), Value is the function number.
A list of tokens for a formula is always terminated by a token of type XSLP_EOF.

The loop writes each term in order, and adds terms (using the XSLP_PLUS operator) after the first pass through the loop.

  for (i=1; i<nSide-1; i++) {
    for (j=i+1; j<nSide; j++) {
      RowIndex[nForm] = iRow++;
      FormulaStart[nForm++] = nToken;

      Type[nToken] = XSLP_COL;
      Value[nToken++] = nSide+i;
      Type[nToken] = XSLP_CON;
      Value[nToken++] = 2;
      Type[nToken] = XSLP_OP;
      Value[nToken++] = XSLP_EXPONENT;
      Type[nToken] = XSLP_COL;
      Value[nToken++] = nSide+j;
      Type[nToken] = XSLP_CON;
      Value[nToken++] = 2;
      Type[nToken] = XSLP_OP;
      Value[nToken++] = XSLP_PLUS;
      Type[nToken] = XSLP_CON;
      Value[nToken++] = 2;
      Type[nToken] = XSLP_COL;
      Value[nToken++] = nSide+i;
      Type[nToken] = XSLP_OP;
      Value[nToken++] = XSLP_MULTIPLY;
      Type[nToken] = XSLP_COL;
      Value[nToken++] = nSide+j;
      Type[nToken] = XSLP_OP;
      Value[nToken++] = XSLP_MULTIPLY;
      Type[nToken] = XSLP_RB;
      Value[nToken++] = 0;
      Type[nToken] = XSLP_COL;
      Value[nToken++] = j;
      Type[nToken] = XSLP_COL;
      Value[nToken++] = i;
      Type[nToken] = XSLP_OP;
      Value[nToken++] = XSLP_MINUS;
      Type[nToken] = XSLP_IFUN;
      Value[nToken++] = Cos;
      Type[nToken] = XSLP_OP
      Value[nToken++] = XSLP_MULTIPLY;
      Type[nToken] = XSLP_OP;
      Value[nToken++] = XSLP_MINUS;
      Type[nToken] = XSLP_EOF;
      Value[nToken++] = 0;
    }
  }

This writes the formula for the distances between pairs of vertices. It follows the same principle as the previous formula, writing the formula in parsed form as:
RHOi 2 ^ RHOj 2 ^ + 2 RHOi * RHOj * ) THETAj THETAi - COS * -

XSLPloadformulass(sprob, nForm, RowIndex, FormulaStart, 1, Type, Value);

The XSLPloadformulas is the most efficient way of loading non-linear formulas into a problem. There is an XSLPaddformulas function which is identical except that it does not delete any existing nonlinear formulas first. There is also an XSLPchgformula function, which can be used to change individual nonlinear formulas one at a time. Because we are using internal parsed format, the "Parsed" flag in the argument list is set to 1.

Adding the non-linear part of the problem using character formulae

Provided that all entities – in particular columns and user functions – have explicit and unique names, the non-linear part can be input by writing the formulae as character strings. This is not as efficient as using the XSLPloadformulas() function but is generally easier to understand.

/* Build up nonlinear formulas */
/* Allow space for largest formula - approx 50 characters per side for area */
  FormBuffer = (char *) malloc(50*nSide);

We shall be using large formulae, so we need a character buffer large enough to hold the largest formula we are using. The estimate here is 50 characters per side of the polygon for the area formula, which is the largest we are using.

/* Area */

  BufferPos = 0;
  for (i=1; i<nSide-1; i++) {
    if (i > 1) {
      BufferPos = BufferPos + sprintf(&FormBuffer[BufferPos], " + ");
    }
    BufferPos = BufferPos + sprintf(&FormBuffer[BufferPos], "RHO%d * RHO%d *
                SIN ( THETA%d - THETA%d )", i+1, i, i+1, i);
  }
  XSLPchgformulatext(sprob, 0, nSide, FormBuffer);

The area formula is of the form:
(RHO2*RHO1*SIN(THETA2-THETA1) + RHO3*RHO2*SIN(THETA3-THETA2) + ... ) / 2
The loop writes the product for each consecutive pair of vertices and also puts in the "+" sign after the first one.

The XSLPchgformulatext function is a variation of XSLPchgformula but uses a character string for the formula instead of passing it as arrays of tokens. The arguments to the function are:

RowIndex	the index of the row.
FormBuffer	the formula, written in character form.

In this case, RowIndex is zero.

/* Distances */
  for (i=1; i<nSide-1; i++) {
    for (j=i+1; j<nSide; j++) {
      sprintf(FormBuffer, "RHO%d ^ 2 + RHO%d ^ 2 - 2 * RHO%d * RHO%d *
              COS ( THETA%d - THETA%d )", j, i, j, i, j, i);
      XSLPchgformulatext(sprob, iRow, FormBuffer);
      iRow++;
    }

This creates the formula for the distance between pairs of vertices and writes each into a new row.

Provided you have given names to any user functions in your program, you can use them in a formula in exactly the same way as SIN and COS have been used above.

Checking the data

Xpress NonLinear includes the function XSLPwriteprob which writes out a non-linear problem in text form which can then be checked manually. Indeed, the problem can then be run using the XSLP console program, provided there are no user functions which refer back into your compiled program. In particular, this facility does allow small versions of a problem to be checked before moving on to the full size ones.

XSLPwriteprob(sprob, "testmat", "");

The first argument is the Xpress NonLinear problem pointer; the second is the name of the matrix to be produced (the suffix ".mat" will be added automatically). The last argument allows various different types of output including "scrambled" names – that is, internally-generated names will be used rather than those you have provided. For checking purposes, this is obviously not a good idea.

Solving and printing the solution

XSLPmaxim(sprob, "");

The XSLPmaxim and XSLPminim functions perform a non-linear maximization or minimization on the current problem. The second argument can be used to pass flags as defined in the Xpress NonLinear Reference Manual.

XPRSwriteprtsol(mprob);

The standard Xpress Optimizer solution print can be obtained by using the XPRSwriteprtsol function. The row and column activities and dual values can be obtained using the XPRSgetlpsol function.

In addition, you can use the XSLPgetvar function to obtain the values of SLP variables – that is, of variables which are in nonlinear formulas, or which have nonlinear formulas. If you are using cascading (see the Xpress NonLinear reference manual for more details) so that Xpress NonLinear recalculates the values of the dependent SLP variables at each SLP iteration, then the value from XSLPgetvar will be the recalculated value, whereas the value from XPRSgetlpsol will be the value from the LP solution (before recalculation).

Closing the program

The XSLPdestroyprob function frees any system resources allocated by Xpress NonLinear for the specific problem. The problem pointer is then no longer valid. XPRSdestroyprob performs a similar function for the underlying linear problem mprob. The XSLPfree function frees any system resources allocated by Xpress NonLinear. You must then call XPRSfree to perform a similar operation for the optimizer.

  XSLPdestroyprob(sprob);
  XPRSdestroyprob(mprob);
  XSLPfree();
  XPRSfree();

If these functions are not called, the program may appear to have worked and terminated correctly. However, in such a case there may be areas of memory which are not returned to the system when the program terminates and so repeated executions of the program will result in progressive loss of available memory to the system, which will manifest iself in poorer performance and could ultimately produce a system crash.

Adding initial values

So far, Xpress NonLinear has started by using values which it estimates for itself. Because most of the variables are bounded, these initial values are fairly reasonable, and the model will solve. However, in general, you will need to provide initial values for at least some of the variables. Initial values are provided using the XSLPsetinitval function.