arch/m68k-all/m680x0/fpsp/slog2.sa

   1 *       $NetBSD: slog2.sa,v 1.2 1994/10/26 07:49:52 cgd Exp $
   2
   3 *       MOTOROLA MICROPROCESSOR & MEMORY TECHNOLOGY GROUP
   4 *       M68000 Hi-Performance Microprocessor Division
   5 *       M68040 Software Package
   6 *
   7 *       M68040 Software Package Copyright (c) 1993, 1994 Motorola Inc.
   8 *       All rights reserved.
   9 *
  10 *       THE SOFTWARE is provided on an "AS IS" basis and without warranty.
  11 *       To the maximum extent permitted by applicable law,
  12 *       MOTOROLA DISCLAIMS ALL WARRANTIES WHETHER EXPRESS OR IMPLIED,
  13 *       INCLUDING IMPLIED WARRANTIES OF MERCHANTABILITY OR FITNESS FOR A
  14 *       PARTICULAR PURPOSE and any warranty against infringement with
  15 *       regard to the SOFTWARE (INCLUDING ANY MODIFIED VERSIONS THEREOF)
  16 *       and any accompanying written materials.
  17 *
  18 *       To the maximum extent permitted by applicable law,
  19 *       IN NO EVENT SHALL MOTOROLA BE LIABLE FOR ANY DAMAGES WHATSOEVER
  20 *       (INCLUDING WITHOUT LIMITATION, DAMAGES FOR LOSS OF BUSINESS
  21 *       PROFITS, BUSINESS INTERRUPTION, LOSS OF BUSINESS INFORMATION, OR
  22 *       OTHER PECUNIARY LOSS) ARISING OF THE USE OR INABILITY TO USE THE
  23 *       SOFTWARE.  Motorola assumes no responsibility for the maintenance
  24 *       and support of the SOFTWARE.
  25 *
  26 *       You are hereby granted a copyright license to use, modify, and
  27 *       distribute the SOFTWARE so long as this entire notice is retained
  28 *       without alteration in any modified and/or redistributed versions,
  29 *       and that such modified versions are clearly identified as such.
  30 *       No licenses are granted by implication, estoppel or otherwise
  31 *       under any patents or trademarks of Motorola, Inc.
  32
  33 *
  34 *       slog2.sa 3.1 12/10/90
  35 *
  36 *       The entry point slog10 computes the base-10
  37 *       logarithm of an input argument X.
  38 *       slog10d does the same except the input value is a
  39 *       denormalized number.
  40 *       sLog2 and sLog2d are the base-2 analogues.
  41 *
  42 *       INPUT:  Double-extended value in memory location pointed to
  43 *               by address register a0.
  44 *
  45 *       OUTPUT: log_10(X) or log_2(X) returned in floating-point
  46 *               register fp0.
  47 *
  48 *       ACCURACY and MONOTONICITY: The returned result is within 1.7
  49 *               ulps in 64 significant bit, i.e. within 0.5003 ulp
  50 *               to 53 bits if the result is subsequently rounded
  51 *               to double precision. The result is provably monotonic
  52 *               in double precision.
  53 *
  54 *       SPEED:  Two timings are measured, both in the copy-back mode.
  55 *               The first one is measured when the function is invoked
  56 *               the first time (so the instructions and data are not
  57 *               in cache), and the second one is measured when the
  58 *               function is reinvoked at the same input argument.
  59 *
  60 *       ALGORITHM and IMPLEMENTATION NOTES:
  61 *
  62 *       slog10d:
  63 *
  64 *       Step 0.   If X < 0, create a NaN and raise the invalid operation
  65 *                 flag. Otherwise, save FPCR in D1; set FpCR to default.
  66 *       Notes:    Default means round-to-nearest mode, no floating-point
  67 *                 traps, and precision control = double extended.
  68 *
  69 *       Step 1.   Call slognd to obtain Y = log(X), the natural log of X.
  70 *       Notes:    Even if X is denormalized, log(X) is always normalized.
  71 *
  72 *       Step 2.   Compute log_10(X) = log(X) * (1/log(10)).
  73 *            2.1  Restore the user FPCR
  74 *            2.2  Return ans := Y * INV_L10.
  75 *
  76 *
  77 *       slog10:
  78 *
  79 *       Step 0.   If X < 0, create a NaN and raise the invalid operation
  80 *                 flag. Otherwise, save FPCR in D1; set FpCR to default.
  81 *       Notes:    Default means round-to-nearest mode, no floating-point
  82 *                 traps, and precision control = double extended.
  83 *
  84 *       Step 1.   Call sLogN to obtain Y = log(X), the natural log of X.
  85 *
  86 *       Step 2.   Compute log_10(X) = log(X) * (1/log(10)).
  87 *            2.1  Restore the user FPCR
  88 *            2.2  Return ans := Y * INV_L10.
  89 *
  90 *
  91 *       sLog2d:
  92 *
  93 *       Step 0.   If X < 0, create a NaN and raise the invalid operation
  94 *                 flag. Otherwise, save FPCR in D1; set FpCR to default.
  95 *       Notes:    Default means round-to-nearest mode, no floating-point
  96 *                 traps, and precision control = double extended.
  97 *
  98 *       Step 1.   Call slognd to obtain Y = log(X), the natural log of X.
  99 *       Notes:    Even if X is denormalized, log(X) is always normalized.
 100 *
 101 *       Step 2.   Compute log_10(X) = log(X) * (1/log(2)).
 102 *            2.1  Restore the user FPCR
 103 *            2.2  Return ans := Y * INV_L2.
 104 *
 105 *
 106 *       sLog2:
 107 *
 108 *       Step 0.   If X < 0, create a NaN and raise the invalid operation
 109 *                 flag. Otherwise, save FPCR in D1; set FpCR to default.
 110 *       Notes:    Default means round-to-nearest mode, no floating-point
 111 *                 traps, and precision control = double extended.
 112 *
 113 *       Step 1.   If X is not an integer power of two, i.e., X != 2^k,
 114 *                 go to Step 3.
 115 *
 116 *       Step 2.   Return k.
 117 *            2.1  Get integer k, X = 2^k.
 118 *            2.2  Restore the user FPCR.
 119 *            2.3  Return ans := convert-to-double-extended(k).
 120 *
 121 *       Step 3.   Call sLogN to obtain Y = log(X), the natural log of X.
 122 *
 123 *       Step 4.   Compute log_2(X) = log(X) * (1/log(2)).
 124 *            4.1  Restore the user FPCR
 125 *            4.2  Return ans := Y * INV_L2.
 126 *
 127
 128 SLOG2    IDNT    2,1 Motorola 040 Floating Point Software Package
 129
 130         section 8
 131
 132         xref    t_frcinx
 133         xref    t_operr
 134         xref    slogn
 135         xref    slognd
 136
 137 INV_L10  DC.L $3FFD0000,$DE5BD8A9,$37287195,$00000000
 138
 139 INV_L2   DC.L $3FFF0000,$B8AA3B29,$5C17F0BC,$00000000
 140
 141         xdef    slog10d
 142 slog10d:
 143 *--entry point for Log10(X), X is denormalized
 144         move.l          (a0),d0
 145         blt.w           invalid
 146         move.l          d1,-(sp)
 147         clr.l           d1
 148         bsr             slognd                  ...log(X), X denorm.
 149         fmove.l         (sp)+,fpcr
 150         fmul.x          INV_L10,fp0
 151         bra             t_frcinx
 152
 153         xdef    slog10
 154 slog10:
 155 *--entry point for Log10(X), X is normalized
 156
 157         move.l          (a0),d0
 158         blt.w           invalid
 159         move.l          d1,-(sp)
 160         clr.l           d1
 161         bsr             slogn                   ...log(X), X normal.
 162         fmove.l         (sp)+,fpcr
 163         fmul.x          INV_L10,fp0
 164         bra             t_frcinx
 165
 166
 167         xdef    slog2d
 168 slog2d:
 169 *--entry point for Log2(X), X is denormalized
 170
 171         move.l          (a0),d0
 172         blt.w           invalid
 173         move.l          d1,-(sp)
 174         clr.l           d1
 175         bsr             slognd                  ...log(X), X denorm.
 176         fmove.l         (sp)+,fpcr
 177         fmul.x          INV_L2,fp0
 178         bra             t_frcinx
 179
 180         xdef    slog2
 181 slog2:
 182 *--entry point for Log2(X), X is normalized
 183         move.l          (a0),d0
 184         blt.w           invalid
 185
 186         move.l          8(a0),d0
 187         bne.b           continue                ...X is not 2^k
 188
 189         move.l          4(a0),d0
 190         and.l           #$7FFFFFFF,d0
 191         tst.l           d0
 192         bne.b           continue
 193
 194 *--X = 2^k.
 195         move.w          (a0),d0
 196         and.l           #$00007FFF,d0
 197         sub.l           #$3FFF,d0
 198         fmove.l         d1,fpcr
 199         fmove.l         d0,fp0
 200         bra             t_frcinx
 201
 202 continue:
 203         move.l          d1,-(sp)
 204         clr.l           d1
 205         bsr             slogn                   ...log(X), X normal.
 206         fmove.l         (sp)+,fpcr
 207         fmul.x          INV_L2,fp0
 208         bra             t_frcinx
 209
 210 invalid:
 211         bra             t_operr
 212
 213         end