NAG Library Callback Functions in R (A More Complex Example)

Index

Introduction
Step 1: Wrapping the Library Routine
Step 2: Wrapping the Callback Function
Step 3: Compiling the Code for Use in R
Step 4: Writing R Code to call the New DLL
Amending the Example to Call the Fortran Library
Additional Information

Introduction

This document expands on the information given in NAG Library Callback Functions in R, and highlights the additional information needed to call one of the more complicated NAG optimisation routines, namely e04ucc from within R. As before, we concentrate on the NAG C Library first and then explain what additional steps are required to call the equivalent Fortran routine (E04UCF).

The full source of all the wrappers for calling both the C and Fortran versions of E04UC are available from the NAG website.

Back to Top

Step 1: Wrapping the Library Routine

Writing the wrapper for the library routine can be broken down into a number of tasks:

  1. Develop the interface for the wrapper.

    The only required feature of the interface is that it must return an object of type SEXP. For this example we are using:

    SEXP we04ucc(SEXP psize, SEXP a, SEXP lower, SEXP upper, SEXP objfun, 
             SEXP confun, SEXP initial, SEXP optlist, SEXP rho);
    The NAG routine being wrapped, e04ucc, has 15 arguments; seven input arguments (n,nclin,ncnlin,a,tda,bl and bu), four input / output arguments (x,options,comm and ifail), two output arguments (objf and g) and two callback functions (objfun and confun). In the prototype for we04ucc we have combined n, nclin and ncnlin into an R vector of length 3, called psize. The value of a will be held in an R matrix of the same name. The argument tda is not required as a is being supplied by the user in the form of an R matrix and hence its size can be readily ascertained. The remaining two input arguments, bl and bu hold the lower and upper limits for the constraints. In our prototype we have converted these vectors into the lists: lower and upper. The first element of these lists correspond to the limits for the simple constraints, the second elements the limits for the linear constraints and the third element the limits for the non-linear constraints. The three elements of lower are then combined to form the vector bl in the wrapper prior to calling the NAG routine. Similarly the three elements of upper are combined to form bu.

    Of the input / output arguments, the value of x on entry is held in initial and the optional parameters usually supplied via the Nag_E04_Opt structure options are instead supplied in the list optlist. The output arguments as well as the value of x on exit and the optional output arguments usually returned in options are returned in a list.

    The callback functions objfun and confun can be supplied either as an R function, or the body of an R function.

    When supplying optional parameters to our interface, the optlist list is used, Calling the R wrapper with a value of optlist equal to:

    optlist = list(list=TRUE,"con_deriv"=FALSE,"print_level"=0,"optim_tol"=1.0e-8)
    for example, is equivalent to calling the C Library routine with: 
    options.list = Nag_TRUE;
options.print_level=Nag_NoPrint;
options.optim_tol=1.0e-4;
    The names of the elements of the optlist list correspond to the names in the options structure. In the case of enumerated types (as with print_level above), an integer value should be provided, starting from 0. The ordering of the enumeration is the same as given in the "Constraint" section of the documentation for that particular parameter. Therefore for print_level; 0 = Nag_NoPrint, 1 = Nag_Soln, 2 = Nag_Iter etc.
  2. Convert the R objects into standard C variable types that can be understood by the NAG Library routine.

    The checking and conversion of the R objects into standard C variables types can be carried out in similar manner to that described in the we04abc example. As some tasks occur often, for example, checking and converting an R vector into a C vector, a number of utility functions were used:

    • nag_extractDoubleVector : for checking and extracting a double precision vector from an R object.
    • nag_extractIntegerVector: for checking and extracting an integer vector from an R object.
    • nag_extractDoubleMatrix : for checking and extracting a double precision matrix from an R object.

    all of these utility routines return a numeric code if an error occurs, therefore the code snippet:

    if ((failed = nag_extractDoubleVector(n,initial,x))) {
  we04ucc_clean_up(p,ca,bu,bl,x,g,poptions);
  sprintf(err_msg,
           "'initial' must be a numeric vector of length %ld (code = %ld)",
                 (long) n, (long) failed);
  error(err_msg);
}
    can be used to extract the n initial values from the R object initial into the C array x. The memory for x needs to have been allocated prior to calling nag_extractDoubleVector. If an invalid value for initial was supplied, as indicated by an non-zero value for failed, then an error message is displayed. The routine we04ucc_clean_up ensures that any memory allocated in the wrapper is released before returning control to R. The code for these utilities is included with the other example source for this document.

    No lists were used in the interface of our previous example (e04abc), therefore we will go into a little more detail on how to handle lists in R. To check whether an R object is a list, the following can be used:

    if (!isNewList(lower)) {
  error("'lower' must be a list with three elements");
}
    Note the use of isNewList rather than isList. The function isList is supplied with the R headers, but does not check for an R object of type "list".

    Each element of a list can be accessed using the VECTOR_ELT function, so:

    SEXP lower_s;
