The SinSurface application example

Introduction

Below we describe the SinSurface program example. The actual computation is performed by the functions calcPlot, sinFun2 and their help functions.

The two functions calcPlot and sinFun2 in the SinSurface package will be translated to C++ together with the declaration of the global array xyMatrix.

The array xyMatrix represents a 21x21 grid on which which the numeric function sinFun2 will be computed.

The function calcPlot accepts four arguments which are coordinates that describe a square in the x-y plane, and one counter (iter) to make the function repeat the computation as many times as necessary for measuring execution time. For each point on a 21x21 grid in that square, the numeric function sinFun2 is called to compute a value that is stored as a matrix element in the matrix representing the grid.

The function sinFun2 computes essentially the same values as Sin(x + y), but in a more complicated way using a rather large expression obtained through conversion of the arguments into polar coordinates (through arcTan) and then using series expansion of both Sin and Cos in 10 terms. The resulting large symbolic expression (more than a page) becomes the body of sinFun2, and is then used as input to CompileEvaluateFunction (see Section 7.5 and Section 7.6.3 in the manual) to generate efficient C++ code.

The SinSurface example application computes a function (here sinFun2) over a 2-D grid. The function values are first stored in the matrix xyMatrix before being plotted. The execution of compiled C++ code for the function sinFun2 is approximately 500-1000 times faster than evaluating the same function interpretively within Mathematica.

To run this example, evaluate this notebook cell by cell or the whole at once.

Initialization

Load the MathCode package

First put a Needs statement, to make sure that the MathCode application is loaded:

Needs["MathCode`"]

Check current directory

Check the current directory, since a number of files will be placed there during the code generation process.

Directory[]

You might want to set the directory to some place where all generated files are put, e.g. the one below, or some other directory.

SetDirectory[$MCRoot<>$PathnameSeparator<>"Demos"<>$PathnameSeparator<>"SinSurface"];

Check whether this file is in the current directory.

ExistFileQ["SinSurface.nb"]

Start of the SinSurface package

The SinSurface package starts in the usual way by a BeginPackage declaration which references other packages.The MathCodeContexts are needed in order to call the code generation related functions.

\(\(BeginPackage["SinSurface`", \ {MathCodeContexts}]; \ 

    Off[General::spell]; 
Clear["SinSurface`*"]; \ \)\)

Exported symbols

Define possibly exported symbols. Even though it is not necessary here, we enclose these names within a Begin["SinSurface‘"] ... End[] kind of “context bracket”, since this can be put into a cell in the notebook, and conveniently re-evaluated (only this cell!) when a new name is added to the list below.

Begin["SinSurface`"]\(\ xyMatrix; \ 
\ calcPlot; 
\ sinFun1; 
\ sinFun2; \ 
\ arcTan; \ 

    \ sin; \ 
\ cos; 
\ plot; \ 
\ cplus; \ 
\ plotfile; 
End[]\)

Setting compilation options

This defines how the functions and variables in the package should be compiled to C++. By default, all typed variables and functions are compiled. However, the compilation process can be controlled in a more detailed manner by setting compilation options via SetCompilationOptions. For example, in this package the function sinFun2 should be symbolically evaluated before being translated to code since it contains symbolic operations; the functions sin, cos, and arcTan should not be compiled at all since they are expanded within the body of sinFun2. The remaining typed function, calcPlot, will be compiled in the normal way.

\(\(SetCompilationOptions[
\ \ EvaluateFunctions -> {sinFun2}, 
\ \ 
      UnCompiledFunctions -> {sin, cos, arcTan}, 
\ \ 
      MainFileAndFunction -> "int main(){return 0;}"
]; \)\)

The body of the SinSurface package

Begin the implementation section of the SinSurface package, where functions are defined. This is usually private, to avoid accidental name shadowing due to identical local variables in several packages.

Begin["SinSurface`Private`"];

Global variables

Declare global variables.

Declare[
  Real[21,21]@  xyMatrix 
];

sin, cos

Taylor-expanded sin and cos functions.

sin[Real@ x_ ]->Real := Normal[Series[Sin[z], {z,0,10}]] /. z -> x ;
cos[Real@ x_ ]->Real := Normal[Series[Cos[z], {z,0,10}]] /. z -> x ;

arcTan

Conversion of grid point to an angle, called by sinFun2.

arcTan[Real@ x_,Real@ y_]->Real :=  (
    If[x < 0, Pi, 0] + If[ x == 0, Sign[y]*Pi/2,  ArcTan[y/x] ]
  );

sinFun2

The function sinFun2 is the function to be computed and plotted, called by calcPlot. It provides a more complicated and computationally heavier way (series expansion) to calculate approximately the same as Sin[x+y]. This gives an example of a combination of symbolic and numeric operations as well as a rather standard mix of arithmetic operations. The expanded symbolic expression which comprises the body of SinFun2 is between 1 and 2 pages long when printed.

Note that the types of local variables to sinFun2 need not necessarily be declared since CompileEvaluateFunction will symbolically expand the whole body before translating it to C++ code.

