sysdeps/i386/fpu/s_asinhl.S

   1 /* ix87 specific implementation of arcsinh.
   2    Copyright (C) 1996, 1997 Free Software Foundation, Inc.
   3    This file is part of the GNU C Library.
   4    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1996.
   5
   6    The GNU C Library is free software; you can redistribute it and/or
   7    modify it under the terms of the GNU Lesser General Public
   8    License as published by the Free Software Foundation; either
   9    version 2.1 of the License, or (at your option) any later version.
  10
  11    The GNU C Library is distributed in the hope that it will be useful,
  12    but WITHOUT ANY WARRANTY; without even the implied warranty of
  13    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14    Lesser General Public License for more details.
  15
  16    You should have received a copy of the GNU Lesser General Public
  17    License along with the GNU C Library; if not, write to the Free
  18    Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
  19    02111-1307 USA.  */
  20
  21 #include <machine/asm.h>
  22
  23 #ifdef __ELF__
  24         .section .rodata
  25 #else
  26         .text
  27 #endif
  28
  29         .align ALIGNARG(4)
  30         ASM_TYPE_DIRECTIVE(huge,@object)
  31 huge:   .tfloat 1e+4930
  32         ASM_SIZE_DIRECTIVE(huge)
  33         .align ALIGNARG(4)
  34         /* Please note that we use double value for 1.0.  This number
  35            has an exact representation and so we don't get accuracy
  36            problems.  The advantage is that the code is simpler.  */
  37         ASM_TYPE_DIRECTIVE(one,@object)
  38 one:    .double 1.0
  39         ASM_SIZE_DIRECTIVE(one)
  40         /* It is not important that this constant is precise.  It is only
  41            a value which is known to be on the safe side for using the
  42            fyl2xp1 instruction.  */
  43         ASM_TYPE_DIRECTIVE(limit,@object)
  44 limit:  .double 0.29
  45         ASM_SIZE_DIRECTIVE(limit)
  46
  47 #ifdef PIC
  48 #define MO(op) op##@GOTOFF(%edx)
  49 #else
  50 #define MO(op) op
  51 #endif
  52
  53         .text
  54 ENTRY(__asinhl)
  55         movl    12(%esp), %ecx
  56         movl    $0x7fff, %eax
  57         andl    %ecx, %eax
  58         andl    $0x8000, %ecx
  59         movl    %eax, %edx
  60         orl     $0xffff8000, %edx
  61         incl    %edx
  62         jz      7f                      // x in ±Inf or NaN
  63         xorl    %ecx, 12(%esp)
  64         fldt    4(%esp)                 // |x|
  65         cmpl    $0x3fde, %eax
  66         jb      2f                      // |x| < 2^-34
  67         fldln2                          // log(2) : |x|
  68         cmpl    $0x4020, %eax
  69         fxch                            // |x| : log(2)
  70         ja      3f                      // |x| > 2^34
  71 #ifdef  PIC
  72         call    1f
  73 1:      popl    %edx
  74         addl    $_GLOBAL_OFFSET_TABLE_+[.-1b], %edx
  75 #endif
  76         cmpl    $0x4000, %eax
  77         ja      5f                      // |x| > 2
  78
  79         // 2^-34 <= |x| <= 2 => y = sign(x)*log1p(|x|+|x|^2/(1+sqrt(1+|x|^2)))
  80         fld     %st                     // |x| : |x| : log(2)
  81         fmul    %st(1)                  // |x|^2 : |x| : log(2)
  82         fld     %st                     // |x|^2 : |x|^2 : |x| : log(2)
  83         faddl   MO(one)                 // 1+|x|^2 : |x|^2 : |x| : log(2)
  84         fsqrt                           // sqrt(1+|x|^2) : |x|^2 : |x| : log(2)
  85         faddl   MO(one)                 // 1+sqrt(1+|x|^2) : |x|^2 : |x| : log(2)
  86         fdivrp                          // |x|^2/(1+sqrt(1+|x|^2)) : |x| : log(2)
  87         faddp                           // |x|+|x|^2/(1+sqrt(1+|x|^2)) : log(2)
  88         fcoml   MO(limit)
  89         fnstsw
  90         sahf
  91         ja      6f
  92         fyl2xp1
  93         jecxz   4f
  94         fchs
  95 4:      ret
  96
  97 7:      fldt    4(%esp)
  98         ret
  99
 100 6:      faddl   MO(one)
 101         fyl2x
 102         jecxz   4f
 103         fchs
 104 4:      ret
 105
 106         // |x| < 2^-34 => y = x (inexact iff |x| != 0.0)
 107         .align ALIGNARG(4)
 108 2:
 109 #ifdef  PIC
 110         call    1f
 111 1:      popl    %edx
 112         addl    $_GLOBAL_OFFSET_TABLE_+[.-1b], %edx
 113 #endif
 114         jecxz   4f
 115         fchs                            // x
 116 4:      fld     %st                     // x : x
 117         fldt    MO(huge)                // huge : x : x
 118         faddp                           // huge+x : x
 119         fstp    %st(0)                  // x
 120         ret
 121
 122         // |x| > 2^34 => y = sign(x) * (log(|x|) + log(2))
 123         .align ALIGNARG(4)
 124 3:      fyl2x                           // log(|x|)
 125         fldln2                          // log(2) : log(|x|)
 126         faddp                           // log(|x|)+log(2)
 127         jecxz   4f
 128         fchs
 129 4:      ret
 130
 131         // |x| > 2 => y = sign(x) * log(2*|x| + 1/(|x|+sqrt(x*x+1)))
 132         .align ALIGNARG(4)
 133 5:      fld     %st                     // |x| : |x| : log(2)
 134         fadd    %st, %st(1)             // |x| : 2*|x| : log(2)
 135         fld     %st                     // |x| : |x| : 2*|x| : log(2)
 136         fmul    %st(1)                  // |x|^2 : |x| : 2*|x| : log(2)
 137         faddl   MO(one)                 // 1+|x|^2 : |x| : 2*|x| : log(2)
 138         fsqrt                           // sqrt(1+|x|^2) : |x| : 2*|x| : log(2)
 139         faddp                           // |x|+sqrt(1+|x|^2) : 2*|x| : log(2)
 140         fdivrl  MO(one)                 // 1/(|x|+sqrt(1+|x|^2)) : 2*|x| : log(2)
 141         faddp                           // 2*|x|+1/(|x|+sqrt(1+|x|^2)) : log(2)
 142         fyl2x                           // log(2*|x|+1/(|x|+sqrt(1+|x|^2)))
 143         jecxz   4f
 144         fchs
 145 4:      ret
 146 END(__asinhl)
 147 weak_alias (__asinhl, asinhl)