sysdeps/x86_64/fpu/e_logl.S

   1 /*
   2  * Written by J.T. Conklin <jtc@netbsd.org>.
   3  * Public domain.
   4  *
   5  * Adapted for `long double' by Ulrich Drepper <drepper@cygnus.com>.
   6  * Adapted for x86-64 by Andreas Jaeger <aj@suse.de>.
   7  */
   8
   9 #include <machine/asm.h>
  10
  11
  12         .section .rodata.cst8,"aM",@progbits,8
  13
  14         .p2align 3
  15         .type one,@object
  16 one:    .double 1.0
  17         ASM_SIZE_DIRECTIVE(one)
  18         /* It is not important that this constant is precise.  It is only
  19            a value which is known to be on the safe side for using the
  20            fyl2xp1 instruction.  */
  21         .type limit,@object
  22 limit:  .double 0.29
  23         ASM_SIZE_DIRECTIVE(limit)
  24
  25
  26 #ifdef PIC
  27 # define MO(op) op##(%rip)
  28 #else
  29 # define MO(op) op
  30 #endif
  31
  32         .text
  33 ENTRY(__ieee754_logl)
  34         fldln2                  // log(2)
  35         fldt    8(%rsp)         // x : log(2)
  36         fxam
  37         fnstsw
  38         fld     %st             // x : x : log(2)
  39         testb   $1, %ah
  40         jnz     3f              // in case x is NaN or +-Inf
  41         movzwl  8+8(%rsp), %eax
  42         cmpl    $0xc000, %eax
  43         jae     6f              // x <= -2, avoid overflow from -LDBL_MAX - 1.
  44 4:      fsubl   MO(one)         // x-1 : x : log(2)
  45 6:      fld     %st             // x-1 : x-1 : x : log(2)
  46         fabs                    // |x-1| : x-1 : x : log(2)
  47         fcompl  MO(limit)       // x-1 : x : log(2)
  48         fnstsw                  // x-1 : x : log(2)
  49         andb    $0x45, %ah
  50         jz      2f
  51         fxam
  52         fnstsw
  53         andb    $0x45, %ah
  54         cmpb    $0x40, %ah
  55         jne     5f
  56         fabs                    // log(1) is +0 in all rounding modes.
  57 5:      fstp    %st(1)          // x-1 : log(2)
  58         fyl2xp1                 // log(x)
  59         ret
  60
  61 2:      fstp    %st(0)          // x : log(2)
  62         fyl2x                   // log(x)
  63         ret
  64
  65 3:      testb   $4, %ah
  66         jnz     4b              // in case x is +-Inf
  67         fstp    %st(1)
  68         fstp    %st(1)
  69         fadd    %st(0)
  70         ret
  71 END (__ieee754_logl)
  72
  73
  74 ENTRY(__logl_finite)
  75         fldln2                  // log(2)
  76         fldt    8(%rsp)         // x : log(2)
  77         fld     %st             // x : x : log(2)
  78         fsubl   MO(one)         // x-1 : x : log(2)
  79         fld     %st             // x-1 : x-1 : x : log(2)
  80         fabs                    // |x-1| : x-1 : x : log(2)
  81         fcompl  MO(limit)       // x-1 : x : log(2)
  82         fnstsw                  // x-1 : x : log(2)
  83         andb    $0x45, %ah
  84         jz      2b
  85         fxam
  86         fnstsw
  87         andb    $0x45, %ah
  88         cmpb    $0x40, %ah
  89         jne     7f
  90         fabs                    // log(1) is +0 in all rounding modes.
  91 7:      fstp    %st(1)          // x-1 : log(2)
  92         fyl2xp1                 // log(x)
  93         ret
  94 END (__logl_finite)