[net.micro.amiga] LatFFP program SOURCE

bobp@amiga.UUCP (Robert S. Pariseau) (11/06/85)

TITLE:  LatFFP program SOURCE (LONG!)

The program below shows how to access the Motorola Fast Floating Point
libraries from V1.0 release Lattice C for the Amiga.  As you will see
when you run it, the performance improvement gained by using the
FFP routines (even in this kludgey fashion) is generally around a
factor of 10!

Part of that improvement is just the difference in precision.  The
Lattice stuff provides a 64 bit software implementation of the IEEE
format.  The FFP stuff provides a 32 bit implementation of the
Motorola format.

The variables used for FFP math are defined as ints.  You can't use
FLOATs because V1.0 Lattice C converts FLOAT to DOUBLE during
expression evaluation and when passing arguments.

------------------------Program Notes:

The program will compile and link cleanly using the stuff on the
standard V1.0 Lattice C for Amiga disk.  The Make script in the
examples directory will do all the work for you.  Since my C disk
is rather full, and since I like the dramatic increase in speed,
I usually do my work in ram disk as follows:

  1> cd df1:                               [my C disk]
  1> copy examples/Make to :               [more convenient in root]
  
  1> copy LatFFP.c to ram:
  1> execute Make ram:LatFFP
  1> copy ram:LatFFP to df1:

I've also added a "stack 20000" command in my startup script
(s/startup-sequence) to make sure I don't have to worry about stack
overflows.

The program shows 3 bugs in the V1.0 stuff and includes one kludge.
The bugs and kludge are maked in the source.  Bug (1): V1.0 Lattice
C doesn't properly handle successive assignments of constants.  Use
expressions.  Bug (2): V1.0 FFP doesn't correctly return a result
from SPCmp().  Fake it with subtraction and SPTst() or write your
own compare for now.  Bug (3): V1.0 FFP doesn't correctly return
the cosine result from SPSincos().  Use SPCos() instead.

The kludge is the set of union definitions so that we can use the
same variables for related Lattice and FFP expressions and for
the conversions in and out of IEEE format.  Note that we use
the Lattice IEEE based routines in printf() to get our output.

Note however, that the conversion routines are in the RAM based
MathTrans library.  Therefore, if you only want to use the WCS
based MathFFP library, you are on your own for conversion.

------------------------Program Source Follows:

/***********************************************************************
 *  LatFFP -- Program to show the use of Motorola Fast Floating Point
 *            math libraries with V1.0 Lattice C for the Amiga in
 *            comparison to the DOUBLE precision IEEE floating point
 *            math which is built in to the C.  The FFP format is
 *            a 32 bit format.
 *
 *  Larry Hildenbrand -- Nov. 4, 1985
 *  Bob Pariseau      -- Nov. 4, 1985  (minor editorial corrections)
 *
 ***********************************************************************/

#include <exec/types.h>
#include <math.h>


/* ??? #define   E      2.718281828459045 */  /* V1.0 Lattice C BUG!  See */
#define   PIME   0.423310826                  /* below.  PIME == PI - E.  */


/****  MAY BE BROKEN OUT INTO SEPARATE #include FILE ****/

extern   int     SPFix();
extern   int     SPFlt();
extern   int     SPCmp();
extern   int     SPTst();
extern   int     SPAbs();
extern   int     SPNeg();
extern   int     SPAdd();
extern   int     SPSub();
extern   int     SPMul();
extern   int     SPDiv();

extern   int     SPAtan();
extern   int     SPSin();
extern   int     SPCos();
extern   int     SPTan();
extern   int     SPSincos();
extern   int     SPSinh();
extern   int     SPCosh();
extern   int     SPTanh();
extern   int     SPExp();
extern   int     SPLog();
extern   int     SPPow();
extern   int     SPSqrt();
extern   int     SPTieee();
extern   int     SPFieee();

/********************************************************/


char st1[80] = "3.1415926535897";
char st2[80] = "2.718281828459045";


int MathBase;        /* Basic FFP lib pointer */
int MathTransBase;   /* Transcendental FFP lib pointer */

int dots_good = 0;


union kludge1        /* Can't use FLOAT directly for FFP stuff */
{                    /* because V1.0 Lattice converts FLOAT to */
  FLOAT num1;        /* DOUBLE in expressions and when passing */
  int   i1;          /* parameters.                            */
} k1;

union kludge2
{
  FLOAT num2;
  int   i2;
} k2;

union kludge3
{
  FLOAT num3;
  int   i3;
} k3;

union kludge4
{
  FLOAT num4;
  int   i4;
} k4;

union kludge5
{
  FLOAT num5;
  int   i5;
} k5;

union kludge6
{
  FLOAT num6;
  int   i6;
} k6;



