sysdeps/i386/fpu/s_asinh.S

   1 /* ix87 specific implementation of arcsinh.
   2    Copyright (C) 1996-2019 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, see
  18    <https://www.gnu.org/licenses/>.  */
  19
  20 #include <machine/asm.h>
  21 #include <libm-alias-double.h>
  22
  23         .section .rodata
  24
  25         .align ALIGNARG(4)
  26         .type huge,@object
  27 huge:   .double 1e+300
  28         ASM_SIZE_DIRECTIVE(huge)
  29         .type one,@object
  30 one:    .double 1.0
  31         ASM_SIZE_DIRECTIVE(one)
  32         .type limit,@object
  33 limit:  .double 0.29
  34         ASM_SIZE_DIRECTIVE(limit)
  35
  36 #ifdef PIC
  37 #define MO(op) op##@GOTOFF(%edx)
  38 #else
  39 #define MO(op) op
  40 #endif
  41
  42         .text
  43 ENTRY(__asinh)
  44         movl    8(%esp), %ecx
  45         movl    $0x7fffffff, %eax
  46         andl    %ecx, %eax
  47         andl    $0x80000000, %ecx
  48         movl    %eax, %edx
  49         orl     $0x800fffff, %edx
  50         incl    %edx
  51         jz      7f                      // x in ±Inf or NaN
  52         xorl    %ecx, 8(%esp)
  53         fldl    4(%esp)                 // |x|
  54         cmpl    $0x3e300000, %eax
  55         jb      2f                      // |x| < 2^-28
  56         fldln2                          // log(2) : |x|
  57         cmpl    $0x41b00000, %eax
  58         fxch                            // |x| : log(2)
  59         ja      3f                      // |x| > 2^28
  60 #ifdef  PIC
  61         LOAD_PIC_REG (dx)
  62 #endif
  63         cmpl    $0x40000000, %eax
  64         ja      5f                      // |x| > 2
  65
  66         // 2^-28 <= |x| <= 2 => y = sign(x)*log1p(|x|+|x|^2/(1+sqrt(1+|x|^2)))
  67         fld     %st                     // |x| : |x| : log(2)
  68         fmul    %st(1)                  // |x|^2 : |x| : log(2)
  69         fld     %st                     // |x|^2 : |x|^2 : |x| : log(2)
  70         faddl   MO(one)                 // 1+|x|^2 : |x|^2 : |x| : log(2)
  71         fsqrt                           // sqrt(1+|x|^2) : |x|^2 : |x| : log(2)
  72         faddl   MO(one)                 // 1+sqrt(1+|x|^2) : |x|^2 : |x| : log(2)
  73         fdivrp                          // |x|^2/(1+sqrt(1+|x|^2)) : |x| : log(2)
  74         faddp                           // |x|+|x|^2/(1+sqrt(1+|x|^2)) : log(2)
  75         fcoml   MO(limit)
  76         fnstsw
  77         sahf
  78         ja      6f
  79         fyl2xp1
  80         jecxz   4f
  81         fchs
  82 4:      ret
  83
  84 7:      fldl    4(%esp)
  85         ret
  86
  87 6:      faddl   MO(one)
  88         fyl2x
  89         jecxz   4f
  90         fchs
  91 4:      ret
  92
  93         // |x| < 2^-28 => y = x (inexact iff |x| != 0.0)
  94         .align ALIGNARG(4)
  95 2:
  96 #ifdef  PIC
  97         LOAD_PIC_REG (dx)
  98 #endif
  99         jecxz   4f
 100         fchs                            // x
 101 4:      fld     %st                     // x : x
 102         faddl   MO(huge)                // huge+x : x
 103         fstp    %st(0)                  // x
 104         cmpl    $0x00100000, %eax
 105         jae     8f
 106         subl    $8, %esp
 107         cfi_adjust_cfa_offset (8)
 108         fld     %st(0)
 109         fmul    %st(0)
 110         fstpl   (%esp)
 111         addl    $8, %esp
 112         cfi_adjust_cfa_offset (-8)
 113 8:      ret
 114
 115         // |x| > 2^28 => y = sign(x) * (log(|x|) + log(2))
 116         .align ALIGNARG(4)
 117 3:      fyl2x                           // log(|x|)
 118         fldln2                          // log(2) : log(|x|)
 119         faddp                           // log(|x|)+log(2)
 120         jecxz   4f
 121         fchs
 122 4:      ret
 123
 124         // |x| > 2 => y = sign(x) * log(2*|x| + 1/(|x|+sqrt(x*x+1)))
 125         .align ALIGNARG(4)
 126 5:      fld     %st                     // |x| : |x| : log(2)
 127         fadd    %st, %st(1)             // |x| : 2*|x| : log(2)
 128         fld     %st                     // |x| : |x| : 2*|x| : log(2)
 129         fmul    %st(1)                  // |x|^2 : |x| : 2*|x| : log(2)
 130         faddl   MO(one)                 // 1+|x|^2 : |x| : 2*|x| : log(2)
 131         fsqrt                           // sqrt(1+|x|^2) : |x| : 2*|x| : log(2)
 132         faddp                           // |x|+sqrt(1+|x|^2) : 2*|x| : log(2)
 133         fdivrl  MO(one)                 // 1/(|x|+sqrt(1+|x|^2)) : 2*|x| : log(2)
 134         faddp                           // 2*|x|+1/(|x|+sqrt(1+|x|^2)) : log(2)
 135         fyl2x                           // log(2*|x|+1/(|x|+sqrt(1+|x|^2)))
 136         jecxz   4f
 137         fchs
 138 4:      ret
 139 END(__asinh)
 140 libm_alias_double (__asinh, asinh)