Ticket #9962: 0001-libroot.so-update-glibc-s-e_pow.S-on-x86.-Fixes-9962.patch

File 0001-libroot.so-update-glibc-s-e_pow.S-on-x86.-Fixes-9962.patch, 9.9 KB (added by jessicah, 10 years ago)
  • src/system/libroot/posix/glibc/arch/x86/e_pow.S

    From 96b0b1bd078020834648b5ab83b23248e4217ef6 Mon Sep 17 00:00:00 2001
    From: Jessica Hamilton <jessica.l.hamilton@gmail.com>
    Date: Mon, 2 Jun 2014 12:57:35 +1200
    Subject: [PATCH] libroot.so: update glibc's e_pow.S on x86. Fixes #9962
    
    ---
     src/system/libroot/posix/glibc/arch/x86/e_pow.S | 213 +++++++++++++++++-------
     1 file changed, 150 insertions(+), 63 deletions(-)
    
    diff --git a/src/system/libroot/posix/glibc/arch/x86/e_pow.S b/src/system/libroot/posix/glibc/arch/x86/e_pow.S
    index 997dd30..76be3d1 100644
    a b  
    11/* ix87 specific implementation of pow function.
    2    Copyright (C) 1996, 1997, 1998, 1999, 2001 Free Software Foundation, Inc.
     2   Copyright (C) 1996-2014 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/* ix87 specific implementation of pow function.
     10   Copyright (C) 1996-2014 Free Software Foundation, Inc.
    311   This file is part of the GNU C Library.
    412   Contributed by Ulrich Drepper <drepper@cygnus.com>, 1996.
    513
     
    1422   Lesser General Public License for more details.
    1523
    1624   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.  */
     25   License along with the GNU C Library; if not, see
     26   <http://www.gnu.org/licenses/>.  */
    2027
    2128#include <machine/asm.h>
    2229
    23 #ifdef __ELF__
    24     .section .rodata
    25 #else
    26     .text
    27 #endif
     30    .section .rodata.cst8,"aM",@progbits,8
    2831
    29     .align ALIGNARG(4)
    30     ASM_TYPE_DIRECTIVE(infinity,@object)
     32    .p2align 3
     33    .type one,@object
     34one:    .double 1.0
     35    ASM_SIZE_DIRECTIVE(one)
     36    .type limit,@object
     37limit:  .double 0.29
     38    ASM_SIZE_DIRECTIVE(limit)
     39    .type p63,@object
     40p63:    .byte 0, 0, 0, 0, 0, 0, 0xe0, 0x43
     41    ASM_SIZE_DIRECTIVE(p63)
     42    .type p10,@object
     43p10:    .byte 0, 0, 0, 0, 0, 0, 0x90, 0x40
     44    ASM_SIZE_DIRECTIVE(p10)
     45
     46    .section .rodata.cst16,"aM",@progbits,16
     47
     48    .p2align 3
     49    .type infinity,@object
    3150inf_zero:
    3251infinity:
    3352    .byte 0, 0, 0, 0, 0, 0, 0xf0, 0x7f
    3453    ASM_SIZE_DIRECTIVE(infinity)
    35     ASM_TYPE_DIRECTIVE(zero,@object)
     54    .type zero,@object
    3655zero:   .double 0.0
    3756    ASM_SIZE_DIRECTIVE(zero)
    38     ASM_TYPE_DIRECTIVE(minf_mzero,@object)
     57    .type minf_mzero,@object
    3958minf_mzero:
    4059minfinity:
    4160    .byte 0, 0, 0, 0, 0, 0, 0xf0, 0xff
    4261mzero:
    4362    .byte 0, 0, 0, 0, 0, 0, 0, 0x80
    4463    ASM_SIZE_DIRECTIVE(minf_mzero)
    45     ASM_TYPE_DIRECTIVE(one,@object)
    46 one:    .double 1.0
    47     ASM_SIZE_DIRECTIVE(one)
    48     ASM_TYPE_DIRECTIVE(limit,@object)
    49 limit:  .double 0.29
    50     ASM_SIZE_DIRECTIVE(limit)
    5164
    5265#ifdef PIC
    53 #define MO(op) op##@GOTOFF(%ecx)
    54 #define MOX(op,x,f) op##@GOTOFF(%ecx,x,f)
     66# define MO(op) op##@GOTOFF(%ecx)
     67# define MOX(op,x,f) op##@GOTOFF(%ecx,x,f)
    5568#else
    56 #define MO(op) op
    57 #define MOX(op,x,f) op(,x,f)
     69# define MO(op) op
     70# define MOX(op,x,f) op(,x,f)
    5871#endif
    5972
    6073    .text
    ENTRY(__ieee754_pow)  
    6376    fxam
    6477
    6578#ifdef  PIC
    66     call    1f
    67 1:  popl    %ecx
    68     addl    $_GLOBAL_OFFSET_TABLE_+[.-1b], %ecx
     79    LOAD_PIC_REG (cx)
    6980#endif
    7081
    7182    fnstsw
    ENTRY(__ieee754_pow)  
    7485    cmpb    $0x40, %ah  // is y == 0 ?
    7586    je  11f
    7687
    77     cmpb    $0x05, %ah  // is y == ±inf ?
     88    cmpb    $0x05, %ah  // is y == ±inf ?
    7889    je  12f
    7990
    8091    cmpb    $0x01, %ah  // is y == NaN ?
    ENTRY(__ieee754_pow)  
    89100    movb    %ah, %dh
    90101    andb    $0x45, %ah
    91102    cmpb    $0x40, %ah
    92     je  20f     // x is ±0
     103    je  20f     // x is ±0
    93104
    94105    cmpb    $0x05, %ah
    95     je  15f     // x is ±inf
     106    je  15f     // x is ±inf
     107
     108    cmpb    $0x01, %ah
     109    je  32f     // x is NaN
    96110
    97111    fxch            // y : x
    98112
     113    /* fistpll raises invalid exception for |y| >= 1L<<63.  */
     114    fld %st     // y : y : x
     115    fabs            // |y| : y : x
     116    fcompl  MO(p63)     // y : x
     117    fnstsw
     118    sahf
     119    jnc 2f
     120
    99121    /* First see whether `y' is a natural number.  In this case we
    100122       can use a more precise algorithm.  */
    101123    fld %st     // y : y : x
    ENTRY(__ieee754_pow)  
    104126    fucomp  %st(1)      // y : x
    105127    fnstsw
    106128    sahf
    107     jne 2f
     129    jne 3f
    108130
    109     /* OK, we have an integer value for y.  */
     131    /* OK, we have an integer value for y.  If large enough that
     132       errors may propagate out of the 11 bits excess precision, use
     133       the algorithm for real exponent instead.  */
     134    fld %st     // y : y : x
     135    fabs            // |y| : y : x
     136    fcompl  MO(p10)     // y : x
     137    fnstsw
     138    sahf
     139    jnc 2f
    110140    popl    %eax
    111141    popl    %edx
    112142    orl $0, %edx
    ENTRY(__ieee754_pow)  
    132162    fstp    %st(0)      // ST*x
    133163    ret
    134164
    135     /* y is ±NAN */
     165    /* y is ±NAN */
    13616630: fldl    4(%esp)     // x : y
    137167    fldl    MO(one)     // 1.0 : x : y
    138168    fucomp  %st(1)      // x : y
    ENTRY(__ieee754_pow)  
    14317331: fstp    %st(1)
    144174    ret
    145175
     17632: addl    $8, %esp
     177    fstp    %st(1)
     178    ret
     179
     180    .align ALIGNARG(4)
     1812:  // y is a large integer (absolute value at least 1L<<10), but
     182    // may be odd unless at least 1L<<64.  So it may be necessary
     183    // to adjust the sign of a negative result afterwards.
     184    fxch            // x : y
     185    fabs            // |x| : y
     186    fxch            // y : x
    146187    .align ALIGNARG(4)
    147 2:  /* y is a real number.  */
     1883:  /* y is a real number.  */
    148189    fxch            // x : y
    149190    fldl    MO(one)     // 1.0 : x : y
    150     fld %st(1)      // x : 1.0 : x : y
    151     fsub    %st(1)      // x-1 : 1.0 : x : y
    152     fabs            // |x-1| : 1.0 : x : y
    153     fcompl  MO(limit)   // 1.0 : x : y
     191    fldl    MO(limit)   // 0.29 : 1.0 : x : y
     192    fld %st(2)      // x : 0.29 : 1.0 : x : y
     193    fsub    %st(2)      // x-1 : 0.29 : 1.0 : x : y
     194    fabs            // |x-1| : 0.29 : 1.0 : x : y
     195    fucompp         // 1.0 : x : y
    154196    fnstsw
    155197    fxch            // x : 1.0 : y
    156198    sahf
    ENTRY(__ieee754_pow)  
    168210    f2xm1           // 2^fract(y*log2(x))-1 : int(y*log2(x))
    169211    faddl   MO(one)     // 2^fract(y*log2(x)) : int(y*log2(x))
    170212    fscale          // 2^fract(y*log2(x))*2^int(y*log2(x)) : int(y*log2(x))
    171     addl    $8, %esp
    172213    fstp    %st(1)      // 2^fract(y*log2(x))*2^int(y*log2(x))
     214    testb   $2, %dh
     215    jz  292f
     216    // x is negative.  If y is an odd integer, negate the result.
     217    fldl    20(%esp)    // y : abs(result)
     218    fld %st     // y : y : abs(result)
     219    fabs            // |y| : y : abs(result)
     220    fcompl  MO(p63)     // y : abs(result)
     221    fnstsw
     222    sahf
     223    jnc 291f
     224
     225    // We must find out whether y is an odd integer.
     226    fld %st     // y : y : abs(result)
     227    fistpll (%esp)      // y : abs(result)
     228    fildll  (%esp)      // int(y) : y : abs(result)
     229    fucompp         // abs(result)
     230    fnstsw
     231    sahf
     232    jne 292f
     233
     234    // OK, the value is an integer, but is it odd?
     235    popl    %eax
     236    popl    %edx
     237    andb    $1, %al
     238    jz  290f        // jump if not odd
     239    // It's an odd integer.
     240    fchs
     241290:    ret
     242291:    fstp    %st(0)      // abs(result)
     243292:    addl    $8, %esp
    173244    ret
    174245
    175246
    176     // pow(x,±0) = 1
     247    // pow(x,±0) = 1
    177248    .align ALIGNARG(4)
    17824911: fstp    %st(0)      // pop y
    179250    fldl    MO(one)
    180251    ret
    181252
    182     // y == ±inf
     253    // y == ±inf
    183254    .align ALIGNARG(4)
    18425512: fstp    %st(0)      // pop y
    185     fldl    4(%esp)     // x
    186     fabs
    187     fcompl  MO(one)     // < 1, == 1, or > 1
     256    fldl    MO(one)     // 1
     257    fldl    4(%esp)     // x : 1
     258    fabs            // abs(x) : 1
     259    fucompp         // < 1, == 1, or > 1
    188260    fnstsw
    189261    andb    $0x45, %ah
    190262    cmpb    $0x45, %ah
    ENTRY(__ieee754_pow)  
    208280    ret
    209281
    210282    .align ALIGNARG(4)
    211     // x is ±inf
     283    // x is ±inf
    21228415: fstp    %st(0)      // y
    213285    testb   $2, %dh
    214286    jz  16f     // jump if x == +inf
    215287
     288    // fistpll raises invalid exception for |y| >= 1L<<63, so test
     289    // that (in which case y is certainly even) before testing
     290    // whether y is odd.
     291    fld %st     // y : y
     292    fabs            // |y| : y
     293    fcompl  MO(p63)     // y
     294    fnstsw
     295    sahf
     296    jnc 16f
     297
    216298    // We must find out whether y is an odd integer.
    217299    fld %st     // y : y
    218300    fistpll (%esp)      // y
    ENTRY(__ieee754_pow)  
    222304    sahf
    223305    jne 17f
    224306
    225     // OK, the value is an integer, but is the number of bits small
    226     // enough so that all are coming from the mantissa?
     307    // OK, the value is an integer.
    227308    popl    %eax
    228309    popl    %edx
    229310    andb    $1, %al
    230311    jz  18f     // jump if not odd
    231     movl    %edx, %eax
    232     orl %edx, %edx
    233     jns 155f
    234     negl    %eax
    235 155:    cmpl    $0x00200000, %eax
    236     ja  18f     // does not fit in mantissa bits
    237312    // It's an odd integer.
    238313    shrl    $31, %edx
    239314    fldl    MOX(minf_mzero, %edx, 8)
    ENTRY(__ieee754_pow)  
    256331    ret
    257332
    258333    .align ALIGNARG(4)
    259     // x is ±0
     334    // x is ±0
    26033520: fstp    %st(0)      // y
    261336    testb   $2, %dl
    262337    jz  21f     // y > 0
    263338
    264     // x is ±0 and y is < 0.  We must find out whether y is an odd integer.
     339    // x is ±0 and y is < 0.  We must find out whether y is an odd integer.
    265340    testb   $2, %dh
    266341    jz  25f
    267342
     343    // fistpll raises invalid exception for |y| >= 1L<<63, so test
     344    // that (in which case y is certainly even) before testing
     345    // whether y is odd.
     346    fld %st     // y : y
     347    fabs            // |y| : y
     348    fcompl  MO(p63)     // y
     349    fnstsw
     350    sahf
     351    jnc 25f
     352
    268353    fld %st     // y : y
    269354    fistpll (%esp)      // y
    270355    fildll  (%esp)      // int(y) : y
    ENTRY(__ieee754_pow)  
    273358    sahf
    274359    jne 26f
    275360
    276     // OK, the value is an integer, but is the number of bits small
    277     // enough so that all are coming from the mantissa?
     361    // OK, the value is an integer.
    278362    popl    %eax
    279363    popl    %edx
    280364    andb    $1, %al
    281365    jz  27f     // jump if not odd
    282     cmpl    $0xffe00000, %edx
    283     jbe 27f     // does not fit in mantissa bits
    284366    // It's an odd integer.
    285367    // Raise divide-by-zero exception and get minus infinity value.
    286368    fldl    MO(one)
    ENTRY(__ieee754_pow)  
    296378    ret
    297379
    298380    .align ALIGNARG(4)
    299     // x is ±0 and y is > 0.  We must find out whether y is an odd integer.
     381    // x is ±0 and y is > 0.  We must find out whether y is an odd integer.
    30038221: testb   $2, %dh
    301383    jz  22f
    302384
     385    // fistpll raises invalid exception for |y| >= 1L<<63, so test
     386    // that (in which case y is certainly even) before testing
     387    // whether y is odd.
     388    fcoml   MO(p63)     // y
     389    fnstsw
     390    sahf
     391    jnc 22f
     392
    303393    fld %st     // y : y
    304394    fistpll (%esp)      // y
    305395    fildll  (%esp)      // int(y) : y
    ENTRY(__ieee754_pow)  
    308398    sahf
    309399    jne 23f
    310400
    311     // OK, the value is an integer, but is the number of bits small
    312     // enough so that all are coming from the mantissa?
     401    // OK, the value is an integer.
    313402    popl    %eax
    314403    popl    %edx
    315404    andb    $1, %al
    316405    jz  24f     // jump if not odd
    317     cmpl    $0xffe00000, %edx
    318     jae 24f     // does not fit in mantissa bits
    319406    // It's an odd integer.
    320407    fldl    MO(mzero)
    321408    ret