src/float/transcendental/cl_LF_ln10.cc

   1 // cl_ln10().
   2
   3 // General includes.
   4 #include "cl_sysdep.h"
   5
   6 // Specification.
   7 #include "cl_F_tran.h"
   8
   9
  10 // Implementation.
  11
  12 #include "cl_lfloat.h"
  13 #include "cl_LF.h"
  14
  15 static inline const cl_LF compute_ln10_old (uintC len)
  16 {
  17         return ln(cl_I_to_LF(10,len));
  18 }
  19
  20 // ln 10 =
  21 // = 46 atanh(1/31) + 34 atanh(1/49) + 20 atanh(1/161)
  22 // = 478 atanh(1/251) + 180 atanh(1/449) - 126 atanh(1/4801) + 206 atanh(1/8749)
  23
  24 static inline const cl_LF compute_ln10_p235 (uintC len)
  25 {
  26         var uintC actuallen = len+1;
  27         return shorten(  The(cl_LF)(46 * cl_atanh_recip(31,actuallen))
  28                        + The(cl_LF)(34 * cl_atanh_recip(49,actuallen))
  29                        + The(cl_LF)(20 * cl_atanh_recip(161,actuallen)),
  30                        len
  31                       );
  32 }
  33
  34 static inline const cl_LF compute_ln10_p2357 (uintC len)
  35 {
  36         var uintC actuallen = len+1;
  37         return shorten(  The(cl_LF)(478 * cl_atanh_recip(251,actuallen))
  38                        + The(cl_LF)(180 * cl_atanh_recip(449,actuallen))
  39                        - The(cl_LF)(126 * cl_atanh_recip(4801,actuallen))
  40                        + The(cl_LF)(206 * cl_atanh_recip(8749,actuallen)),
  41                        len
  42                       );
  43 }
  44
  45 #define compute_ln10 compute_ln10_p2357
  46
  47 const cl_LF cl_ln10 (uintC len)
  48 {
  49         var uintC oldlen = TheLfloat(cl_LF_ln10)->len; // vorhandene Länge
  50         if (len < oldlen)
  51                 return shorten(cl_LF_ln10,len);
  52         if (len == oldlen)
  53                 return cl_LF_ln10;
  54
  55         // TheLfloat(cl_LF_ln10)->len um mindestens einen konstanten Faktor
  56         // > 1 wachsen lassen, damit es nicht zu häufig nachberechnet wird:
  57         var uintC newlen = len;
  58         oldlen += floor(oldlen,2); // oldlen * 3/2
  59         if (newlen < oldlen)
  60                 newlen = oldlen;
  61
  62         // gewünschte > vorhandene Länge -> muß nachberechnen:
  63         cl_LF_ln10 = compute_ln10(newlen);
  64         return (len < newlen ? shorten(cl_LF_ln10,len) : cl_LF_ln10);
  65 }