src/float/transcendental/cl_F_exp.cc

   1 // exp().
   2
   3 // General includes.
   4 #include "base/cl_sysdep.h"
   5
   6 // Specification.
   7 #include "cln/float.h"
   8
   9
  10 // Implementation.
  11
  12 #include "float/transcendental/cl_F_tran.h"
  13 #include "cln/float.h"
  14 #include "float/cl_F.h"
  15 #include "cln/lfloat.h"
  16 #include "float/lfloat/cl_LF.h"
  17
  18 #include "base/cl_inline.h"
  19 #include "float/lfloat/elem/cl_LF_minusp.cc"
  20 #include "float/lfloat/misc/cl_LF_exponent.cc"
  21
  22 namespace cln {
  23
  24 // Division durch ln(2).
  25 inline const cl_F_div_t cl_floor_ln2 (const cl_F& x)
  26 {
  27         // Bei 0<=x<1/2 kann man sofort q:=0 setzen.
  28         if (!minusp(x) && (float_exponent(x) < 0))
  29                 return cl_F_div_t(0,x);
  30         else
  31                 return floor2(x,cl_ln2(x));
  32 }
  33 inline const cl_LF_div_t cl_floor_ln2 (const cl_LF& x)
  34 {
  35         // Bei 0<=x<1/2 kann man sofort q:=0 setzen.
  36         if (!minusp_inline(x) && (float_exponent_inline(x) < 0))
  37                 return cl_LF_div_t(0,x);
  38         else
  39                 return floor2(x,The(cl_LF)(cl_ln2(x)));
  40 }
  41
  42 const cl_F exp (const cl_F& x)
  43 {
  44 // Methode:
  45 // d := (float-digits x),
  46 // Genauigkeit um sqrt(d)+max(integer-length(e)) Bits erhöhen,
  47 // (q,r) := (floor x ln(2))
  48 // Ergebnis ist exp(q*ln(2)+r) = (scale-float exp(r) q).
  49
  50         // Rechengenauigkeit erhöhen und durch ln(2) dividieren:
  51         if (longfloatp(x) && (TheLfloat(x)->len >= 84)) {
  52                 DeclareType(cl_LF,x);
  53                 var cl_LF_div_t q_r = cl_floor_ln2(extend(x,TheLfloat(x)->len+1));
  54                 var cl_I& q = q_r.quotient;
  55                 var cl_LF& r = q_r.remainder;
  56                 return cl_float(scale_float(expx_ratseries(r),q),x);
  57         } else {
  58                 var cl_F_div_t q_r = cl_floor_ln2(cl_F_extendsqrtx(x));
  59                 var cl_I& q = q_r.quotient;
  60                 var cl_F& r = q_r.remainder;
  61                 return cl_float(scale_float(expx_naive(r),q),x);
  62         }
  63 }
  64
  65 }  // namespace cln