NAG Technical Report 2/2009

Calling NAG Library Routines from Java


Start of report   Skip to Example 1   Skip to Example 2   Skip to Example 3   Skip to Example 5   Skip to Example 6

4.4. Example 4

A nonlinear minimization routine, E04UCF

Here we show how to call a NAG Fortran Library routine which takes two user-supplied evaluation routines as arguments: subroutine E04UCF. The example demonstrates how to create Java arrays in C interface code. Thanks go to Sebastien Girard of Imperial College, London, for inspiration for this example.

Contents

  1. Subroutine specification from the NAG Fortran Library Manual
  2. Declaring the native function in our Java program
  3. Compiling the Java program
  4. Generating a header file for use by C
  5. Implementing the native function in C code which interfaces to Fortran
  6. Building the shareable library or DLL
  7. Running the program
  8. Quick summary of how to build the minimization solver example

  1. Subroutine specification from the NAG Fortran Library Manual
  2. According to the Fortran Library Manual, the prototype for subroutine E04UCF looks like this:

          SUBROUTINE E04UCF(N, NCLIN, NCNLN, LDA, LDCJ, LDR, A, BL, BU, CONFUN,
         1                  OBJFUN, ITER, ISTATE, C, CJAC, CLAMDA, OBJF, OBJGRD,
         2                  R, X, IWORK, LIWORK, WORK, LWORK, IUSER, USER, IFAIL)
                            INTEGER N, NCLIN, NCNLN, LDA, LDCJ, LDR, ITER,
         1                  ISTATE(N+NCLIN+NCNLN), IWORK(LIWORK), LIWORK, LWORK,
         2                  IUSER(*), IFAIL
                            DOUBLE PRECISION  A(LDA,*), BL(N+NCLIN+NCNLN),
         1                  BU(N+NCLIN+NCNLN), C(*), CJAC(LDCJ,*),
         2                  CLAMDA(N+NCLIN+NCNLN), OBJF, OBJGRD(N), R(LDR,N),
         3                  X(N), WORK(LWORK), USER(*)
                            EXTERNAL CONFUN, OBJFUN
    
    The subroutine E04UCF is designed to minimize an arbitrary smooth function f(x) subject to constraints (which may include simple bounds on the variables, linear constraints and smooth nonlinear constraints).

    The function f(x) to be minimized is evaluated by user-supplied routine argument OBJFUN to E04UCF, which is declared as

          SUBROUTINE OBJFUN(MODE, N, X, OBJF, OBJGRD, NSTATE, IUSER, USER)
                            INTEGER MODE, N, NSTATE, IUSER(*)
                            DOUBLE PRECISION X(N), OBJF, OBJGRD(N), USER(*)
    
    The routine OBJFUN may also need to evaluate the gradient of the objective function. In addition, a user-supplied routine argument CONFUN is needed to evaluate nonlinear constraints of the minimization problem:
          SUBROUTINE CONFUN(MODE, NCNLN, N, LDCJ, NEEDC, X, C, CJAC, NSTATE,
         1                  IUSER, USER)
                            INTEGER MODE, NCNLN, N, LDCJ, NEEDC(NCNLN), NSTATE,
         1                  IUSER(*)
                            DOUBLE PRECISION X(N), C(NCNLN), CJAC(LDCJ,N), USER(*)
    

    For full description of the roles of all routine arguments consult the E04UCF routine document in the NAG Fortran Library Manual.

  3. Declaring the native function in our Java program
  4. The NAG Fortran Library uses a simple integer, IFAIL, to return an error code (compare this with the NagError structure type used by the C Library).

    Although we are going to call a routine from the NAG Fortran Library, we still implement our interface code in C for convenience. Thus the full program will be a mixture of three languages - Java, C and Fortran.

      // Declaration of the Native (C) function
      private native int e04ucf(int n, int nclin, int ncnln,
                                double[] a, int lda, double[] bl, double[] bu,
                                String confun, String objfun,
                                double[] objf, double[] objgrd, double[] x);
    
    i.e. a method with return type int. Note that we choose not to pass all possible arguments - for example, the arrays cjac, clamda and istate are missing. We could include these arguments if we wanted the information they contain to be returned to Java; here we don't. Since we are also not using the ifail argument, we will use the int return value to send back any error code.

    As with Example 3, note that we cannot pass subroutine arguments directly from Java to C, and so here we just pass the names of methods via the String arguments confun and objfun.

  5. Compiling the Java program
  6. Here is the complete source code of our Java program Minimization.java.
    public class Minimization
    {
    
      // Declaration of C native function, NAG routine e04ucf
      private native int e04ucf(int n, int nclin, int ncnln,
                                double[] a, int lda, double[] bl, double[] bu,
                                String confun, String objfun,
                                double[] objf, double[] objgrd, double[] x);
    
      // An interface to e04uef, an option setting routine for e04ucf
      private native void e04uef(String option);
    
      static
      {
        System.loadLibrary("nagCJavaInterface");
      }
    
      /* A routine to evaluate the nonlinear constraint functions and
         Jacobian. This gets called from NAG routine e04ucf via the
         Java Native Interface. N.B. cjac is stored as a 1D array
         rather than 2D array for convenience. */
      private void confun(int mode, int ncnln, int n, int ldcj,
                          int[] needc, double[] x, double[] c,
                          double[] cjac, int nstate)
      {
        if (nstate == 1)
          {
            // First call to confun. Set all Jacobian elements to zero.
            // Note that this will only work when 'Derivative Level = 3'
            // (the default (see Section 11.2).
            for (int j=0; j<n; j++)
              {
                for (int i=0; i<ncnln; i++)
                  {
                    // Notice how we address the array cjac so that contents
                    // are in the order required by a 2D Fortran array.
                    cjac[i+j*ldcj] = 0.0;
                  }
              }
          }
    
        if (needc[0] > 0)
          {
            if (mode == 0 || mode == 2)
              {
                c[0] = x[0]*x[0] + x[1]*x[1] + x[2]*x[2] + x[3]*x[3];
              }
            if (mode == 1 || mode == 2)
              {
                cjac[0+0*ldcj] = 2.0e0*x[0];
                cjac[0+1*ldcj] = 2.0e0*x[1];
                cjac[0+2*ldcj] = 2.0e0*x[2];
                cjac[0+3*ldcj] = 2.0e0*x[3];
              }
          }
    
        if (needc[1] > 0)
          {
            if (mode == 0 || mode == 2)
              {
                c[1] = x[0]*x[1]*x[2]*x[3];
              }
            if (mode == 1 || mode == 2)
              {
                cjac[1+0*ldcj] = x[1]*x[2]*x[3];
                cjac[1+1*ldcj] = x[0]*x[2]*x[3];
                cjac[1+2*ldcj] = x[0]*x[1]*x[3];
                cjac[1+3*ldcj] = x[0]*x[1]*x[2];
              }
          }
      }
    
      /* A routine to evaluate the objective function and its gradient.
         This gets called from NAG routine e04ucf via the Java Native
         Interface */
      private double objfun(int mode, int n, double[] x,
                            double[] objgrd, int nstate)
      {
        double objf = 0.0;
        if (mode == 0 || mode == 2)
          {
            objf = x[0]*x[3]*(x[0]+x[1]+x[2]) + x[2];
          }
    
        if (mode == 1 || mode == 2)
          {
            objgrd[0] = x[3]*(2.0e0*x[0]+x[1]+x[2]);
            objgrd[1] = x[0]*x[3];
            objgrd[2] = x[0]*x[3] + 1.0e0;
            objgrd[3] = x[0]*(x[0]+x[1]+x[2]);
          }
        return objf;
      }
    
      // Main program
      public static void main(String args[])
      {
    
        Minimization nlp = new Minimization();
    
        // Pass the names of the constraint function and the objective
        // function evaluation routines.
        nlp.Solve("confun", "objfun");
      }
    
      private void Solve(String confunction, String objfunction)
      {
    
        // n -- the number of variables (excluding slacks)
        int n = 4;
    
        // nclin -- the number of linear constraints
        int nclin = 1;
    
        // ncnln -- the number of nonlinear constraints
        int ncnln = 2;
    
        // a[lda*n] -- array of linear constraints, where
        // lda = max(1, nclin). Although the NAG routine e04ucf
        // has A as a two dimensional matrix, for ease of access via the
        // Java Native Interface (JNI) it is much easier to store it as
        // a one dimensional array in Java. We still require the
        // value lda which Fortran will be told is the leading dimension
        // of its 2D array.
        int lda = java.lang.Math.max(1,nclin);
        double[] a;
        a = new double[lda*n];
        // a[i+j*lda] references array element a[i,j] in Fortran order.
        a[0+0*lda] = 1.0;
        a[0+1*lda] = 1.0;
        a[0+2*lda] = 1.0;
        a[0+3*lda] = 1.0;
    
        // bl[n+nclin+ncnln] -- lower bounds for all the variables and general constraints
        double[] bl = {1.0, 1.0, 1.0, 1.0, -1.0e+25, -1.0e+25, 25.0};
    
        // bu[n+nclin+ncnln] -- upper bounds for all the variables and general constraints
        double[] bu = {5.0, 5.0, 5.0, 5.0, 20.0, 40.0, 1.0e+25};
    
        // x[n] -- initial estimate of the solution
        double[] x = {1.0, 5.0, 5.0, 1.0};
    
        // objf[1] -- an array of length 1 to hold the final objective
        // function value computed by e04ucf
        double[] objf = new double[1];
    
        // objgrd[n] -- an array to hold the gradient of the objectve function,
        // computed by e04ucf
        double[] objgrd = new double[n];
    
        // ifail -- output error variable.
        int ifail;
    
        int i;
    
        // Set some options for e04ucf
        e04uef("Nolist");          // Turn off echoing of options by e04uef
        e04uef("Print Level = 0"); // Turn off e04ucf internal monitoring information
    
        System.out.println(" Running e04ucf example program from Java");
        System.out.println(" ----------------------------------------");
        System.out.println(" Problem:");
        System.out.println("");
        System.out.println("   Minimize F(x) = x0*x3*(x0 + x1 + x2) + x2");
        System.out.println("   Subject to bounds");
        for (i = 0; i < n; i++)
          System.out.println("     " + bl[i] + " <= x" + i + " <= " + bu[i]);
        System.out.println("   General linear constraint");
        System.out.println("     x0 + x1 + x2 + x3 <= 20");
        System.out.println("   Nonlinear constraints");
        System.out.println("     x0^2 + x1^2 + x2^2 + x3^2 <= 40");
        System.out.println("     x0*x1*x2*x3 >= 25");
    
        // Call the NAG Library routine e04ucf via the Java Native Interface
        ifail = e04ucf(n, nclin, ncnln, a, lda, bl, bu, confunction, objfunction,
                       objf, objgrd, x);
    
        // Output some results
        System.out.println("");
        System.out.println(" Results returned by NAG nonlinear minimization routine e04ucf");
        System.out.println(" -------------------------------------------------------------");
        System.out.println(" Fail code ifail = " + ifail);
        if (ifail == 0)
          {
            System.out.println(" Final objective function value = " + objf[0]);
            System.out.println(" Solution vector x:");
            for (i = 0; i < n; i++)
              System.out.println("   x[" + i + "] = " + x[i]);
          }
      }
    
    }
    
    Some points to note about this program:

