Calling NAG Library Routines from Scilab

NAG Technical Report 5/2009

Start of report Skip to Example 1 Skip to Example 2 Skip to Example 3 Skip to Example 4 Reference Page

3.5. Example 5

Pseudorandom Numbers taken from the Normal Distribution

We show how to call a routine from the NAG Fortran Library: G05LAF Although we are calling a NAG Fortran Library routine, we must use C code for the glue between Scilab and NAG. In this example we will demonstrate calling a random number generator which returns a vector of pseudorandom numbers generated from the Normal distribution. We will include an option to generate either a repeatable or non-repeatable sequence of random numbers.

Function prototype from the NAG Fortran Library Manual
C interface
Compiling with the Builder Script
Calling the function

Function prototype from the NAG Fortran Library Manual
According to the Fortran Library Manual, the prototype for function G05LAF looks like this:
```
  SUBROUTINE G05LAF(XMU, VAR, N, X, IGEN, ISEED, IFAIL)
  integer N, IGEN, ISEED(4), IFAIL
  double precision XMU, VAR, X(N)
```
See the NAG FORTRAN Library Manual [6] for detailed argument descriptions.

C interface

Here is the source code of our C interface nag_intext5.c:

/* Example 5
   =========

   Shows how to implement a Scilab wrapper calling a NAG routine from the
   Fortran library */

#include "stack-c.h"
#undef Complex
#define Complex NagComplex
#define SciComplex doublecomplex

#include "stdio.h"

... (other header file definitions missing here which are in the source code) ...

int nag_intext5(char *fname)
{
  // to call this function in scilab use:
  // [x,ifail] = nag_psdrnd_norm_fun(n,xmu,var,repeatable)
  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 n, i, min, max;
  int repeatable;

  Integer igen, ifail, ndim;
  Integer *iseed=0;

  double xmu, var;
  double *x;

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

  GetRhsVar(1, "i", &m1, &n1, &l1); // input n
  GetRhsVar(2, "d", &m2, &n2, &l2); // input xmu
  GetRhsVar(3, "d", &m3, &n3, &l3); // input var
  GetRhsVar(4, MATRIX_OF_BOOLEAN_DATATYPE, &m4, &n4, &l4); // repeatable

  // there are also other datatypes for integers, real, string etc.
  // MATRIX_OF_INTEGER_DATATYPE on the istk
  // MATRIX_OF_RATIONAL_DATATYPE on the sstk
  // MATRIX OF DOUBLE DATATYPE on the stk
  // MATRIX_OF_COMPLEX_DATATYPE on the zstk
  // MATRIX_OF_BOOLEAN_DATATYPE on the istk
  // STRING_DATATYPE on the cstk


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

  // Read in n
  //*============*
  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 output array x
  // ==================================
  if ( !(x = malloc(n*sizeof(double))))
    {
      sciprint("Allocation failure\n");
      Error(999); goto END;

    }


  // Read in xmu
  //===========
  if (m2!= 1 || n2!=1)
    {
      sciprint("%s: Dimension should be 1x1 for arg 2\r\n",fname);
      sciprint("n = %i", n);
      Error(999); goto END;
    }
  else
    {
      xmu = *stk(l2);
    }

  // Read in var
  //==============
  if (m3!=1 || n3!=1)
    {
      sciprint("%s: Dimension should be 1x1 for arg 3\r\n",fname);
      Error(999); goto END;
    }
  else
    {
      var = *stk(l3);
    }

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


  //***********************************
  // End of reading in input variables
  //***********************************


  if ( !(iseed = malloc(4*sizeof(Integer))))
    {
      sciprint("Allocation failure\n");
      Error(999); goto END;

    }

  // Call initialisation routine for this rng
  // ========================================
  igen = 0; // Basic generator

  iseed[0] = 23;
  iseed[1] = 4;
  iseed[2] = 10;
  iseed[3] = 1985;

  // Initialise the rng either repeatably or not repeatably
  if (repeatable)
    {
      g05kbf_(&igen, iseed);
    }
  else
    {
      g05kcf_(&igen, iseed);
    }

  // Call random number generator based on Normal PDF
  //* ================================================
  ifail = -1;

  g05laf_(&xmu, &var, &ndim, x, &igen, iseed, &ifail);

  // Print NAG error message if routine fails
  //=========================================
  if (ifail != 0)
    sciprint("g05laf failed");

  // Create output variables for Scilab
  //====================================
  CreateVar(5, "d", &n, &m1, &l5); // x
  CreateVar(6, "i", &n1, &m1, &l6); // ifail

  // Assign scilab output variables
  //================================
  for (i=0;i<n;i++)
    {
      *stk(l5+i) = x[i];
    }

  *istk(l6) = (int) ifail;

  // Assign order for output variables
  //===================================
  LhsVar(1) = 5;
  LhsVar(2) = 6;

 END:

  // Deallocate memory
  //===================
  if (x) free(x);
  if (iseed) free(iseed);

  return(0);
}

Points to note about this code:

Omitted from the code shown above, but included in the source code, is the header file source code for the G05 routines used in the above interface, copied from the Nag Fortran header file nagmk21_cf_us.h found on the nag website www.nag.com under the 'Downloads -> C header files' section.
As is typical on Linux/UNIX systems, we need to add an underscore to the name of a Fortran routine called from C, e.g. g05laf_
In order to interface properly with the Fortran library all variables must be passed by reference. Any scalar variables in the routine call are preceded by an ampersand, but arrays are passed without ampersand as they are already pointers.
We have also included an option in the Scilab interface to choose whether or not the random sequence is repeatable or not, demonstrating the flexibility of writing one's own interface to Scilab.

Compiling with the Builder Script

To compile and link this function, you need to run the build script nag_builder5.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.

This builder script looks slightly different from the previous ones, as we are linking to the NAG Fortran library.
```
    exec nag_builder5.sce
    exec loader.sce
```
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_example5.sce.
```
    --> exec nag_example5.sce
    ifail  =
           0.

      x  =

           0.9074485
           1.1109788
           0.8802789
           1.4147883
           2.1821836
           3.5388501
           1.3572214
           3.0131266
           2.5570866
           2.1172563
```
Tip: If you get any error messages in Scilab check the troubleshooting section.

Start of report Skip to Example 1 Skip to Example 2 Skip to Example 3 Skip to Example 4 Reference Page