workbench/libs/mathieeesingtrans/ieeespexp.c

   1 /*
   2     Copyright © 1995-2003, The AROS Development Team. All rights reserved.
   3     $Id$
   4 */
   5
   6 #include "mathieeesingtrans_intern.h"
   7
   8 /*****************************************************************************
   9
  10     NAME */
  11
  12         AROS_LH1(float, IEEESPExp,
  13
  14 /*  SYNOPSIS */
  15         AROS_LHA(float, y, D0),
  16
  17 /*  LOCATION */
  18         struct Library *, MathIeeeSingTransBase, 13, MathIeeeSingTrans)
  19
  20 /*  FUNCTION
  21         Calculate e^x
  22
  23     INPUTS
  24
  25     RESULT
  26
  27     BUGS
  28         IEEE single precision number
  29
  30         flags:
  31         zero     : result is zero
  32         negative : 0
  33         overflow : the result was out of range for the IEEE single precision
  34                    format
  35
  36     INTERNALS
  37         e^(>= 89): return 0x7f800000;
  38         e^(2^(<=-24)): return one;
  39
  40 *****************************************************************************/
  41 {
  42     AROS_LIBFUNC_INIT
  43
  44     const LONG ExpTable[] =
  45     {
  46         0x6da12cc2,  /* e^64  */
  47         0x568fa1fe,  /* e^32  */
  48         0x4b07975e,  /* e^16  */
  49         0x453a4f53,  /* e^8   */
  50         0x425a6481,  /* e^4   */
  51         0x40ec7325,  /* e^2   */
  52         0x402df854,  /* e^1   */
  53         0x3fd3094c,  /* e^(1/2) */
  54         0x3fa45af2,  /* e^(1/4) */
  55         0x3f910b02,  /* e^(1/8) */
  56         0x3f88415b,  /* e^(1/16) */
  57         0x3f84102b,  /* e^(1/32) */
  58         0x3f820405,  /* e^(1/64) */
  59         0x3f810101,  /* e^(1/128) */
  60         0x3f808040,  /* e^(1/256) */
  61         0x3f804010,  /* e^(1/512) */
  62         0x3f802004,  /* e^(1/1024) */
  63         0x3f801001,  /* e^(1/2048) */
  64         0x3f800800,  /* e^(1/4096) */
  65         0x3f800400,  /* e^(1/8192) */
  66         0x3f800200,  /* e^(1/16384) */
  67         0x3f800100,  /* e^(1/32768) */
  68         0x3f800080,  /* e^(1/65536) */
  69         0x3f800040,  /* e^(1/131072) */
  70         0x3f800020,  /* e^(1/262144) */
  71         0x3f800010,  /* e^(1/524288) */
  72         0x3f800008,  /* e^(1/1048576) */
  73         0x3f800004,  /* e^(1/2097152) */
  74         0x3f800002,  /* e^(1/4194304) */
  75         0x3f800001,  /* e^(1/8388608) */
  76     };
  77     ULONG Res, i;
  78     LONG Mantisse;
  79     char Exponent;
  80
  81     Exponent = ((y & IEEESPExponent_Mask) >> 23) -0x7f;
  82
  83     /* e^0 = 1, e^(2^(<=-24)) = 1 */
  84     if ( 0 == y || Exponent <= -24 ) return one;
  85
  86     /* e^(>= 89) = overflow) */
  87     if (Exponent > 6)
  88     {
  89         SetSR(Overflow_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
  90         return IEEESP_Pinfty;
  91     }
  92
  93     i = 6 - Exponent;
  94
  95     Mantisse = (y & IEEESPMantisse_Mask) << 9;
  96     Res = ExpTable[i++];
  97
  98     while ( 0 != Mantisse && i <= 29 )
  99     {
 100         /* is the highest bit set? */
 101         if ( Mantisse < 0 )
 102         {
 103             Res = IEEESPMul(Res, ExpTable[i]);
 104             if (0 == Res)
 105             {
 106                 SetSR(Zero_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
 107                 return Res;
 108             }
 109             if (IEEESP_Pinfty == Res)
 110             {
 111                 SetSR(Overflow_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
 112                 return Res;
 113             }
 114         }
 115         i++;
 116         Mantisse <<= 1;
 117     }
 118
 119     if ( y < 0) return IEEESPDiv(one, Res);
 120
 121     return Res;
 122
 123     AROS_LIBFUNC_EXIT
 124 }