show_dot() { if( ++dots_good == 1000) { dots_good = 0; printf(".");} }



show_result( num ) FLOAT num; { printf("\nResult = %f", num); }



show_result_ffp(in_val)   /* Convert to IEEE and display */
  int in_val;
{
  union kludge_sr
  {
    FLOAT new_iv_f;
    int   new_iv_i;
  } k;

  k.new_iv_i = SPTieee(in_val);
  show_result(k.new_iv_f);
}




main()                  /* Lattice/FFP test code */
{
   UWORD i;
   int i3;

   printf("C-ROM & Shared RAM Library Test Facility for the Amiga-Lattice Basic/Trans Math Libraries");

/*****************************************/
/*****************************************/
/***                                   ***/
/***  OPEN ROM AND RAM FFP LIBRARIES   ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   if((MathBase = OpenLibrary("mathffp.library", 0)) < 1 ) {
       printf("\n\n*** ERROR ***  Can't open mathffp.library: vector = %lx\n", MathBase);
       exit();
   }
   else {
       printf("\n\nSuccessfully opened mathffp.library: vector = %lx\n", MathBase);
   }

   if((MathTransBase = OpenLibrary("mathtrans.library", 0)) < 1 ) {
       printf("\n\n*** ERROR ***  Can't open mathtrans.library: vector = %lx\n", MathTransBase);
       CloseLibrary(MathBase);
       exit();
   }
   else {
       printf("\n\nSuccessfully opened mathtrans.library: vector = %lx\n", MathTransBase);
   }


/*****************************************/
/*****************************************/
/***                                   ***/
/***  COMPILER & FFP S.P. ADDITION     ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;                  /* V1.0 Lattice C BUG!  Can't have two  */
   k2.num2 = k1.num1 - PIME;      /* constant assignments in a row.  Fake */
                                  /* it by making the second be an        */
                                  /* expression!                          */


   printf("\n\n50,000 additions of %s to %s (Compiler Interface)\n", st1, st2);
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k3.num3 = k1.num1 + k2.num2;
      show_dot();
   }
   show_result( k3.num3 );


   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   printf("\n\n50,000 additions of %s to %s (Function Interface)\n", st1, st2 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k3.i3 = SPAdd(k2.i2, k1.i1);
      show_dot();
   }
   show_result_ffp( k3.i3 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  COMPILER & FFP S.P. SUBTRACTION  ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;
   k2.num2 = k1.num1 - PIME;

   printf("\n\n50,000 subtractions of %s from %s (Compiler Interface)\n", st1, st2 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k3.num3 = k2.num2 - k1.num1;
      show_dot();
   }
   show_result( k3.num3 );


   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   printf("\n\n50,000 subtractions of %s from %s (Function Interface)\n", st1, st2 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k3.i3 = SPSub(k1.i1, k2.i2);
      show_dot();
   }
   show_result_ffp( k3.i3 );

   
/*****************************************/
/*****************************************/
/***                                   ***/
/***  COMPILER & FFP S.P. MULTIPLYS    ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;
   k2.num2 = k1.num1 - PIME;

   printf("\n\n50,000 multiplies of %s by %s (Compiler Interface)\n", st1, st2 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k3.num3 = k1.num1 * k2.num2;
      show_dot();
   }
   show_result( k3.num3 );
   

   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   printf("\n\n50,000 multiplies of %s by %s (Function Interface)\n", st1, st2 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k3.i3 = SPMul(k2.i2, k1.i1);
      show_dot();
   }
   show_result_ffp( k3.i3 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  COMPILER & FFP S.P. DIVISION     ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;
   k2.num2 = k1.num1 - PIME;

   printf("\n\n50,000 divides of %s by %s (Compiler Interface)\n", st1, st2 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k3.num3 = k1.num1 / k2.num2;
      show_dot();
   }
   show_result( k3.num3 );


   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   printf("\n\n50,000 divides of %s by %s (Function Interface)\n", st1, st2 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k3.i3 = SPDiv(k2.i2, k1.i1);
      show_dot();
   }
   show_result_ffp( k3.i3 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  COMPILER & FFP S.P. TRUNCATION   ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;
   k2.num2 = k1.num1 - PIME;

   printf("\n\n50,000 fixes of %s (Compiler Interface)\n", st1 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      i3 = (int) k1.num1;
      if( ++dots_good == 1000) { dots_good = 0; printf(".");}
   }
   printf("\nResult = %d", i3 );


   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   printf("\n\n50,000 fixes of %s (Function Interface)\n", st1 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      i3 = SPFix(k1.i1);
      if( ++dots_good == 1000) { dots_good = 0; printf(".");}
   }
   printf("\nResult = %d", i3);


/*****************************************/
/*****************************************/
/***                                   ***/
/***  COMPILER & FFP S.P. FLOATATION   ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   i3 = 5;

   printf("\n\n50,000 floats of %d (Compiler Interface)\n", i3 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
       k4.num4 = (FLOAT) i3;
      show_dot();
   }
   show_result( k4.num4 );


   printf("\n\n50,000 floats of %d (Function Interface)\n", i3 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
       k4.i4 = SPFlt(i3);
      show_dot();
   }
   show_result_ffp( k4.i4 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  COMPILER & FFP S.P. NEGATION     ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;
   k2.num2 = k1.num1 - PIME;

   printf("\n\n50,000 negates of %s (Compiler Interface)\n", st2 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k4.num4 = -k2.num2;
      show_dot();
   }
   show_result( k4.num4 );


   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   printf("\n\n50,000 negates of %s (Function Interface)\n", st2 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k4.i4 = SPNeg(k2.i2);
      show_dot();
   }
   show_result_ffp( k4.i4 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  COMPILER & FFP S.P. ABSOLUTE VAL ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;
   k4.num4 = PIME - k1.num1;

   printf("\n\n50,000 absolute values of %f (Compiler Interface)\n", k4.num4 );
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k5.num5 = fabs(k4.num4);
      show_dot();
   }
   show_result( k5.num5 );


   printf("\n\n50,000 absolute values of %f (Function Interface)\n", k4.num4 );
   k4.i4 = SPFieee(k4.i4);
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k5.i5 = SPAbs(k4.i4);
      show_dot();
   }
   show_result_ffp( k5.i5 );


   printf("\n\n*** HIT RETURN TO CONTINUE ***");   getch();


/*****************************************/
/*****************************************/
/***                                   ***/
/***  COMPILER & FFP S.P. COMPARE      ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;
   k2.num2 = k1.num1 - PIME;

   if (k2.num2 >= k1.num1)
      printf("\n\n%f is greater than or equal to %f (Compiler Interface)\n", k2.num2, k1.num1);
   else
      printf("\n\n%f is less than %f (Compiler Interface)\n", k2.num2, k1.num1);


   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   printf("\n*** SPCmp(k2.i2, k1.i1) returned %d ***", SPCmp(k2.i2, k1.i1));


   if (SPCmp(k2.i2, k1.i1))     /* V1.0 FFP Bug.  SPCmp() broken. */

   {  k1.num1 = PI;
      k2.num2 = k1.num1 - PIME;
      printf("\n\n%f is greater than or equal to %f (Function Interface)\n", k2.num2, k1.num1);
   }
   else
   {  k1.num1 = PI;
      k2.num2 = k1.num1 - PIME;
      printf("\n\n%f is less than %f (Function Interface)\n", k2.num2, k1.num1);
   }


   if (k1.num1 >= k2.num2)
      printf("\n\n%f is greater than or equal to %f (Compiler Interface)\n", k1.num1, k2.num2);
   else
      printf("\n\n%f is less than %f (Compiler Interface)\n", k1.num1, k2.num2);


   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   if (SPCmp(k1.i1, k2.i2))
   {  k1.num1 = PI;
      k2.num2 = k1.num1 - PIME;
      printf("\n\n%f is greater than or equal to %f (Function Interface)\n", k1.num1, k2.num2);
   }
   else
   {  k1.num1 = PI;
      k2.num2 = k1.num1 - PIME;
      printf("\n\n%f is less than %f (Function Interface)\n", k1.num1, k2.num2);
   }


/*****************************************/
/*****************************************/
/***                                   ***/
/***  COMPILER & FFP S.P. TEST         ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   if (k2.num2)
      printf("\n\n%f is not equal to 0.0 (Compiler Interface)\n", k2.num2);
   else
      printf("\n\n%f is equal to 0.0 (Compiler Interface)\n", k2.num2);


   k2.i2 = SPFieee(k2.i2);

   if (SPTst(k2.i2))
   {  k1.num1 = PI;
      k2.num2 = k1.num1 - PIME;
      printf("\n\n%f is not equal to 0.0 (Function Interface)\n", k2.num2);
   }
   else
   {  k1.num1 = PI;
      k2.num2 = k1.num1 - PIME;
      printf("\n\n%f is equal to 0.0 (Function Interface)\n", k2.num2);
   }


   k2.num2 = 0.0;

   if (k2.num2)
      printf("\n\n%f is not equal to 0.0 (Compiler Interface)\n", k2.num2);
   else
      printf("\n\n%f is equal to 0.0 (Compiler Interface)\n", k2.num2);


   k2.i2 = SPFieee(k2.i2);

   if (SPTst(k2.i2))
   {  k2.num2 = 0.0;
      printf("\n\n%f is not equal to 0.0 (Function Interface)\n", k2.i2);
   }
   else
   {  k2.num2 = 0.0;
      printf("\n\n%f is equal to 0.0 (Function Interface)\n", k2.i2);
   }


/*****************************************/
/*****************************************/
/***                                   ***/
/***  FFP S.P. SQUARE ROOT             ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;
   k2.num2 = k1.num1 - PIME;

   printf("\n\n50,000 square roots of %f\n", k2.num2);
   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);
   for( dots_good = 0, i= 1; i < 50000; i++ )
   {
      k3.i3 = SPSqrt( k2.i2 );
      show_dot();
   }
   show_result_ffp( k3.i3 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  FFP S.P. NATURAL LOGARITHM       ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   printf("\n\n40,000 logarithms of %s\n", st1 );
   for( dots_good = 0, i= 1; i < 40000; i++ )
   {
      k3.i3 = SPLog( k1.i1 );
      show_dot();
   }
   show_result_ffp( k3.i3 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  FFP S.P. EXPONENT (BASE e)       ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   printf("\n\n40,000 exponents of %s\n", st1 );
   for( dots_good = 0, i= 1; i < 40000; i++ )
   {
      k3.i3 = SPExp( k1.i1 );
      show_dot();
   }
   show_result_ffp( k3.i3 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  FFP S.P. SINE, COSINE & TANGENT  ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;
   k2.num2 = k1.num1 - PIME;
   k1.num1 = k1.num1 / 6.0;
   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   printf("\n\n20,000 sines, cosines and tangents of %s / 6 radians\n", st1 );
   for( dots_good = 0, i= 1; i < 20000; i++ )
   {
      k2.i2 = SPSin(k1.i1);
      k3.i3 = SPCos(k1.i1);
      k4.i4 = SPTan(k1.i1);
      k5.i5 = SPSincos(&k6.i6, k1.i1);  /* V1.0 FFP BUG!  Cosine return  */
                                        /* value (k6.i6) of SPSincos is  */
                                        /* broken.  Function result      */
                                        /* (sine -- k5.i5) is OK.        */

      show_dot();
   }
   show_result_ffp( k2.i2 );   show_result_ffp( k3.i3 );
   show_result_ffp( k4.i4 );
   show_result_ffp( k5.i5 );   show_result_ffp( k6.i6 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  FFP S.P. ARCTANGENT              ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   printf("\n\n20,000 arctangents of the tangent of %s / 6 radians\n", st1 );
   for( dots_good = 0, i= 1; i < 20000; i++ )
   {
      k2.i2 = SPAtan(k4.i4);
      show_dot();
   }
   show_result_ffp( k2.i2 );


/***************************************************/
/***************************************************/
/***                                             ***/
/***  FFP S.P. HYPERBOLIC SINE, COSINE & TANGENT ***/
/***                                             ***/
/***************************************************/
/***************************************************/

   k1.num1 = PI;
   k2.num2 = k1.num1 - PIME;
   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   printf("\n\n20,000 hyperbolic sines, cosines and tangents of %s radians\n", st1 );
   for( dots_good = 0, i= 1; i < 20000; i++ )
   {
      k2.i2 = SPSinh( k1.i1 );
      k3.i3 = SPCosh( k1.i1 );
      k4.i4 = SPTanh( k1.i1 );
      show_dot();
   }
   show_result_ffp( k2.i2 );   show_result_ffp( k3.i3 );   show_result_ffp( k4.i4 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  FFP S.P. POWER FUNCTION          ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   k1.num1 = PI;
   k2.num2 = k1.num1 - PIME;
   k1.i1 = SPFieee(k1.i1);
   k2.i2 = SPFieee(k2.i2);

   printf("\n\n10,000 %s raised to the %s power\n", st1, st2 );
   for( dots_good = 0, i= 1; i < 10000; i++ )
   {
      k3.i3 = SPPow( k2.i2, k1.i1 );
      show_dot();
   }
   show_result_ffp( k3.i3 );


/*****************************************/
/*****************************************/
/***                                   ***/
/***  CLOSE ROM AND RAM FFP LIBRARIES  ***/
/***                                   ***/
/*****************************************/
/*****************************************/

   RemLibrary(MathTransBase);    /* Mark lib for Expunge() from RAM upon  */
                                 /* CloseLibrary() by last opener.  Else  */
                                 /* lib stays in RAM for others to use    */
                                 /* quickly (no need to go to disk) until */
                                 /* AllocMem() finds it needs the memory  */
                                 /* for some other purpose.               */
  
   CloseLibrary(MathTransBase);  /* Close transcendental math RAM library */



   CloseLibrary(MathBase);       /* Close basic math ROM library */


   printf("\n\nEnd C-ROM & Shared RAM Library Test (LATTICE) \n");
}

-----------------------That's all for now!