External Calling
(Using Compiled Code in Maple)
|
Introduction
|
|
External calling is a feature in Maple that allows you to seamlessly integrate your compiled C, Fortran, or Java code into Maple. Functions written in these languages can be linked and used as if they were native Maple procedures. There are many clear benefits to using this feature.
Reuse Code.
Many users have invested a lot of time writing fast algorithms in other languages. These algorithms do not need to be rewritten. They can usually be called as is from Maple.
Extend the Power of Maple.
Numerous routines for doing a variety of mathematical and non-math operations are available on the Web. You can use these routines in Maple to improve the range of capabilities of Maple, or write your own. You can add a database object, link to controlled hardware via a serial port, or interface with another program. The possibilities are endless.
Improve Performance.
It is hard to beat the speed of compiled C and Fortran code. You can use external calling to access fast numerical libraries. In some cases, Maple already does this. Access to the NAG library routines and other numerical algorithms is already built into Maple using the external calling mechanism.
|
|
Calling an External Function
|
|
|
Simple C Example
|
|
Consider the following C function.
int __stdcall mult( int a, int b )
{
return a * b;
}
This function was saved in a file, mult.c, and compiled into the shared library, mult.dll.
To run this example, mult.dll must be created and placed in the c:/shared directory.
mult.dll can be built using the following Microsoft Visual C/C++ compiler command line.
C:\shared> cl -Gz mult.c -link -dll -export:mult -out:mult.dll
Alternately, it can be built using the Watcom compiler that comes with Maple with a compile and link command as follows. Note that the include path here isn't necessary, but is needed with more complicated examples that use the external API.
C:\shared> wcl386 -c -bd -DM_DECL=__stdcall -DM_CDECL=__cdecl -i="c:\watcom-1.3\h;c:\watcom-1.3\h\nt;c:\watcom-1.3\maple11.00\include" mult.c -fo=mult.obj
C:\shared> wlink system nt_dll libpath 'c:\Program Files\Maple 11\bin.win' name mult.dll file mult.obj library 'maplec.lib' export '_mult@8'
A similar function (without the __stdcall declaration) can be compiled on MacOS, Linux, or any of the Unix platforms.
% gcc-shared mult.c -o libmult.so
To access this function from Maple, provide a description of the external function, including:
Name of the function in the shared library (mult)
Name of the shared library (mult.dll)
Parameter and return types (all ints, which are written integer[4] in Maple)
>
|
|
The function, my_mult, can be used like any other maple procedure now.
>
|
|
| (2.1.1) |
|
|
Simple Fortran Example
|
|
Consider the following Fortran function.
INTEGER mult( a, b )
INTEGER a
INTEGER b
mult = a * b
END
It has been compiled into the shared library, mult.dll.
To access this function from Maple, provide a description of the external function, including:
Name of the function in the shared library (mult)
Name of the shared library (mult.dll)
Keyword FORTRAN to indicate which language the target is written in
Parameter and return types (all INTEGERs, which are written integer[4] in Maple)
>
|
|
The function, fmult, can be used like any other maple procedure now.
>
|
|
| (2.2.1) |
|
|
Simple Java Example
|
|
Consider the following Java function.
public class mult {
static int jmult( int a, int b )
{
return a * b;
}
}
It has been saved in the file, mult.java, and compiled into byte code in the class library, mult.class using the command javac mult.java.
To access this function from Maple, provide a description of the external function, including:
Name of the function in the shared library (mult)
Name of the class library (mult.class)
Path to this class and any other associated class files
Keyword JAVA to indicate which language the target is written in
Parameter and return types (all ints, which are written integer[4] in Maple)
>
|
|
The function, jmult, can be used like any other maple procedure now.
>
|
|
| (2.3.1) |
|
|
Type Declarations
|
|
The basic type conversions are as follows.
Maple C Java Fortran
integer[1] char byte BYTE
integer[2] short short INTEGER*2
integer[4] int int INTEGER
integer[8] long long* long INTEGER*8
float[4] float float REAL
float[8] double double DOUBLE PRECISION
*C's long size is compiler dependent, usually 32-bits, but often 64-bits when using a 64-bit compiler.
Compound types can be built out of the basic types.
Maple C
Array(datatype=integer[4],order=..., etc.) int x[size];
string[n] char s[n];
complex[4] struct { float r, i };
complex[8] struct { double r, i };
REF(integer[4]) int *x;
|
|
Matrix Example
|
|
This example uses external C code to multiply two matrices. The source code compiled into mat_mult.dll is as follows.
void mat_mult( double *A, double *B, double *C, int I, int J, int K )
{
int i, j, k;
double t;
for( i = 0; i < I; i++ ) {
for( k = 0; k < K; k++ ) {
t = 0.0;
for( j = 0; j < J; j++ ) {
t += A[i*J+j] * B[j*K+k];
}
C[i*K+k] = t;
}
}
}
The Maple declaration including the name of the external function and library, plus the parameter types is as follows.
>
|
|
The function is now ready to use. Create a pair of C_order Matrices with datatype=float[8] to test your function.
>
|
|
| (2.5.1) |
>
|
|
| (2.5.2) |
Create an empty Matrix of the right size to store the result.
>
|
|
| (2.5.3) |
Call the external function.
>
|
|
>
|
|
| (2.5.4) |
Compare the result with the one the built-in Maple routine gives.
>
|
|
| (2.5.5) |
|
|
Compiler Notes
|
|
|
_stdcall on Windows
|
|
Most UNIX C compilers use one and only one calling convention. Windows C compilers can use any one of three conventions, __cdecl, __fastcall, or __stdcall. Maple only recognizes functions compiled with __stdcall. This is not the default used in most Windows compilers. When using Microsoft Visual C/C++, compile with the flag /Gz, or declare functions __stdcall in the code.
|
|
Functions need to be exported
|
|
Some compilers will export, or make public, all functions in your program. This means programs using the shared library can find all functions in order to make use of them. Other compilers need explicit decorations on the function declaration, and/or a file that lists all the functions that should be exported. The MSVC compiler can export functions in a variety of ways, including listing them in a .def file, putting the __declspec(dllexport) decoration on the function definition, or by using the /export:name option to the compiler.
|
|
MSVC command line
|
|
To build and link the Matrix Example shown above with the MSVC compiler, do one of the following.
Two step method (separate compile and link):
cl -c -Gz mat_mult.c
link /dll /export:mat_mult mat_mult.obj
One step method (compile and link on same line):
cl -Gz mat_mult.c -link -dll -export:mat_mult -out:mat_mult.dll
|
|
gcc command line
|
|
gcc on Linux and other UNIX platforms use __stdcall and export all functions by default. Use the -shared option to create a shared library.
gcc -shared mat_mult.c -olibmat_mult.so
|
|
Fortran Compilers
|
|
Fortran compilers use several different calling conventions. Fortran external call has only been tested with the GNU g77 compiler on various platforms. There are known incompatibilities with other compilers, especially when using string arguments. Try to get a simple example working before moving to a more complicated one.
|
|
Calling Java
|
|
To call Java functions, jvm.dll needs to be in your PATH. You may need to install a Java Developer's Kit (JDK) to get this dll. UNIX installations will need at least libjvm.so and libhpi.so in the PATH. Typically the following directories need to be added to your PATH.
Windows: $JAVA/jre/bin/classic, $JAVA/jre/bin
Solaris: $JAVA/jre/lib/sparc, $JAVA/jre/lib/sparc/native_threads
Linux: $JAVA/jre/lib/i386/classic, $JAVA/jre/lib/i386/native_threads, $JAVA/jre/lib/i386
Note: These paths may differ depending on the version and supplier of the JDK you install. $JAVA is the base directory where the JDK was installed.
|
|
|
|
Advanced Features
|
|
|
Wrapper Generation
|
|
So far, all the examples given leave most of the behind the scenes work up to Maple. Maple's internal representation of some data types is not necessarily compatible with the representation used by other languages. As a result sometimes conversions are needed. Often Maple can automatically do the conversions. For example, in the mult examples above, Maple needs to turn its INTPOS data structure into a 4-byte hardware integer before calling the external functions. In some cases the conversions are non-trivial and cannot be done automatically. As an intermediate step before reading the section on writing custom wrappers, understand that Maple can generate and compile wrappers on the fly. A wrapper is a bit of code that makes use of Maple's external API for type conversions before calling the real external function.
Use the keyword, WRAPPER, to tell Maple to do translations via a compiled wrapper. Consider the same Matrix example presented earlier. This time specify the keyword, WRAPPER.
>
|
|
The ultimate result is the same (provided you have a compiler that is properly set up). The matmult function is defined and can be used like any other Maple function. This time, however, all translations are done by executing the code in the newly created and compiled mwrap_mat_mult.c file.
|
mwrap_mat_mult.c
|
|
/* MWRAP_mat_mult Wrapper
Generated automatically by Maple 8
Do not edit this file. */
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <mplshlib.h>
#include <maplec.h>
MKernelVector mapleKernelVec;
typedef void *MaplePointer;
ALGEB *args;
static ALGEB Array[3];
static FLOAT64 *MWRAPARRAY_A0_ALLOC( ALGEB m ) {
MapleRaiseError(mapleKernelVec,"cannot create ARRAY with unspecified bounds");
}
static ALGEB MWRAPARRAY_A0_TO_MAPLE( FLOAT64 *A ) {
RTableSettings s;
if( Array[0] && (void*)A == RTableDataBlock(mapleKernelVec,Array[0]) )
return(Array[0]);
MapleRaiseError(mapleKernelVec,"cannot create ARRAY with unspecified bounds");
}
static void MWRAPARRAY_A0_TO_C( FLOAT64 **A, ALGEB m ) {
RTableSettings s;
if( Array[0] && (void*)(*A) == RTableDataBlock(mapleKernelVec,Array[0]) )
return;
RTableGetSettings(mapleKernelVec,&s,m);
if( s.data_type != RTABLE_FLOAT64
|| !IsMapleNULL(mapleKernelVec,s.index_functions) ) {
s.data_type = RTABLE_FLOAT64;
Array[0] = RTableCopy(mapleKernelVec,&s,m);
}
else Array[0] = m;
*A = (FLOAT64 *)RTableDataBlock(mapleKernelVec,Array[0]);
}
static FLOAT64 *MWRAPARRAY_A1_ALLOC( ALGEB m ) {
MapleRaiseError(mapleKernelVec,"cannot create ARRAY with unspecified bounds");
}
static ALGEB MWRAPARRAY_A1_TO_MAPLE( FLOAT64 *A ) {
RTableSettings s;
if( Array[1] && (void*)A == RTableDataBlock(mapleKernelVec,Array[1]) )
return(Array[1]);
MapleRaiseError(mapleKernelVec,"cannot create ARRAY with unspecified bounds");
}
static void MWRAPARRAY_A1_TO_C( FLOAT64 **A, ALGEB m ) {
RTableSettings s;
if( Array[1] && (void*)(*A) == RTableDataBlock(mapleKernelVec,Array[1]) )
return;
RTableGetSettings(mapleKernelVec,&s,m);
if( s.data_type != RTABLE_FLOAT64
|| !IsMapleNULL(mapleKernelVec,s.index_functions) ) {
s.data_type = RTABLE_FLOAT64;
Array[1] = RTableCopy(mapleKernelVec,&s,m);
}
else Array[1] = m;
*A = (FLOAT64 *)RTableDataBlock(mapleKernelVec,Array[1]);
}
static FLOAT64 *MWRAPARRAY_A2_ALLOC( ALGEB m ) {
MapleRaiseError(mapleKernelVec,"cannot create ARRAY with unspecified bounds");
}
static ALGEB MWRAPARRAY_A2_TO_MAPLE( FLOAT64 *A ) {
RTableSettings s;
if( Array[2] && (void*)A == RTableDataBlock(mapleKernelVec,Array[2]) )
return(Array[2]);
MapleRaiseError(mapleKernelVec,"cannot create ARRAY with unspecified bounds");
}
static void MWRAPARRAY_A2_TO_C( FLOAT64 **A, ALGEB m ) {
RTableSettings s;
if( Array[2] && (void*)(*A) == RTableDataBlock(mapleKernelVec,Array[2]) )
return;
RTableGetSettings(mapleKernelVec,&s,m);
if( s.data_type != RTABLE_FLOAT64
|| !IsMapleNULL(mapleKernelVec,s.index_functions) ) {
s.data_type = RTABLE_FLOAT64;
Array[2] = RTableCopy(mapleKernelVec,&s,m);
}
else Array[2] = m;
*A = (FLOAT64 *)RTableDataBlock(mapleKernelVec,Array[2]);
}
/* main - MWRAP_mat_mult */
ALGEB MWRAP_mat_mult( MKernelVector kv,
void (*fn) ( FLOAT64 *a1, FLOAT64 *a2, FLOAT64 *a3, INTEGER32 a4, INTEGER32 a5, INTEGER32 a6 ),
ALGEB fn_args )
{
FLOAT64 *a1;
FLOAT64 *a2;
FLOAT64 *a3;
INTEGER32 a4;
INTEGER32 a5;
INTEGER32 a6;
int i;
mapleKernelVec = kv;
args = (ALGEB*) fn_args;
if( MapleNumArgs(mapleKernelVec,(ALGEB)args) != 6 )
MapleRaiseError(mapleKernelVec,"Incorrect number of arguments");
for( i=0; i<3; ++i ) Array[i] = NULL;
/* integer[4] */
a4 = MapleToInteger32(mapleKernelVec,args[4]);
/* integer[4] */
a5 = MapleToInteger32(mapleKernelVec,args[5]);
/* integer[4] */
a6 = MapleToInteger32(mapleKernelVec,args[6]);
/* ARRAY(datatype = float[8]) */
if( IsMaplePointerNULL(mapleKernelVec,args[1]) )
a1 = NULL;
else if( IsMapleUnassignedName(mapleKernelVec,args[1]) )
a1 = MWRAPARRAY_A0_ALLOC( args[1]);
else {
MWRAPARRAY_A0_TO_C(&a1,args[1]);
}
/* ARRAY(datatype = float[8]) */
if( IsMaplePointerNULL(mapleKernelVec,args[2]) )
a2 = NULL;
else if( IsMapleUnassignedName(mapleKernelVec,args[2]) )
a2 = MWRAPARRAY_A1_ALLOC( args[2]);
else {
MWRAPARRAY_A1_TO_C(&a2,args[2]);
}
/* ARRAY(datatype = float[8]) */
if( IsMaplePointerNULL(mapleKernelVec,args[3]) )
a3 = NULL;
else if( IsMapleUnassignedName(mapleKernelVec,args[3]) )
a3 = MWRAPARRAY_A2_ALLOC( args[3]);
else {
MWRAPARRAY_A2_TO_C(&a3,args[3]);
}
(*fn)(a1, a2, a3, a4, a5, a6);
return( ToMapleNULL(mapleKernelVec) );
}
|
The generated wrapper has some ability to convert Arrays with the wrong data type, and handle one level of procedure parameters (callbacks). Most importantly, the generated wrapper gives a good template example of how to write your own wrapper (after cleaning up the extra auto-generated layers in the file).
|
Fortran External Call with Wrapper
|
|
|
Introduction
|
|
The Fortran external API (given in libmaplefortran.so) provides routines that translate Maple data types to Fortran data types and vice versa. This allows you to invoke Fortran routines from the Maple environment. To create a Fortran external call from Maple:
1. Create a Fortran wrapper library.
In the wrapper library, translate all the Maple input parameters into Fortran data types using the routines provided by the Fortran external API. Since Fortran 77 does not provide mechanism for dynamic memory allocation, a static array is to be declared for each array parameter. This static array acts as a container and must be large enough to cover all the input array sizes. Every array parameter has to be copied to its corresponding container. After the translation is done and all the parameters have been converted to native Fortran data types, the core computation can proceed. When finished, the result must be converted back to a Maple data type.
2. Compile the wrapper library.
The wrapper must be linked with maplefortran.dll, which provides the Fortran external API routines.
3. Declare the function in Maple.
Use define_external to create a link to the Fortran wrapper shared library. This link allows invocation of the Fortran routine from Maple.
|
|
FFT Example with Wrapper
|
|
The following code invokes a Fortran external call to fftwrap that computes the Fast Fourier Transform (FFT) of a rectangular real signal. FFT is the fundamental operation of signal processing applications. FFT transforms a signal from the time domain to its spectrum in the frequency domain. The Fortran routines to compute the FFT are taken from SLATEC Common Mathematical Library, Version 4.1. The source code can be downloaded from the website http://netlib2.cs.utk.edu/slatec
|
fftwrap.f (source of SLATEC fftwrap function)
|
|
C***Fortran wrapper fftwrap for the FFT routine
C***Input must be a real signal
C***Every wrapper routine accepts two parameters
C*** kv - address to the Maple kernel vector
C*** args - address to the array of input arguments
INTEGER FUNCTION fftwrap(kv, args)
INTEGER kv
INTEGER args
C Header file for the Fortran external API
INCLUDE "maplefortran.hf"
PARAMETER( MAXFFT = 4096 )
CHARACTER ERRMSG*100
INTEGER ERRMSG_LEN
REAL W( 2*MAXFFT + 15 )
C***Array used as "container" to the input array parameters
C X - real input signal
C RY - real part of FFT(X)
C IY - imaginary part of FFT(X)
C N - length of input signal
REAL X( MAXFFT )
DOUBLEPRECISION RY( MAXFFT )
DOUBLEPRECISION IY( MAXFFT )
C length of input signal
INTEGER N
C SX, SRY, SIY, and SN are addresses of Maple DAGS that correspond to
C X, RY, IY, and N, respectively
INTEGER SX, SRY, SIY, SN
C Bounds to Maple rtables
INTEGER RBOUNDS(2)
C Bounds to Fortran arrays
INTEGER FBOUNDS(2)
ERRMSG_LEN = 100
IF ( num_args(kv, args) .NE. 4 ) THEN
ERRMSG = 'Four parameters expected - X, RY, IY, N'
CALL error( kv, ERRMSG, ERRMSG_LEN )
fftwrap = 0
RETURN
ENDIF
C Extract each argument from the array args
SX = extract_arg( kv, args, 1 )
SRY = extract_arg( kv, args, 2 )
SIY = extract_arg( kv, args, 3 )
SN = extract_arg( kv, args, 4 )
C Argument checking
IF ( is_rtable(kv, SX) .EQ. 0 .OR.
* rtable_is_real(kv, SX) .EQ. 0 ) THEN
ERRMSG = 'input X has to be a real vector'
CALL error( kv, ERRMSG, ERRMSG_LEN )
fftwrap = 0
RETURN
ENDIF
IF ( is_rtable(kv, SRY) .EQ. 0 .OR.
* rtable_is_real(kv, SRY) .EQ. 0 ) THEN
ERRMSG = 'output RY has to be a real vector'
CALL error( kv, ERRMSG, ERRMSG_LEN )
fftwrap = 0
RETURN
ENDIF
IF ( is_rtable(kv, SIY) .EQ. 0 .OR.
* rtable_is_real(kv, SIY) .EQ. 0 ) THEN
ERRMSG = 'output IY has to be a real vector'
CALL error( kv, ERRMSG, ERRMSG_LEN )
fftwrap = 0
RETURN
ENDIF
IF ( is_integer(kv, SN) .EQ. 0 ) THEN
ERRMSG = 'input N has to be an integer'
CALL error( kv, ERRMSG, ERRMSG_LEN )
fftwrap = 0
RETURN
ENDIF
C Translate Maple data type to Fortran data type
N = maple_to_integer32( kv, SN )
RBOUNDS(1) = 1
RBOUNDS(2) = N
FBOUNDS(1) = 1
FBOUNDS(2) = MAXFFT
CALL copy_to_array( kv,
* SX,
* X,
* 1,
* RBOUNDS,
* 1,
* FBOUNDS,
* IRTABLE_FLOAT32
* )
C Invoke the FFT routine to find the real and imaginary part of FFT(X)
CALL realfft( X, RY, IY, N )
C Translate Fortran data type back to Maple data type
RBOUNDS(1) = 1
RBOUNDS(2) = N
FBOUNDS(1) = 1
FBOUNDS(2) = MAXFFT
CALL copy_to_rtable( kv,
* RY,
* SRY,
* 1,
* FBOUNDS,
* 1,
* RBOUNDS,
* IRTABLE_FLOAT64
* )
RBOUNDS(1) = 1
RBOUNDS(2) = N
FBOUNDS(1) = 1
FBOUNDS(2) = MAXFFT
CALL copy_to_rtable( kv,
* IY,
* SIY,
* 1,
* FBOUNDS,
* 1,
* RBOUNDS,
* IRTABLE_FLOAT64
* )
fftwrap = 0
END
|
|
realfft.f (source of SLATEC realfft function)
|
|
C***FFT routine - it is the front end to the routines of the SLATEC library
C X - real input signal
C RY - real part of FFT(X)
C IY - imaginary part of FFT(X)
C N - length of input signal
SUBROUTINE realfft( X, RY, IY, N )
REAL X(*)
DOUBLEPRECISION RY(*)
DOUBLEPRECISION IY(*)
INTEGER N
PARAMETER( MAXFFT = 4096 )
REAL W( 2*MAXFFT + 15 )
C Call the SLATEC routines to compute the FFT of the real signal X
CALL RFFTI( N, W )
CALL RFFTF( N, X, W )
C Since the input is a real signal, the FFT routine only computes the first
C half of the spectrum. The other half can be obtained by taking the
C conjugate of the first half.
C Real part of the FFT
RY(1) = X(1)
DO 10, I = 1, N/2-1
RY(I+1) = X(2*I)
RY(N-I+1) = X(2*I)
10 CONTINUE
RY(N/2+1) = X(N)
C Imaginary part of the FFT
IY(1) = 0
DO 20, I = 1, N/2-1
IY(I+1) = X(2*I+1)
IY(N-I+1) = -X(2*I+1)
20 CONTINUE
IY(N/2+1) = 0
END
|
Setup the parameters.
>
|
|
Build rectangular input signal.
>
|
|
Centralize and plot input signal.
>
|
|
Real and Imaginary part of FFT(X)
>
|
|
EXTERNAL CALL (with wrapper) TO FIND FFT OF X
>
|
|
Centralize and plot the magnitude of FFT(X).
>
|
|
|
|
FFT Example without Wrapper
|
|
|
Introduction
|
|
Alternatively, Maple also allows you to invoke Fortran routines directly without using any wrapper program. In this case, Maple makes all the data type translation automatically and hence eliminates the need for a wrapper program. This method is cleaner and requires less overhead than the wrapper method because array input does not have to be copied to a container.
|
Again, using the FFT example, you can make an external call to the FFT routine realfft directly without the need of a wrapper program. All the data type translation is done automatically by Maple.
The following Maple code directly invokes the FFT routine realfft to compute the FFT of a rectangular signal.
Setup
>
|
|
Build rectangular input signal.
>
|
|
Centralize and plot input signal.
>
|
|
Real and Imaginary part of FFT(X)
>
|
|
EXTERNAL CALL (without wrapper) TO FIND FFT OF X
>
|
|
Centralize and plot the magnitude of FFT(X).
>
|
|
|
|
|
|
Writing Your Own Wrapper
|
|
For complete control over data conversions, Maple allows modification of existing wrappers and creation of custom wrappers. There are numerous C and Fortran functions available for translating and manipulating data structures, using define_external.
Custom wrappers have the following advantages over using the default auto-conversion method.
Efficiency
Because you have complete control in the wrapper, you can select which conversions need to be done right away, and which can be deferred to later (if at all). Sometimes parameters need to be converted on the way out, but not on the way in and vice versa. Avoiding extra conversions can reduce the overhead of calling an external function, especially when working with large objects or using an external function that gets called hundreds of times.
Access to More Objects
Without writing your own wrapper, you cannot manipulate structures such as lists, tables, or Arrays of polynomials. Custom wrappers give you access to almost all of Maple's standard data structures.
Direct Access to Maple Objects
Maple objects can be created and manipulated in-place without ever needing to convert to a native hardware type.
Ability to Call Back into Maple
Maple procedures can be evaluated in external code. Thus symbolic manipulations can be done on intermediate results all in the external code.
Note: The biggest disadvantage is added complexity.
When declaring a wrapper external function in Maple via define_external, there is no need to specify the parameter types. The wrapper always gets the same two arguments -- a kernel handle and an expression sequence (dag) of arguments. The keyword, MAPLE, is used to denote a custom wrapper function.
|
Wrapper Entry Point
|
|
The main entry point in your custom wrapper should accept two parameters, the MKernelVector handle and something of type ALGEB. When Maple calls this main entry point, it will pass the kernel handle and an EXPSEQ dag. Index directly into the EXPSEQ dag to get at the raw Maple arguments. Use the NumArgs function available in maplec.dll to determine the number of arguments. Note that the generated wrapper uses a slightly different argument list than the custom wrapper. The custom wrapper does not get a pointer to the main external function.
|
Generated wrapper main entry point prototype
|
|
ALGEB MWRAP_mult(
MKernelVector kv,
INTEGER32 (*fn) ( INTEGER32 a1, INTEGER32 a2 ),
ALGEB fn_args
);
|
|
Custom wrapper main entry point prototype
|
|
ALGEB MWRAP_mult(
MKernelVector kv,
ALGEB fn_args
);
|
|
|
Custom Matrix Example Wrapper
|
|
The following is the code of the wrapper routine matmultwrap for matrix multiplication.
|
matmult wrapper
|
|
#include <stdio.h>
#include "mplshlib.h"
#include "maplec.h"
void matmult( double *A, double *B, double *C, int I, int J, int K );
void* matmultwrap( MKernelVector kv, ALGEB args )
{
RTableSettings s;
M_INT bounds[4];
double *A, *B, *C;
int I, J, K;
/* Check number of arguments is correct */
if( NumArgs(kv, args) != 6 )
Error( kv, "6 arguments required for matmultwrap" );
/* Argument checking */
if( !IsRTable(kv, (ALGEB) args[1]) ||
!RTableIsReal(kv, (ALGEB) args[1]) )
Error( kv, "A must be a FLOAT64 matrix" );
if( !IsRTable(kv, (ALGEB) args[2]) ||
!RTableIsReal(kv, (ALGEB) args[2]) )
Error( kv, "B must be a FLOAT64 matrix" );
if( !IsRTable(kv, (ALGEB) args[3]) ||
!RTableIsReal(kv, (ALGEB) args[3]) )
Error( kv, "C must be a FLOAT64 matrix" );
if( !IsInteger(kv, (ALGEB) args[4]) )
Error( kv, "I must be an integer" );
if( !IsInteger(kv, (ALGEB) args[5]) )
Error( kv, "J must be an integer" );
if( !IsInteger(kv, (ALGEB) args[6]) )
Error( kv, "K must be an integer" );
A = (FLOAT64 *) RTableDataReturn( kv, (ALGEB) args[1] );
B = (FLOAT64 *) RTableDataReturn( kv, (ALGEB) args[2] );
C = (FLOAT64 *) RTableDataReturn( kv, (ALGEB) args[3] );
I = MapleToInteger32( kv, (ALGEB) args[4] );
J = MapleToInteger32( kv, (ALGEB) args[5] );
K = MapleToInteger32( kv, (ALGEB) args[6] );
bounds[0] = 1;
bounds[1] = I;
bounds[2] = 1;
bounds[3] = K;
RTableGetDefaults( kv, &s );
s.num_dimensions = 2;
s.data_type = 6;
s.storage = 4;
s.order = 1;
matmult( A, B, C, I, J, K );
return NULL;
}
|
>
|
|
The Maple definition of this function would look like the following.
|
|
C Example with Callbacks
|
|
This example uses a C implementation of Newton's method of finding roots. The equation evaluation is done in Maple and the root finding loop is implemented in C.
Solve the following equation.
>
|
|
| (3.2.3.1) |
|
newton.c
|
|
The external routine requires the above function, its derivative, a tolerance, and an initial guess. The code is as follows:
/* newton.c*/
#include <stdio.h>
#include <math.h>
/* An implementation of Newton's method for finding a root
of the given function, f. The derivative of f, fprime
must also be provided.
The initial guess is continuously improved until
the absolute value of f(guess) is less than or equal to
the given tolerance. The improved guess is returned.
*/
double newton(
double (*f) (double),
double (*fprime) (double),
double guess,
double tolerance )
{
while (fabs(f(guess)) > tolerance ) {
guess = guess - f(guess)/fprime(guess);
}
return(guess);
}
On Solaris, this was compiled into a shared library using the command,
cc -G newton.c -o libnewton.so
On Windows, this was compiled using MSVC with the command,
cl newton.c -Gz -Fenewton.dll -link -dll -export:newton
|
|
Method 1: Generated wrapper
|
|
The external function must be described in Maple so the given parameters can be properly translated to C hardware format. The arguments needed by define_external are:
1. Name of the external function
2. Parameter types passed to the external function including the return value
3. Keyword 'WRAPPER' to indicate a wrapper file should be generated
4. Name of the shared library (dll) containing the external function
>
|
|
|
|
Method 2: Custom Wrapper
|
|
When customizing your wrapper it is often a good idea to start with the generated wrapper.
|
mwrap_newton.c
|
|
/* MWRAP_newton Wrapper */
/* this is needed in order to use Maple's external API */
#include <maplec.h>
/* prototype for external function */
extern double newton(
double (*f) (double),
double (*fprime) (double),
double guess,
double tolerance );
/* some global variables need to be used in the callbacks */
MKernelVector mapleKernelVec;
static ALGEB Proc[2];
/* callback - f */
static FLOAT64 wrap_f( FLOAT64 x )
{
ALGEB arg1, mr;
/* convert float[8] from C to Maple */
arg1 = ToMapleFloat(mapleKernelVec,(double) x);
/* evaluate Proc[0] */
mr = EvalMapleProc(mapleKernelVec,Proc[0], 1, arg1);
/* convert the result back to C */
return( MapleToFloat64(mapleKernelVec,mr) );
}
/* callback - fprime */
static FLOAT64 wrap_fprime( FLOAT64 x )
{
ALGEB arg1, mr;
/* convert float[8] from C to Maple */
arg1 = ToMapleFloat(mapleKernelVec,(double) x);
/* evaluate Proc[1] */
mr = EvalMapleProc(mapleKernelVec,Proc[1], 1, arg1);
/* convert the result back to C */
return( MapleToFloat64(mapleKernelVec,mr) );
}
/* main - MWRAP_newton */
/* When writing your own wrapper there is no function pointer passed here */
ALGEB MWRAP_newton( MKernelVector kv, ALGEB fn_args )
{
FLOAT64 guess, tolerance;
FLOAT64 r;
ALGEB *args;
mapleKernelVec = kv;
args = (ALGEB*) fn_args;
/* argument checking and conversion */
if( MapleNumArgs(mapleKernelVec,(ALGEB)args) != 4 )
MapleRaiseError(mapleKernelVec,"Incorrect number of arguments");
/* PROC(x::float[8],RETURN::float[8]) */
if( !IsMapleProcedure(mapleKernelVec,args[1]) )
MapleRaiseError(mapleKernelVec,"procedure expected");
Proc[0] = args[1];
/* PROC(x::float[8],RETURN::float[8]) */
if( !IsMapleProcedure(mapleKernelVec,args[2]) )
MapleRaiseError(mapleKernelVec,"procedure expected");
Proc[1] = args[2];
/* float[8] */
guess = MapleToFloat64(mapleKernelVec,args[3]);
/* float[8] */
tolerance = MapleToFloat64(mapleKernelVec,args[4]);
/* call the external function */
r = newton(wrap_f, wrap_fprime, guess, tolerance);
/* convert the result back to Maple format */
return( ToMapleFloat(mapleKernelVec,r) );
}
This was compiled on Solaris with the following command.
cc -G mwrap_newton.c -o libwrapnewton.so -lnewton -lmaplec -I $MAPLE/extern/include
This can be compiled using MSVC with the following command.
cl -Gz newton.c mwrap_newton.c -Fe:mwrap_newton.dll -Gz -I$MAPLE_ROOT/exterrn/include -link $MAPLE_ROOT/bin.win/maplec.lib -dll -export:MWRAP_newton -out:mwrap_newton.dll
|
The external function does not need to be described in Maple because the parameter conversion is done in the manually written wrapper. The arguments needed by define_external are:
1. Name of the external function
2. Keyword 'MAPLE' to indicate that the wrapper was pregenerated
3. Name of the shared library (dll) containing the external function
>
|
|
|
|
Using Methods 1 and 2
|
|
Now that the newton procedure is defined, you can use it to solve the initial question. That is, you can find a root of the following equation.
| (3.2.3.4.1) |
>
|
|
You need the derivative of this equation.
>
|
|
| (3.2.3.4.2) |
Both f and fprime must be procedures. Turn the above polynomials into procedures.
>
|
|
| (3.2.3.4.3) |
Using an initial guess of x = 1, find a root to a tolerance of 1e-8.
>
|
|
| (3.2.3.4.4) |
Check to make sure the answer above works.
>
|
|
| (3.2.3.4.5) |
Plot this function to help you pick another guess.
>
|
|
>
|
|
| (3.2.3.4.6) |
>
|
|
| (3.2.3.4.7) |
Note: For the last root, the tolerance was reduced to avoid an infinite loop. The external source code should be changed to exit if it loops too many times; otherwise,Maple will not be able to recover since it cannot interrupt an external function. (Interrupting an external function that operates on Maple data structures could leave Maple in a bad state producing unexpected results, or crashing the program.)
>
|
|
| (3.2.3.4.8) |
|
|
|
|
Combine CodeGeneration and External Call
|
|
It is possible to combine External Calling with the CodeGeneration package to generate external functions, compile and access them on the fly. This may be useful for large problems, or moderate size functions that are called many times. In this example the series command is used to expand x^x about 0, up to order 3. The codegen[makeproc] command is used to make a procedure with parameter x from the result.
>
|
|
| (3.3.1) |
Use CodeGeneration[C] to generate C code. Optimize the Maple code first, and put the result in file test1.c. (Make sure you do not already have a file test1.c, because CodeGeneration appends to existing files.)
>
|
|
>
|
|
Create a shared library.
>
|
|
Create a link to the external function using the define_external command.
>
|
|
Test to see if you get the same result from functions f1 and f2.
>
|
|
| (3.3.2) |
In this case, the external function executes much faster than the Maple function using software floats. Timings of 50000 executions are as follows (on a 1.7GHz Pentium IV ).
>
|
|
| (3.3.3) |
| (3.3.4) |
For comparison, using Maple's built-in hardware float evaluator produces slightly faster results.
| (3.3.5) |
|
|
|
Exercises
|
|
|
Extend the Simple Example(s)
|
|
1. Choose one of the source files in the simple examples. Change the mult function to multiply double precision floating point numbers instead of integers.
2. Create an add function in the same file. Export both the add and mult functions in the new shared library. Note that the name add is acceptable to use in your dll, but the result of define_external should be assigned to myadd since add is already used in Maple. Similarly, you must put single quotes around the name add inside define_external.
|
|
Extend the Matrix Example
|
|
1. Write an external function that adds two matrices.
2. Write an external function that finds the biggest element inside a Matrix.
|
|
Pass Arguments By-Reference
|
|
Consider the following C function that doubles a number in-place.
void doubleme( int *n ) {
*n *= 2;
}
1. How would this be declared using define_external?
2. Write a custom wrapper to implement this function.
|
|
Search the Web
|
|
1. Find some code on the Web. There are many free libraries for doing all kinds of numerical operations available on the Web. Find code that handles FFTs, linear-optimization, statistics, or basic linear algebra. Use external calling to integrate this code into Maple.
|
|
Return to Index for Example Worksheets
|