lower_s = VECTOR_ELT(lower,i);
    sets the R object lower_s equal to the ith element of the list lower. The elements of the list are referenced from 0. In a similar manner the following loops through each element of the list optlist and prints out the name of the element. 
    SEXP names;
int i;

PROTECT(names = getAttrib(optlist,R_NamesSymbol));
for (i = 0; i < length(optlist); i++) {
  this_name = CHAR(STRING_ELT(names,i));
  printf("%s\n", this_name);
}
UNPROTECT(names);

The remaining steps in writing the wrapper for the library routine, namely:

  • Prepare any relevant user supplied information into a form that can be passed through the NAG Library routine to the callback function wrapper.
  • Call the NAG routine.
  • Check the error returns from the NAG Library.
  • Return the results from the NAG function to the user.

are carried out in a similar manner to that described in the example for e04abc.

Back to Top

Step 2: Wrapping the Callback Function

As described in NAG Library Callback Functions in R there are two ways that the callback function can be supplied to the C wrapper; as an R function, or as the body of an R function (via the R command body). The basic approach for doing this is the same for the two callback functions required by e04ucc (objfun and confun) as it was for the single callback function used by e04abc. The only difference is the number of parameters the callback functions require.

  1. Using the body of an R function

    The simplest method of interfacing to an R function is when the body of the function is used. In this case increasing the number of parameters can be achieved by increasing the number of calls to the defineVar function, for example:

    defineVar(install("n"), nag_i2r(1,&n), rho);
defineVar(install("x"), nag_d2r(n,x), rho);
defineVar(install("objgrd"), nag_d2r(n,objgrd), rho);
defineVar(install("flag"), nag_i2r(1,&comm->flag), rho);
defineVar(install("first"), nag_l2l(1,&comm->first), rho);

PROTECT(ans = eval(objfun,rho));
    sets up the five arguments required by objfun (n, x, objgrd, flag and first) before calling objfun and saving the result in an R object called ans. As in the previous example, we make use of the utility functions nag_i2r, nag_d2r and nag_l2l to convert C variables into R objects.

  2. Using the function itself

    In the example for e04abc, when the callback function had been supplied as an R function (as opposed to its body) we made use of the lang2 function. The lang2 function is supplied with the R header files, however it is restricted to calling routines with a single argument. As well as lang2 the R headers include routines for calling functions with two and three arguments, called lang3 and lang4 respectively. If more arguments are required the user must write their own. Fortunately it is straight forward to extend the code for lang4 to cope with more arguments. Utility routines for up to seven arguments, based on lang4, called lang5 up to lang8, have been included in the example source for this document. To call the objective function, objfun, when objfun is an actual function, you could therefore use:

    PROTECT(R_fcall = lang6(objfun,
                        nag_i2r(1,&n),
                        nag_d2r(n,x), 
                        nag_d2r(n,objgrd), 
                        nag_i2r(1,&comm->flag), 
                        nag_l2l(1,&comm->first)));

PROTECT(ans = eval(R_fcall,rho));
  3. Using a callback function that returns a list

    Both the objfun and confun callback functions must return more than one object. In the case of objfun this is the value of the objective function and the value of the derivative, and for confun it is the value of the non-linear constraints and their derivative. In order to facilitate this, both objfun and confun return a list.

    As we can see from above, it makes no difference to the eval function whether the callback function returns a scalar, vector, matrix or list, either way the result is returned and stored in an R object of type SEXP. It is how we process this object in the wrapper that needs to change. One way of dealing with the returned list would be to assume that the first element was, say, the value of the objective function and the second the value of the derivatives. The list could therefore be accessed via the VECTOR_ELT function as before, for example:

    SEXP tobj, tobjgrd;

PROTECT(tobj = VECTOR_ELT(lower,0));
PROTECT(tobjgrd = VECTOR_ELT(lower,1));
    However, in this example, we have used a utility function, called nag_getListElement, which extracts an element from the list based on its name. The code snippet: 
    SEXP tobj, tobjgrd;

PROTECT(tobj = nag_getListElement(ans, "objf"));
PROTECT(tobjgrd = nag_getListElement(ans, "objgrd"));
    takes the value of the objective function from the element of the list called "objf", and the derivatives from the one called "objgrd", irrespective of the order they appear in the list. The nag_getListElement function searches the names of the elements until it finds the specified value (str in the following code snippet). If the required name is not found then a null value (R_NilValue) is returned: 
    SEXP elmt, names;

elmt = R_NilValue;
names = getAttrib(list, R_NamesSymbol);

for (i = 0; i < length(list); i++) {
  if (!strcmp(CHAR(STRING_ELT(names,i)),str)) {
    elmt = VECTOR_ELT(list,i);
    break;
  }
}
return elmt;
    This function is based on the "getListElement" function described in Writing R Extensions.

