root/arch/m68k/math-emu/fp_log.c

DEFINITIONS

1. fp_fsqrt
2. fp_fetoxm1
3. fp_fetox
4. fp_ftwotox
5. fp_ftentox
6. fp_flogn
7. fp_flognp1
8. fp_flog10
9. fp_flog2
10. fp_fgetexp
11. fp_fgetman

   1 /*
   2 
   3   fp_trig.c: floating-point math routines for the Linux-m68k
   4   floating point emulator.
   5 
   6   Copyright (c) 1998-1999 David Huggins-Daines / Roman Zippel.
   7 
   8   I hereby give permission, free of charge, to copy, modify, and
   9   redistribute this software, in source or binary form, provided that
  10   the above copyright notice and the following disclaimer are included
  11   in all such copies.
  12 
  13   THIS SOFTWARE IS PROVIDED "AS IS", WITH ABSOLUTELY NO WARRANTY, REAL
  14   OR IMPLIED.
  15 
  16 */
  17 
  18 #include "fp_emu.h"
  19 
  20 static const struct fp_ext fp_one =
  21 {
  22         .exp = 0x3fff,
  23 };
  24 
  25 extern struct fp_ext *fp_fadd(struct fp_ext *dest, const struct fp_ext *src);
  26 extern struct fp_ext *fp_fdiv(struct fp_ext *dest, const struct fp_ext *src);
  27 
  28 struct fp_ext *
  29 fp_fsqrt(struct fp_ext *dest, struct fp_ext *src)
  30 {
  31         struct fp_ext tmp, src2;
  32         int i, exp;
  33 
  34         dprint(PINSTR, "fsqrt\n");
  35 
  36         fp_monadic_check(dest, src);
  37 
  38         if (IS_ZERO(dest))
  39                 return dest;
  40 
  41         if (dest->sign) {
  42                 fp_set_nan(dest);
  43                 return dest;
  44         }
  45         if (IS_INF(dest))
  46                 return dest;
  47 
  48         /*
  49          *               sqrt(m) * 2^(p)        , if e = 2*p
  50          * sqrt(m*2^e) =
  51          *               sqrt(2*m) * 2^(p)      , if e = 2*p + 1
  52          *
  53          * So we use the last bit of the exponent to decide whether to
  54          * use the m or 2*m.
  55          *
  56          * Since only the fractional part of the mantissa is stored and
  57          * the integer part is assumed to be one, we place a 1 or 2 into
  58          * the fixed point representation.
  59          */
  60         exp = dest->exp;
  61         dest->exp = 0x3FFF;
  62         if (!(exp & 1))         /* lowest bit of exponent is set */
  63                 dest->exp++;
  64         fp_copy_ext(&src2, dest);
  65 
  66         /*
  67          * The taylor row around a for sqrt(x) is:
  68          *      sqrt(x) = sqrt(a) + 1/(2*sqrt(a))*(x-a) + R
  69          * With a=1 this gives:
  70          *      sqrt(x) = 1 + 1/2*(x-1)
  71          *              = 1/2*(1+x)
  72          */
  73         fp_fadd(dest, &fp_one);
  74         dest->exp--;            /* * 1/2 */
  75 
  76         /*
  77          * We now apply the newton rule to the function
  78          *      f(x) := x^2 - r
  79          * which has a null point on x = sqrt(r).
  80          *
  81          * It gives:
  82          *      x' := x - f(x)/f'(x)
  83          *          = x - (x^2 -r)/(2*x)
  84          *          = x - (x - r/x)/2
  85          *          = (2*x - x + r/x)/2
  86          *          = (x + r/x)/2
  87          */
  88         for (i = 0; i < 9; i++) {
  89                 fp_copy_ext(&tmp, &src2);
  90 
  91                 fp_fdiv(&tmp, dest);
  92                 fp_fadd(dest, &tmp);
  93                 dest->exp--;
  94         }
  95 
  96         dest->exp += (exp - 0x3FFF) / 2;
  97 
  98         return dest;
  99 }
 100 
 101 struct fp_ext *
 102 fp_fetoxm1(struct fp_ext *dest, struct fp_ext *src)
 103 {
 104         uprint("fetoxm1\n");
 105 
 106         fp_monadic_check(dest, src);
 107 
 108         return dest;
 109 }
 110 
 111 struct fp_ext *
 112 fp_fetox(struct fp_ext *dest, struct fp_ext *src)
 113 {
 114         uprint("fetox\n");
 115 
 116         fp_monadic_check(dest, src);
 117 
 118         return dest;
 119 }
 120 
 121 struct fp_ext *
 122 fp_ftwotox(struct fp_ext *dest, struct fp_ext *src)
 123 {
 124         uprint("ftwotox\n");
 125 
 126         fp_monadic_check(dest, src);
 127 
 128         return dest;
 129 }
 130 
 131 struct fp_ext *
 132 fp_ftentox(struct fp_ext *dest, struct fp_ext *src)
 133 {
 134         uprint("ftentox\n");
 135 
 136         fp_monadic_check(dest, src);
 137 
 138         return dest;
 139 }
 140 
 141 struct fp_ext *
 142 fp_flogn(struct fp_ext *dest, struct fp_ext *src)
 143 {
 144         uprint("flogn\n");
 145 
 146         fp_monadic_check(dest, src);
 147 
 148         return dest;
 149 }
 150 
 151 struct fp_ext *
 152 fp_flognp1(struct fp_ext *dest, struct fp_ext *src)
 153 {
 154         uprint("flognp1\n");
 155 
 156         fp_monadic_check(dest, src);
 157 
 158         return dest;
 159 }
 160 
 161 struct fp_ext *
 162 fp_flog10(struct fp_ext *dest, struct fp_ext *src)
 163 {
 164         uprint("flog10\n");
 165 
 166         fp_monadic_check(dest, src);
 167 
 168         return dest;
 169 }
 170 
 171 struct fp_ext *
 172 fp_flog2(struct fp_ext *dest, struct fp_ext *src)
 173 {
 174         uprint("flog2\n");
 175 
 176         fp_monadic_check(dest, src);
 177 
 178         return dest;
 179 }
 180 
 181 struct fp_ext *
 182 fp_fgetexp(struct fp_ext *dest, struct fp_ext *src)
 183 {
 184         dprint(PINSTR, "fgetexp\n");
 185 
 186         fp_monadic_check(dest, src);
 187 
 188         if (IS_INF(dest)) {
 189                 fp_set_nan(dest);
 190                 return dest;
 191         }
 192         if (IS_ZERO(dest))
 193                 return dest;
 194 
 195         fp_conv_long2ext(dest, (int)dest->exp - 0x3FFF);
 196 
 197         fp_normalize_ext(dest);
 198 
 199         return dest;
 200 }
 201 
 202 struct fp_ext *
 203 fp_fgetman(struct fp_ext *dest, struct fp_ext *src)
 204 {
 205         dprint(PINSTR, "fgetman\n");
 206 
 207         fp_monadic_check(dest, src);
 208 
 209         if (IS_ZERO(dest))
 210                 return dest;
 211 
 212         if (IS_INF(dest))
 213                 return dest;
 214 
 215         dest->exp = 0x3FFF;
 216 
 217         return dest;
 218 }
 219