## 3.4. Example 4

### Evaluating a multi-dimensional integral (d01wcc)

In this example we will demonstrate two ways of providing a user-defined function to a NAG routine in Scilab. The function we shall use is d01wcc, which evaluates a multidimensional integral of the function f(x) where x is a multidimensional vector. The user-supplied function will evaluate this function at any x.

The two ways of providing this function to the NAG routine are either to write the function in C, or to write it as a Scilab function. Both approaches have their advantages and disadvantages. The C code is a lot easier to write and requires fewer lines, but in order to change it, each time you would need to rebuild the interface. Writing the function in Scilab code is slightly harder, as you need to write a second interface for the Scilab function in C, but it does give you the flexibility to change your function without having to rebuild the Scilab-NAG interface.

### Contents

1. Function prototype from the NAG C Library Manual
2. Writing the User-supplied Function in C
3. Writing the User-Supplied Function in Scilab
1. ### Function prototype from the NAG C Library Manual

According to the C Library Manual, the prototype for function d01wcc looks like this:

```    #include <nag.h>
#include <nagd02.h>

const double a[], const double b[], Integer *minpts, Integer maxpts, double eps,
double *finval, double *acc, Nag_User *comm, NagError *fail)
```
Note that the user-supplied function is named f, and is of type double. For more information about the argument list of nag_multid_quad_adapt_1, refer to the NAG C Library manual.

1. ### C Interface

Here is the source code of our interface written in C nag_intext4a.c with the user-supplied function also written in C:

```    /* Example 4a
=========

Shows how to implement a Scilab wrapper to a NAG function using a User
defined function written in C */

#include "stack-c.h"
#undef Complex
#define Complex NagComplex
#include <nag.h>
#include <nag_stdlib.h>

#include <nagd01.h>
#define SciComplex doublecomplex

#include "stdio.h"

static double f(Integer n, double z[], Nag_User *comm);

int nag_intext4a(char *fname)
{
// to call this function in scilab use:
int m1,n1,l1;
int m2,n2,l2;
int m3,n3,l3;
int m4,n4,l4;
int m5,n5,l5;
int m6,n6,l6;
int l7;
int l8;
int l9;
int l10;

int n, i, min, max;

Integer ndim, minpts, maxpts;
double eps, finval, acc;
double *a=0, *b=0;
Nag_User comm;
NagError ifail;

// define minimum and maximum left and right hand arguments
int minlhs=1, minrhs=6, maxlhs=4, maxrhs=6;
Nbvars = 0;
CheckRhs(minrhs, maxrhs);
CheckLhs(minlhs,maxlhs);

// Initialise ifail
INIT_FAIL(ifail);

GetRhsVar(1, "i", &m1, &n1, &l1); // input ndim
GetRhsVar(2, "d", &m2, &n2, &l2); // input a
GetRhsVar(3, "d", &m3, &n3, &l3); // input b
GetRhsVar(4, "i", &m4, &n4, &l4); // input minpts
GetRhsVar(5, "i", &m5, &n5, &l5); // input maxpts
GetRhsVar(6, "d", &m6, &n6, &l6); // input eps

//******************************************************************
// Read in input variables from Scilab and convert to NAG variables
//******************************************************************

//*============*
if (m1!=1 || n1!=1)
{
sciprint("%s: Dimension should be 1x1 character for arg 1\r\n",fname);
Error(999); goto END;
}
else
{
n = *istk(l1);
ndim = (Integer) n;
}

// Allocate memory for a(ndim) and b(ndim)
//========================================
if ( !(a = NAG_ALLOC(n,double) ) ||
!(b = NAG_ALLOC(n,double)))
{
sciprint("Allocation failure\n");
Error(999); goto END;
}

//===========
if (m2!= n || n2!=1)
{
sciprint("%s: Dimension should be ndim by 1 for arg 2\r\n",fname);
sciprint("n = %i", n);
Error(999); goto END;
}
else
{
for (i=0; i<n; i++)
{
a[i] = *stk(l2+i);
}
}

//==============
if (m3!=n || n3!=1)
{
sciprint("%s: Dimension should be ndim by 1 for arg 3\r\n",fname);
Error(999); goto END;
}
else
{
for(i=0; i<n; i++)
{
b[i] = *stk(l3+i);
}
}

//============
if (m4!=1 || n4!=1)
{
sciprint("%s: Dimension should be 1x1 for arg 4\r\n",fname);
Error(999); goto END;
}
else
{
min = *istk(l4);
minpts = (Integer)min;
}

//============
if (m5!=1 || n5!=1)
{
sciprint("%s: Dimension should be 1x1 for arg 5\r\n",fname);
Error(999); goto END;
}
else
{
max = *istk(l5);
maxpts = (Integer) max;
}

//==============
if (m6!=1 || n6!=1)
{
sciprint("%s: Dimension should be 1x1 for arg 6\r\n",fname);
Error(999); goto END;
}
else
{
eps = *stk(l6);
}
//***********************************
// End of reading in input variables
//***********************************

//*====================================================
&comm, &ifail);

// Print NAG error message if routine fails
//=========================================
if (ifail.code != NE_NOERROR)
sciprint("Ifail message is %s", ifail.message);

// Create output variables for Scilab
//====================================
CreateVar(7, "d", &n1, &m1, &l7); // finval
CreateVar(8, "d", &n1, &m1, &l8); // acc
CreateVar(9, "i", &n1, &m1, &l9); // minpts
CreateVar(10, "i", &n1, &m1, &l10); // ifail

// Assign scilab output variables
//================================
*stk(l7) = finval;
*stk(l8) = acc;
*istk(l9) = (int) minpts;
*istk(l10) = (int) ifail.code;

// Assign order for output variables
//===================================
LhsVar(1) = 7;
LhsVar(2) = 8;
LhsVar(3) = 9;
LhsVar(4) = 10;

END:

// Deallocate memory
//===================
if (a) NAG_FREE(a);
if (b) NAG_FREE(b);

return(0);
}

static double f(Integer n, double z[], Nag_User *comm)
{
// Function describing the integrand for the NAG routine used above
//==================================================================
double tmp_pwr;
tmp_pwr = z[1]+1.0+z[3];
return z[0]*4.0*z[2]*z[2]*exp(z[0]*2.0*z[2])/(tmp_pwr*tmp_pwr);
}

```