A full copy of both versions of these wrappers is available in the example source.

Back to Top

Step 3: Compiling the Code for Use in R

The code for this example can be compiled in exactly the same way as described in step 3 of NAG Library Callback Functions in R.

Back to Top

Step 4: Writing R Code to Call the New DLL

The R front end for we04ucc can be written in the same way as described for we04abc. In this example, when supplying the body of the function, we have included a small amount of code to ensure that the list returned by the user has the correct names. The following function takes a list (l) and a vector of names (l.names) and checks whether there exists entries in l called l.names[i], for each i.

add.missing.names = function(l,l.names) {
  ## Adds any names from l.names that are missing from l, in order 
  ## (so first missing value is given the first unused name from l.names)
  ll.names = names(l)
  id.to.replace = (1:length(l))[!(ll.names%in%l.names)]
  replace.with = l.names[!(l.names%in%ll.names)]
  nmin = min(length(id.to.replace),length(replace.with))
  ll.names[id.to.replace[1:nmin]] = replace.with[1:nmin]
  names(l) = ll.names
  l
}

So, given a list ll:

add.missing.names(ll,c("objf","objgrd"))

will check that there are elements called "objf" and "objgrd" in the list. If not, assuming it does not exist, it will assign the name "objf" to the first element not called "objf" or "objgrd" and "objgrd" to the second (again, assuming it does not already exist).

Back to Top

Amending the Example to Call the Fortran Library

As before, calling the E04UCF routine from the NAG Fortran Library is very similar to calling e04ucc from the NAG C Library. The main difference, other than those already described in NAG Library Callback Functions in R, is the way that the optional parameters are specified.

Rather than passing optional parameters through a structure (as in the case of the options argument to e04ucc) the E04UC routine from the Fortran Library requires the user to call an option setting routine. Which routine is called depends on whether the thread safe routine E04UCA is being used, or the thread unsafe E04UCF is being used. The only difference between these two routines is that E04UCF holds the optional parameters in a Fortran common block, whereas E04UCA holds them in the routine arguments, LWSAV, IWSAV and RWSAV. In the example supplied with this document we call the thread safe routine E04UCA.

Prior to calling E04UCA the default values for the optional parameters must be initialised through a call to E04WBF (see the documentation for E04UCF for more details on this), for example:

E04WBF(rname,6,cwsav,80,&lcwsav,lwsav,&llwsav,iwsav,&liwsav,rwsav,&lrwsav,&ifail);

If you compare this function call to the documented interface for E04WBF you will notice that it has 12 arguments, as opposed to the documented 10. The two additional parameters (the second and fourth in the above call, taking values of 6 and 80), are hidden arguments and relate to the lengths of the character arguments rname and cwsav. In the C header files supplied with the NAG Fortran Library these variables are included in the prototypes and usually start with length_. So, for example, the prototype for E04WBF is:

extern void __stdcall E04WBF(
  CONST char rname[], 
  int length_rname,
  CONST char cwsav[], 
  int length_cwsav,
  CONST int *lcwsav,
  int lwsav[],
  CONST int *llwsav,
  CONST int iwsav[],
  CONST int *liwsav,
  CONST double rwsav[],
  CONST int *lrwsav,
  int *ifail
);

where length_rname is the length associated with the documented argument RNAME and length_cwsav the length associated with the documented argument CWSAV. These hidden arguments containing the string lengths are required for any Fortran routine that is being called from C and has character arguments, although their location in the prototype depends on the compiler and options used when compiling the Fortran code. Another NAG Fortran routine used in this example which is affected by this is X04ACF. This routine, which opens a file for reading or writing, has a prototype of

extern void __stdcall X04ACF(
  CONST int *nout,
  CONST char file[], 
  int length_file,
  CONST int *mode,
  int *ifail
);

compared to a documented interface of

X04ACF(NOUT,FILE,MODE,IFAIL)

Once again, the additional parameter length_file must contain the length of the filename file. When calling NAG Fortran routines with character arguments from C a couple of things should be noted:

  • The hidden variable should related to the number of characters being used and not the length of the character variable. For example, if the filename file was defined as 
    char *file[100];
strncpy(file,"test.txt",100);
    then the length of file, length_file should be set to 8 (as "test.txt" has 8 characters) and not 100 (even though file has been allocated a size of 100).
  • Some of the NAG Fortran Library routines take vectors of strings as arguments. One such example is E04WBF which expects CWSAV to be a vector of LCWSAV strings, with each string being 80 characters long, i.e.: 
    INTEGER      LCWSAV
CHARACTER*80 CWSAV(LCWSAV)
    If calling this routine from C, you would need to give the vector cwsav a length of 80*lcwsav, i.e. 
    char *cwsav;
cwsav = (char *) malloc(80*lcwsav*sizeof(char));
    and the hidden variable length_cwsav would have a value of 80.

Back to Top

Additional Information

Back to Top