We can compile our Java program with the following command:

  % javac Minimization.java

  • Generating a header file for use by C
  • Having compiled Minimization.java, we can use javah to create a C header file:

      % javah -jni Minimization
    
    The generated header file, Minimization.h, contains these two function prototypes for ther two JNI functions:
      JNIEXPORT jint JNICALL Java_Minimization_e04ucf
        (JNIEnv *, jobject, jint, jint, jint, jdoubleArray, jint, jdoubleArray,
        jdoubleArray, jstring, jstring, jdoubleArray, jdoubleArray, jdoubleArray);
    
      JNIEXPORT void JNICALL Java_Minimization_e04uef
        (JNIEnv *, jobject, jstring);
    

  • Implementing the native function in C code
  • Now that we have created the header file Minimization.h, we can write our C code implementation of Java_Minimization_e04ucf.

    5.1 Source code for the C interface library

    Here is the C source code, from file MinimizationImp.c:

    #include "Minimization.h"
    #include <math.h>
    #include <stdio.h>
    #include <stdlib.h>
    #include <string.h>
    #include <jni.h>
    
    /* Nasty global variables to store pointers to Java environment
       and methods so that we can use them in different parts of this
       C code. */
    JNIEnv *globalJavaEnv;
    jobject globaljavaobject;
    jmethodID globalConfunID;
    jmethodID globalObjfunID;
    
    /* This routine has the interface required by NAG routine e04ucf for
       argument confun. It makes calls back to the Java version of confun */
    void confunFun(int *mode, int *ncnln, int *n,
                   int *ldcj, int needc[], double x[], double c[],
                   double cjac[], int *nstate, int iuser[], double user[])
    {
      int i;
      int *jneedcpt;
      double *jxpt, *jcpt, *jcjacpt;
    
      /* First create some Java arrays to pass to confun */
      jintArray jneedc = (*globalJavaEnv)->NewIntArray(globalJavaEnv, *ncnln);
      jdoubleArray jx = (*globalJavaEnv)->NewDoubleArray(globalJavaEnv, *n);
      jdoubleArray jc = (*globalJavaEnv)->NewDoubleArray(globalJavaEnv, *ncnln);
      jdoubleArray jcjac = (*globalJavaEnv)->NewDoubleArray(globalJavaEnv, ((*ldcj)*(*n)));
    
      /* Copy input arguments to Java arrays needc and x */
      jneedcpt = (*globalJavaEnv)->GetIntArrayElements(globalJavaEnv, jneedc, 0);
      jxpt = (*globalJavaEnv)->GetDoubleArrayElements(globalJavaEnv, jx, 0);
      for (i = 0; i < *ncnln; i++)
        jneedcpt[i] = needc[i];
      for (i = 0; i < *n; i++)
        jxpt[i] = x[i];
      /* Release array elements back to Java (this puts the values
         back into Java arrays jneedc and jx) */
      (*globalJavaEnv)->ReleaseIntArrayElements(globalJavaEnv, jneedc, jneedcpt, 0);
      (*globalJavaEnv)->ReleaseDoubleArrayElements(globalJavaEnv, jx, jxpt, 0);
    
      /* Call the Java method via its method ID */
      (*globalJavaEnv)->CallVoidMethod(globalJavaEnv, globaljavaobject,
                                       globalConfunID, (jint)(*mode),
                                       (jint) (*ncnln), (jint)(*n),
                                       (jint)(*ldcj), jneedc, jx, jc,
                                       jcjac, (jint)(*nstate));
    
      /* Copy results from Java arrays back to C arrays */
      jcpt = (*globalJavaEnv)->GetDoubleArrayElements(globalJavaEnv, jc, 0);
      jcjacpt = (*globalJavaEnv)->GetDoubleArrayElements(globalJavaEnv, jcjac, 0);
      for (i = 0; i < *ncnln; i++)
        c[i] = jcpt[i];
      for (i = 0; i < (*ldcj)*(*n); i++)
        cjac[i] = jcjacpt[i];
    
      /* Release array elements back to Java to free memory */
      (*globalJavaEnv)->ReleaseDoubleArrayElements(globalJavaEnv, jc, jcpt, 0);
      (*globalJavaEnv)->ReleaseDoubleArrayElements(globalJavaEnv, jcjac, jcjacpt, 0);
    }
    
    /* This routine has the interface required by NAG routine e04ucf for
       argument objfun. It makes calls back to the Java version of objfun */
    void objfunFun(int *mode, int *n, double x[], double *objf, double objgrd[],
                   int *nstate, int iuser[], double user[])
    {
      int i;
      double *jobjgrdpt;
    
      /* First create some Java arrays to pass to objfun */
      jdoubleArray jx = (*globalJavaEnv)->NewDoubleArray(globalJavaEnv, *n);
      jdoubleArray jobjgrd = (*globalJavaEnv)->NewDoubleArray(globalJavaEnv, *n);
    
      /* Copy input array x to Java array jx */
      double *jxpt = (*globalJavaEnv)->GetDoubleArrayElements(globalJavaEnv, jx, 0);
      for (i = 0; i < *n; i++)
        jxpt[i] = x[i];
      /* Release array elements back to Java (this puts the values into jx) */
      (*globalJavaEnv)->ReleaseDoubleArrayElements(globalJavaEnv, jx, jxpt, 0);
    
      /* Call Java objfun which fills in array objgrd and returns objf */
      *objf = (*globalJavaEnv)->CallDoubleMethod(globalJavaEnv, globaljavaobject,
                                                 globalObjfunID, (jint) (*mode),
                                                 (jint) (*n), jx, jobjgrd,
                                                 (jint) (*nstate));
    
      /* Get results back from Java to C array objgrd */
      jobjgrdpt = (*globalJavaEnv)->GetDoubleArrayElements(globalJavaEnv,
                                                           jobjgrd, 0);
      for (i = 0; i < *n; i++)
        objgrd[i] = jobjgrdpt[i];
    
      /* Release array elements back to Java to free memory */
      (*globalJavaEnv)->ReleaseDoubleArrayElements(globalJavaEnv, jobjgrd,
                                                   jobjgrdpt, 0);
    }
    
    JNIEXPORT jint JNICALL Java_Minimization_e04ucf
    (
     JNIEnv *env,
     jobject object,
     jint n,
     jint nclin,
     jint ncnln,
     jdoubleArray a,
     jint lda,
     jdoubleArray bl,
     jdoubleArray bu,
     jstring confun,
     jstring objfun,
     jdoubleArray objf,
     jdoubleArray objgrd,
     jdoubleArray x
     )
    {
      /* Local variables and arrays */
      int ldcj;
      int ldr;
      int iter;
      int *istate;
      double *c;
      double *cjac;
      double *clamda;
      double *r;
      int liwork;
      int *iwork;
      int lwork;
      double *work;
      int ifail;
    
      /* N.B. we choose not to use iuser and user arrays in our evaluation
         functions, so these are empty arrays. */
      int iuser[1];
      double user[1];
    
      jclass cls;
      const char *confunction;
      const char *objfunction;
    
      jdouble *a_pt, *bl_pt, *bu_pt, *x_pt, *objf_pt, *objgrd_pt;
    
      /* Array leading dimension information required by the Fortran routine */
      if (ncnln > 0)
        ldcj = ncnln;
      else
        ldcj = 1;
      ldr = n;
    
      /* Compute the amount of workspace we need to supply to e04ucf */
      liwork = 3 * n + nclin + 2 * ncnln;
      if (ncnln == 0 && nclin == 0)
        lwork = 20 * n;
      else if (ncnln == 0 && nclin > 0)
        lwork = 2 * n*n + 20 * n + 11 * nclin;
      else
        lwork = 2 * n*n + n*nclin + 2*n*ncnln + 20*n + 11*nclin + 21*ncnln;
    
      /* Allocate arrays of appropriate size. */
      /* Note that we store cjac as a one dimensional array rather than
         a 2D array as in Fortran, for convenience in communication with
         Java. */
      istate = (int *)malloc((n+nclin+ncnln)*sizeof(int));
      c = (double *)malloc((ncnln)*sizeof(double));
      cjac = (double *)malloc((ldcj*n)*sizeof(double));
      clamda = (double *)malloc((n+nclin+ncnln)*sizeof(double));
      r = (double *)malloc((ldr*n)*sizeof(double));
      iwork = (int *)malloc((liwork)*sizeof(int));
      work = (double *)malloc((lwork)*sizeof(double));
    
      /* Copy the Java env pointers to global space
         so that confunFun and objfunFun can access them. */
      globalJavaEnv = env;
      globaljavaobject = object;
    
      /* Get hold of the name of the user's Java evaluation functions. */
      confunction = (*env)->GetStringUTFChars(env, confun, 0);
      objfunction = (*env)->GetStringUTFChars(env, objfun, 0);
    
      /* Now we have the Java evaluation function names we can use
         them to get hold of handles (method IDs) to the functions.
         Once more, the method IDs are stored globally so that confunFun
         and objfunFun can use them. Note that the Java function signatures
         must be correct. You can find out the signatures after compiling
         the Java program Minimization.java by using the command
           % javap -private -s Minimization
       */
      cls = (*env)->GetObjectClass(env, object);
      globalConfunID = (*env)->GetMethodID(env, cls, confunction, "(IIII[I[D[D[DI)V");
      globalObjfunID = (*env)->GetMethodID(env, cls, objfunction, "(II[D[DI)D");
    
      /* Free up the Java string argument so we don't leak memory. */
      (*env)->ReleaseStringUTFChars(env, confun, confunction);
      (*env)->ReleaseStringUTFChars(env, objfun, objfunction);
    
      if (globalConfunID == 0)
        {
          printf("Cannot find confun method \"%s\" with signature \"(IIII[I[D[D[DI)V\"\n",
                 confunction);
          return -1;
        }
      if (globalObjfunID == 0)
        {
          printf("Cannot find objfun method \"%s\" with signature \"(II[D[DI)D\"\n",
                 objfunction);
          return -1;
        }
    
      /* Extract the arrays from Java */
      a_pt = (*env)->GetDoubleArrayElements(env, a, 0);
      bl_pt = (*env)->GetDoubleArrayElements(env, bl, 0);
      bu_pt = (*env)->GetDoubleArrayElements(env, bu, 0);
      objf_pt = (*env)->GetDoubleArrayElements(env, objf, 0);
      objgrd_pt = (*env)->GetDoubleArrayElements(env, objgrd, 0);
      x_pt = (*env)->GetDoubleArrayElements(env, x, 0);
    
      /* Call to main NAG Library routine e04ucf */
      ifail = -1;
    #ifdef WINDOWS
        E04UCF
    #else
        e04ucf_
    #endif
        (&n, &nclin, &ncnln, &lda, &ldcj, &ldr, a_pt, bl_pt, bu_pt,
         confunFun, objfunFun, &iter, istate, c, cjac, clamda,
         objf_pt, objgrd_pt, r, x_pt, iwork, &liwork, work, &lwork,
         iuser, user, &ifail);
    
      /* Release the array elements back to Java and free memory. */
      (*env)->ReleaseDoubleArrayElements(env, a, a_pt, 0);
      (*env)->ReleaseDoubleArrayElements(env, bl, bl_pt, 0);
      (*env)->ReleaseDoubleArrayElements(env, bu, bu_pt, 0);
      (*env)->ReleaseDoubleArrayElements(env, objf, objf_pt, 0);
      (*env)->ReleaseDoubleArrayElements(env, objgrd, objgrd_pt, 0);
      (*env)->ReleaseDoubleArrayElements(env, x, x_pt, 0);
    
      return ifail;
    
    }
    
    // Interface to option setting routine e04uef
    JNIEXPORT void JNICALL Java_Minimization_e04uef
      (JNIEnv *env, jobject object, jstring option)
    {
      const char *coption;
    
      /* Get hold of the Java option string. */
      coption = (*env)->GetStringUTFChars(env, option, 0);
    
      /* Call the option setting routine */
    #ifdef WINDOWS
      E04UEF(coption, strlen(coption));
    #else
      e04uef_(coption, strlen(coption));
    #endif
    
      /* Free up the Java string argument so we don't leak memory. */
      (*env)->ReleaseStringUTFChars(env, option, coption);
    }
    

    5.2 Description of the C code

    The function named Java_Minimization_e04ucf is our C implementation of the Java-declared method e04ucf.

    We cannot pass the Java methods objfun and confun which evaluate the objective function and nonlinear constraints directly to the NAG Fortran Library routine E04UCF, so we need to wrap them in C functions. These C functions we name objfunFun and confunFun respectively:

       void objfunFun(int *mode, int *n, double x[], double *objf, double objgrd[],
                      int *nstate, int iuser[], double user[]);
    
       void confunFun(int *mode, int *ncnln, int *n,
                      int *ldcj, int needc[], double x[], double c[],
                      double cjac[], int *nstate, int iuser[], double user[]);
    
    These functions have the argument types and return types required by the NAG Library routine E04UCF. Notice in particular that all scalar arguments (such as mode) are passed as pointers, as required by Fortran. Further, notice that argument cjac is declared as a one-dimensional array rather than the 2-D array specified by the Fortran routine. We do this because 2-D Java arrays do not map easily onto Fortran 2-D arrays.

    Inside objfunFun and confunFun we do nothing but call the equivalent Java methods. Once again, the trick is in knowing how to make these calls to Java. We do this using the JNI functions CallDoubleMethod for objfun (because in Java we defined that method to have return type double) and CallVoidMethod for confun.

    The Java methods objfun and confun both need to be passed array arguments, elements of which the methods need to fill in. As with Example 3, we obtain the method ID of the two methods and store them in global variables in our C code so that they can be accessed from inside the C evaluation functions as well as the main JNI function. These IDs are obtained by passing the appropriate names and signatures to calls of JNI function GetMethodID. Note that a good way to discover Java method signatures is to use the command

      % javap -private -s Minimization
    
    after compiling the Java program Minimization.java.

    A complication is that our C functions objfunFun and confunFun need to pass Java arrays to the Java methods, but themselves have only C (or Fortran) style arrays. Therefore the C code needs to create Java arrays, then copy the contents of the C arrays to the Java arrays. To create a Java double array, we use the JNI function NewDoubleArray followed by a call of GetDoubleArrayElements to get a C pointer to the array elements, then ReleaseDoubleArrayElements to put the C contents into the Java array, before calling the Java method. On return from the Java method we again call GetDoubleArrayElements to obtain the results. For integer arrays, JNI functions NewIntArray, GetIntArrayElements and ReleaseIntArrayElements are appropriate. It is very important to get the order of calls to these JNI functions exactly right to ensure that data is in the right place at the right time. Follow the example in MinimizationImp.c.

    5.3 Returning results to Java

    In this example all results are returned to Java via the array arguments which came from the Java call to the native method - apart from the error code IFAIL which is returned via the function name. Notice that we declared the argument objf, which contains the optimal function value, as an array of length 1. This is because Java methods cannot update the contents of scalar arguments, but can update array contents.

  • Building the shareable library or DLL
  • This step is operating-system dependent.

    The compiler flags used were described in Section 7 of Example 1. Note that when building under Microsoft Windows we also add the C compiler switch -DWINDOWS to tell the C code that the name of the Fortran Library routines E04UCF and E04UEF must be given in upper case, as required by Windows versions of the NAG Fortran Library. For UNIX machines the name will typically be in lower case with an appended underscore, when called from C, i.e. "e04ucf_".

  • Running the program
  • Assuming that all has gone well, we can run the program using the command

      % java Minimization
    
    The expected output looks like this:

      Running e04ucf example program from Java
      ----------------------------------------
      Problem:
    
        Minimize F(x) = x0*x3*(x0 + x1 + x2) + x2
        Subject to bounds
          1.0 <= x0 <= 5.0
          1.0 <= x1 <= 5.0
          1.0 <= x2 <= 5.0
          1.0 <= x3 <= 5.0
        General linear constraint
          x0 + x1 + x2 + x3 <= 20
        Nonlinear constraints
          x0^2 + x1^2 + x2^2 + x3^2 <= 40
          x0*x1*x2*x3 >= 25
    
      Results returned by NAG nonlinear minimization routine e04ucf
      -------------------------------------------------------------
      Fail code ifail = 0
      Final objective function value = 17.014017289134703
      Solution vector x:
        x[0] = 1.0
        x[1] = 4.742999642848296
        x[2] = 3.821149976895378
        x[3] = 1.379408294178579
    

    (If you get an error message saying that a library cannot be located, see the tip given in Example 1). For this example, since we linked to the NAG FOrtran Library rather than the NAG C Library, we need to set the relevant environment variable to point to the directory containing the NAG Fortran Library.

  • Quick summary of how to build the minimization example
  • Given the two source files Minimization.java and MinimizationImp.c, issue the following commands:
    Start of report   Skip to Example 1   Skip to Example 2   Skip to Example 3   Skip to Example 5   Skip to Example 6
    Copyright 2009 Numerical Algorithms Group
    Page last updated 2011-01-25 annak
    [NP3671]