Note also that a function to be symbolically expanded before final code generation should be without side effects, e.g. assignment to global variables or input/output, since otherwise these actions will be performed during the symbolic evaluation, which would be problematic.

sinFun2[Real@ x_,Real@ y_]->Real := Block[
   {
   	Real@ {r,xx,yy}
   },
    r = Sqrt[x^2+y^2];
    xx =  r*cos[arcTan[x,y]];
    yy = r*sin[arcTan[x,y]]; 
    sin[xx+yy]
  ];

calcPlot

The function calcPlot calculates data for a plot of sinFun2 over a 21x21 grid, which is returned as a 21x21 array.

calcPlot[Real@ xmin_, Real@ xmax_, Real@ ymin_,
 Real@ ymax_, Integer@ iter_]->Real[_,_] := 
 Module[{  
    Integer@  n = 20,
    Real@ {x,y}, 
    Integer@ {i, j, count}
 },
   For[count=1,count<=iter,count=count+1,
   	For[i=1, i<=(n+1),  i=i+1, 
      For[j=1, j<=(n+1), j=j+1,
	    x = xmin+(xmax-xmin)*(i-1)/n; 
	    y = ymin+(ymax-ymin)*(j-1)/n;
	    xyMatrix[[i,j]] = sinFun2[x,y]
       ]
    ]
   ];
   xyMatrix
 ];

End of the SinSurface package

End[]
EndPackage[];

Execution

We first execute the application interpretively within Mathematica, and then compile the key function and execute the application again. Then we compile the application to C++, build an executable, and call the same functions from Mathematica via MathLink.

Mathematica evaluation

Let Mathematica calculate a plot.

meval = Timing[plot = calcPlot[-2., 2., -2., 2., 1] ][[1]]

Perform the plot.

ListPlot3D[plot];

Using Mathematica standard Compile[]

Redefine sinFun2 to become a compiled version, using Mathematica standard Compile[].

sinFun2 = Compile[{x, y}, Evaluate[sinFun2[x, y]]];

Timing[plot = calcPlot[-2., 2., -2., 2., 1]; ]

compeval = %[[1]];

sinFun2 =.

Fully Expanded Symbolic expression to be compiled to C++

This shows the size of sinFun2 fully expanded. Be prepared for a big output.

sinFun2[a, b]

Using the MathCode code generator

Compile the SinSurface package:

CompilePackage["SinSurface"]

The warnings produced by CompilePackage is due to the lack of type information for the local variables of SinFun2. Since we have chosen to perform symbolic expansion SinFun2 will be evaluated and the local variables therefore disappear.

The generated C++ code

Below follows the generated C++ code from the SinSurface program. Notice that the package name becomes a prefix of the name of the generated C++ function. Thus the Mathematica function SinSurface‘sinFun2 becomes SinSurface_TsinFun2 in C++.

The generated code from SinSurface‘sinFun2 is produced from a large expression by CompileEvaluateFunction. Common subexpression elimination is performed by the code generator, producing many temporary variables and subexpressions which can be seen in the body of the C++ function SinSurface_TsinFun2.

By contrast, the C++ code in the body of function SinSurface_TcalcPlot produced by CompileFunction from the Mathematica function SinSurface‘calcPlot, follows the structure of the original code quite closely.

We give a command to type out the text of the generated C++ file:

!! SinSurface.cc

Compiling and linking the C++ code

The command MakeBinary compiles the generated C++ code using some C++ compiler (in this particular case the MicroSoft Visual C++, version 4.0, for Windows95). The object code is by default linked into two executables: SinSurfaceml.exe for calling the compiled code via MathLink, and SinSurfacemain.exe as a stand-alone executable.

MakeBinary[];

Installing and connecting to Mathematica

The command InstallCode installs and connects the external process containing the compiled and linked SinSurface code.

installedLink = InstallCode["SinSurface"]

Define a time measurement function called AbsTime:

\(SetAttributes[AbsTime, HoldFirst]; 


    AbsTime[x_]\  := \ 
      Module[{
\ \ start, res
\ }, 
\ \ \ start\  = \ AbsoluteTime[]; 

        \ \ \ res = x; 

        \ \ \ {\((AbsoluteTime[] - start)\)\ Second, res}
]; \)

Execution of generated C++ code

Perform the computation and the plot:

(plot  =  calcPlot[-2.0, 2.0, -2.0, 2.0, 500]; )  //  
      AbsTime

Since the external computation was performed 500 times, the time needed for one external computation is:

externaleval  =  %[[1]]/500

Just check that we get the right result

ListPlot3D[plot];

Performance Comparison

The performance between the three forms of execution are compared in the table below. The generated C++ code for this example is roughly 500-1000 times faster than standard interpreted Mathematica code, and almost 50-100 times faster than expression code compiled by the internal Mathematica Compile[] command.

TableForm[(({{meval, meval/externaleval}, {compeval, compeval/externaleval}, {externaleval, 1}})) ,TableHeadings -> {{"Standard Mathematica", "Compile[]", 
              "External C++ via MathLink"}, {"Time consumed", 
              "Relative"}}]

Cleanup

Uninstall[installedLink];
CleanMathCodeFiles[Confirm->False]