• The input function for the NAG routine has a specific interface that has to be adhered to, which is documented in the NAG C Library manual [5]
• The interface doesn't take in the input arguments of the function as its arguments, but only takes fname the name of the function in Scilab. In this example we have chosen fname to be the same as the NAG C library name.
• The Scilab function GetRhsVar is used to retrieve the input arguments of fname
• CreateVar creates a variable into which the NAG output is placed
• The interface doesn't actually return the outputs of the function fname - instead the Scilab array LhsVar stores the variable numbers that are stored on Scilab's stack stk()
2. ### Compiling with the Builder Script

To compile and link this function, you can run the build script nag_builder4a.sce which links this interface to the function name fname that will be seen in Scilab. This must be run in the same directory as the interface file, and when run successfully, generates a loader script loader.sce that needs to be executed in order to dynamically link the newly created library into Scilab.

```    exec nag_builder4a.sce
```
3. ### Calling the function

If the function built and the loader.sce file loaded the function without problem, then we can use the function within Scilab, either on the command-line or in a script. Here we run the example script, nag_example4a.sce.

```    --> exec nag_example4a.sce
finval = 0.5753623
```

Tip: If you get any error messages in Scilab check the troubleshooting section.

3. ### Writing the user-supplied functions in Scilab

1. #### C Interface for the Scilab code

The file nag_intext4b.c contains similar code to the example above.

The only thing that is different between the two source codes nag_intext4a.c and nag_intext4b.c is that in the second one, instead of the user-supplied function being written in C, we have an interface to the user-supplied function written in Scilab. This is so that any variable types can be converted and so that the interface to the NAG routine is correct. Below we just show the part of the file containing the new interface, where we have named the function fSci.

```    static double fSci(Integer n, double z[], Nag_User *comm);

static double fSci(Integer n, double z[], Nag_User *comm)
{
int mlhs, mrhs, ibegin;
int m_u, n_u, l_u;
int m_z, n_z, l_z;
int i;
double retval;
char name[] = "fUser";

// Put input variables for scilab function onto stack in the right order
// On input u (stored in space 7) is a dummy input to be overwritten with
// the output from fUser
// The inputs for fUser must be in the correct order in the last places of
// the stack.  They are then scrapped and the outputs are created starting
// at the same place as the first input variable.
m_u = 1;
n_u = 1;
m_z = (int) n;
n_z = 1;

CreateVar(7, "d", &m_u, &n_u, &l_u); // first input
CreateVar(8, "d", &m_z, &n_z, &l_z); // second input

// Write the input vector z onto the stack in the second input
for (i=0; i<m_z; i++)
{
*stk(l_z+i) = z[i];
}

ibegin = 7; // the first input variable on the stack
mlhs = 1;  // number of left hand side outputs to write onto the stack
mrhs = 2; // number of right hand side inputs to read off the stack
SciString(&ibegin,name,&mlhs,&mrhs); // Call scilab function 'name'

// output from scilab function is now on the stack so return the result
retval = *stk(l_u);

return retval;
}

```

• The User-defined function must still have the same interface as specified in the C Library manual and therefore we must write an interface to the Scilab function.
• The Scilab function fUser is not called directly, but through a Scilab function called SciString. The arguments for SciString are
• ibegin: the variable number of the first input argument. Note that the input arguments must be next to each other and in the order they are to passed into fUser with the first input argument being at ibegin.
• name: the name of the Scilab function.
• mlhs and mrhs: the number of input and output variables of the Scilab function. The input arguments are destroyed after use, and the output arguments are then put on the stack in order starting at ibegin.
2. ### User-supplied functions in Scilab

Source code of our Scilab function is in the example script nag_example4b.sce but could equally have been in a separate script. The source also demonstrates two different ways of writing the function.
```    deff('u=fUser(uu,z)',
['tmp_pwr = z(2)+1.0+z(4)';
'u = z(1)*4.0*z(3)*z(3)*exp(z(1)*2.0*z(3))/(tmp_pwr*tmp_pwr)']);
```
3. ### Compiling with the Builder Script

To compile and link this function, you need to run the build script nag_builder4b.sce which links this interface to the function name fname that will be seen in Scilab. This must be run in the same directory as the interface file, and when run successfully, generates a loader script loader.sce that needs to be executed in order to dynamically link the newly created library into Scilab.

```
exec nag_builder4b.sce
```      --> exec nag_example4b.sce