workbench/libs/mathtrans/spexp.c

   1 /*
   2     Copyright © 1995-2003, The AROS Development Team. All rights reserved.
   3     $Id$
   4 */
   5
   6 #include "mathtrans_intern.h"
   7
   8 /*****************************************************************************
   9
  10     NAME */
  11
  12         AROS_LH1(float, SPExp,
  13
  14 /*  SYNOPSIS */
  15         AROS_LHA(float, fnum1, D0),
  16
  17 /*  LOCATION */
  18         struct Library *, MathTransBase, 13, MathTrans)
  19
  20 /*  FUNCTION
  21         Calculate e^x
  22
  23     INPUTS
  24
  25     RESULT
  26         Motorola fast floating point number
  27
  28         flags:
  29         zero     : result is zero
  30         negative : 0
  31         overflow : the result was out of range for the ffp-format
  32
  33     BUGS
  34
  35     INTERNALS
  36         e^(>= 44): return FFP_Pinfty;
  37
  38 *****************************************************************************/
  39 {
  40     AROS_LIBFUNC_INIT
  41
  42     const LONG ExpTable[] =
  43     {
  44         0x8fa1fe6f,  /* e^32  */
  45         0x87975e58,  /* e^16  */
  46         0xba4f534c,  /* e^8   */
  47         0xda648146,  /* e^4   */
  48         0xec732543,  /* e^2   */
  49         0xadf85442,  /* e^1   */
  50         0xD3094C41,  /* e^(1/2) */
  51         0xA45aF241,  /* e^(1/4) */
  52         0x910b0241,  /* e^(1/8) */
  53         0x88415b41,  /* e^(1/16) */
  54         0x84102b41,  /* e^(1/32) */
  55         0x82040541,  /* e^(1/64) */
  56         0x81010141,  /* e^(1/128) */
  57         0x80804041,  /* e^(1/256) */
  58         0x80401041,  /* e^(1/512) */
  59         0x80200441,  /* e^(1/1024) */
  60         0x80100141,  /* e^(1/2048) */
  61         0x80080041,  /* e^(1/4096) */
  62         0x80040041,  /* e^(1/8192) */
  63         0x80020041,  /* e^(1/16384) */
  64         0x80010041,  /* e^(1/32768) */
  65         0x80008041,  /* e^(1/65536) */
  66         0x80004041,  /* e^(1/131072) */
  67         0x80002041,  /* e^(1/262144) */
  68         0x80001041,  /* e^(1/524288) */
  69         0x80000841,  /* e^(1/1048576) */
  70         0x80000441,  /* e^(1/2097152) */
  71         0x80000241,  /* e^(1/4194304) */
  72         0x80000141   /* e^(1/8388608) */
  73     };
  74
  75     ULONG Res, i;
  76     LONG Mantisse;
  77     char Exponent;
  78
  79     Exponent = (fnum1 & FFPExponent_Mask) - 0x41;
  80
  81     /* e^0 = 1, e^(2^(<=-24)) = 1 */
  82     if ( 0 == fnum1 || Exponent <= -24 ) return one;
  83
  84     /* e^(>= 44) = overflow = infinity) */
  85     if (Exponent > 5)
  86     {
  87         SetSR(Overflow_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
  88         return FFP_Pinfty;
  89     }
  90
  91     i = 5 - Exponent;
  92
  93     Mantisse = (fnum1 & FFPMantisse_Mask);
  94     Res = ExpTable[i++];
  95     Mantisse <<= 1;
  96
  97     while ( 0 != Mantisse && i <= 28 )
  98     {
  99         /* is the highest bit set? */
 100         if ( Mantisse < 0 )
 101         {
 102             Res = SPMul(Res, ExpTable[i]);
 103             if (0 == Res)
 104             {
 105                 SetSR(Zero_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
 106                 return Res;
 107             }
 108             if (FFP_Pinfty == Res)
 109             {
 110                 SetSR(Overflow_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
 111                 return Res;
 112             }
 113         }
 114         i++;
 115         Mantisse <<= 1;
 116     }
 117
 118     if ( (char) fnum1 < 0) return SPDiv(one, Res);
 119
 120     return Res;
 121
 122     AROS_LIBFUNC_EXIT
 123 }