NAG Technical Report 5/2009

Calling NAG Library Routines from Scilab


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

3.3. Example 3

eigvals

Complex Eigenvalalues and Eigenvectors of a Complex General Matrix (f02gcc)

Here we show how to write an interface to a NAG C Library function with complex input and output matrix variables. The function we will demonstrate with is f02gcc which computes the eigenvalues and eigenvectors of a complex system of equations.

Contents

  1. Function prototype from the NAG C Library Manual
  2. C Interface
  3. Compiling with the builder script
  4. Example of calling the function
  1. Function prototype from the NAG C Library Manual

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

         #include <nag.h>
    
         #include <nagf02.h>
    
         void nag_complex_eigensystem_sel(Nag_Select_Eigenvalues crit, Integer n, Complex a[], Integer tda, double wl,
                                          double wu, Integer mest, Integer *m, Complex w[], Complex v[], Integer tdv, 
                                          NagError *fail)
    
    This function takes an input matrix of type Complex amongst other inputs and computes the selected eigenvalues and their corresponding eigenvectors.
  2. C Interface

    Here is the source code of our interface written in C nag_intext3.c:

        /* Example 3
           =========
    
          Shows how to implement a Scilab wrapper calling in Complex matrix
          arguments and characters */
    
        #include "stack-c.h"
        #undef Complex
        #define Complex NagComplex
        #include <nag.h>
        #include <nag_stdlib.h>
        #include <nagf02.h>
    
        #define SciComplex doublecomplex
    
        #include "stdio.h"
    
        // Macros to allow matrix notation in the code.
        // Here the matrices are accessed in row order due to working in C
        #define A_NAG(I,J) a_nag[(I)*n + J]
        #define AMAT(I,J) amat[I*n + J]
        #define V_NAG(I,J) v_nag[(I)*n + J]
        #define VMAT(I,J) vmat[(I)*n + J]
    
        int nag_intext3(char *fname)
        {
          // to call this function in scilab use:
          // [w,v,m,ifail] = nag_cmplx_eisys_sel_fun(crit,n,a,wl,wu,mest)
          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 m7,n7,l7;
          int m8,n8,l8;
          int m9,n9,l9;
          int m10, n10, l10;
    
          int n, tda, mest,i,j, tdv;
          Integer m;
          double wl, wu;
          char crit;
          SciComplex *amat=0, *vmat=0, *wmat=0;
          NagComplex *a_nag=0, *v_nag=0, *w_nag=0;
          Nag_Select_Eigenvalues crit_nag;
          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, "c", &m1, &n1, &l1); // input crit 
          GetRhsVar(2, "i", &m2, &n2, &l2); // input n 
          GetRhsVar(3, "z", &m3, &n3, &l3); // input a 
          GetRhsVar(4, "d", &m4, &n4, &l4); // input wl 
          GetRhsVar(5, "d", &m5, &n5, &l5); // input wu 
          GetRhsVar(6, "i", &m6, &n6, &l6); // input mest 
    
          // Allocate memory for a_nag and amat 
          if ( !(a_nag = NAG_ALLOC(m3*n3,NagComplex) ) ||
               !(amat = NAG_ALLOC(m3*n3,SciComplex)))
            {
              sciprint("Allocation failure\n");
              Error(999); goto END;
            }
    
          //******************************************************************
          // Read in input variables from Scilab and convert to NAG variables 
          //******************************************************************
    
          // Read in crit 
          //*============*
          if (m1!=1 || n1!=1)
            {
              sciprint("%s: Dimension should be 1x1 character for arg 1\r\n",fname);
              Error(999); goto END;
            }
          else
            {
              crit = *cstk(l1);
              if (crit == 'r' || crit == 'R')
                {
                  crit_nag = Nag_Select_RealPart;
                }
              else if( crit == 'm' || crit == 'M')
                {
                  crit_nag = Nag_Select_Modulus;
                }
              else
                {
                  sciprint("%s: Arg 1 should be either 'r' or 'm'. \r\n",fname);
                  Error(999); goto END;
                }
            }
    
          // Read in n 
          //===========
          if (m2!=1 || n2!=1)
            {
              sciprint("%s: Dimension should be 1x1 for arg 2\r\n",fname);
              Error(999); goto END;
            }
          else
            {
              n = *istk(l2);
            }
    
          // Read in amat 
          //==============
          if (m3!=n || n3!=n)
            {
              sciprint("%s: Dimension should be nxn for arg 3\r\n",fname);
              Error(999); goto END;
            }
          else
            {
              // Convert Complex data types Scilab to NAG
              // Matrices are stored in column order in Scilab
              for(i=0; i<n; i++)
                {
                  for(j=0; j<n; j++)
                    {
                      AMAT(i,j) = *zstk(l3 + i + j * n);
                      A_NAG(i,j).re = AMAT(i,j).r;
                      A_NAG(i,j).im = AMAT(i,j).i;
                    }
                }
            }
    
          // Read in wl 
          //============
          if (m4!=1 || n4!=1)
            {
              sciprint("%s: Dimension should be 1x1 for arg 4\r\n",fname);
              Error(999); goto END;
            }
          else
            {
              wl = *stk(l4);
            }
    
          // Read in wu 
          //============
          if (m5!=1 || n5!=1)
            {
              sciprint("%s: Dimension should be 1x1 for arg 5\r\n",fname);
              Error(999); goto END;
            }
          else
            {
              wu = *stk(l5);
            }
    
          // Read in mest 
          //==============
          if (m6!=1 || n6!=1)
            {
              sciprint("%s: Dimension should be 1x1 for arg 6\r\n",fname);
              Error(999); goto END;
            }
          else
            {
              mest = *istk(l6);
            }
          //***********************************
          // End of reading in input variables 
          //***********************************
    
          // Allocate memory for NAG output v_nag(n,mest) and w_nag 
          //========================================================
          if ( !(v_nag = NAG_ALLOC(n*mest,NagComplex) ))
            {
              sciprint("Allocation failure\n");
              Error(999); goto END;
            }
    
          if ( !(w_nag = NAG_ALLOC(MAX(1,n),NagComplex) ))
            {
              sciprint("Allocation failure\n");
              Error(999); goto END;
            }
    
          // Define the other arguments for the NAG routine 
          //================================================
          tda = n;
          tdv = mest;
    
          // Call NAG routine f02gcc 
          //=========================
          nag_complex_eigensystem_sel(crit_nag,n,a_nag,tda,wl,wu,mest,&m,w_nag,v_nag,tdv,&ifail);
    
    
          
          // Create output variables for Scilab 
          //====================================
          CreateVar(7, "z", &n, &mest, &l7); // vmat 
          CreateVar(8, "z", &n, &m1, &l8); // wmat 
          CreateVar(9, "i", &n1, &m1, &l9); // m 
          CreateVar(10, "i", &n1, &m1, &l10); // ifail 
    
          // Allocate memory for output vmat(n,m) and wmat 
          //===============================================
          if ( !(vmat = NAG_ALLOC(n*mest,SciComplex)) ||
               !(wmat = NAG_ALLOC(MAX(1,n),SciComplex)))
            {
              sciprint("Allocation failure\n");
              Error(999); goto END;
            }
    
          // Assign scilab output variables 
          //================================
    
          // Convert complex types NAG to Scilab
          for (i=0; i<n; i++)
            {
              for (j=0; j<mest; j++)
                {
                  VMAT(i,j).r = V_NAG(i,j).re;
                  VMAT(i,j).i = V_NAG(i,j).im;
                  *zstk(l7 +i + j*n) = VMAT(i,j);
                }
            }
          for (i=0; i<n; i++)
            {
              wmat[i].r = w_nag[i].re;
              wmat[i].i = w_nag[i].im;
              *zstk(l8+i) = wmat[i];
            }
          *istk(l9) = m;
          *istk(l10) = ifail.code;
    
          // Print NAG error message if routine fails 
          //==========================================
          if (ifail.code)
            {
              sciprint("ifail message %s",ifail.message);
            }
    
          // Assign order for output variables 
          //===================================
          LhsVar(1) = 8;
          LhsVar(2) = 7;
          LhsVar(3) = 9;
          LhsVar(4) = 10;
    
          // Free allocated memory 
          //=======================
         END:
          if (a_nag) NAG_FREE(a_nag);
          if (amat) NAG_FREE(amat);
          if (v_nag) NAG_FREE(v_nag);
          if (vmat) NAG_FREE(vmat);
          if (w_nag) NAG_FREE(w_nag);
          if (wmat) NAG_FREE(wmat);
    
          return(0);
        }
    
    

    Points to note about this code:

  3. Compiling with the Builder Script

    To compile and link this function, you can run the build script nag_builder3.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_builder3.sce
        exec loader.sce
    
  4. 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_example3.sce.

        --> exec nag_example3.sce
          v  = 
          
          - 0.3865491 + 0.1732346i  - 0.0356136 - 0.1782180i
          - 0.3539288 + 0.4528810i    0.1263743 + 0.2666324i
            0.6123701                 0.0129326 - 0.2965682i
          - 0.0859284 - 0.3283626i    0.8898240
    
          w  =
    
          - 5.0000335 + 2.0060272i  
          3.0022643 - 3.9998187i  
          7.9981945 - 0.9963651i  
          - 6.0004253 - 6.9998434i  
    

    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 4   Skip to Example 5   Reference Page