Return to cbrt.s CVS log

Up to [CSRG BSD Unix] / 43BSDTahoe / usr.lib / libm / tahoe

Annotation of 43BSDTahoe/usr.lib/libm/tahoe/cbrt.s, revision 1.1

1.1     ! root        1: # Copyright (c) 1987 Regents of the University of California.
        !             2: # All rights reserved.
        !             3: #
        !             4: # Redistribution and use in source and binary forms are permitted
        !             5: # provided that the above copyright notice and this paragraph are
        !             6: # duplicated in all such forms and that any documentation,
        !             7: # advertising materials, and other materials related to such
        !             8: # distribution and use acknowledge that the software was developed
        !             9: # by the University of California, Berkeley.  The name of the
        !            10: # University may not be used to endorse or promote products derived
        !            11: # from this software without specific prior written permission.
        !            12: # THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR
        !            13: # IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED
        !            14: # WARRANTIES OF MERCHANTIBILITY AND FITNESS FOR A PARTICULAR PURPOSE.
        !            15: #
        !            16: # All recipients should regard themselves as participants in an ongoing
        !            17: # research project and hence should feel obligated to report their
        !            18: # experiences (good or bad) with these elementary function codes, using
        !            19: # the sendbug(8) program, to the authors.
        !            20: #
        !            21: #      @(#)cbrt.s      5.4 (Berkeley) 6/30/88
        !            22: #
        !            23:        .data
        !            24:        .align  2
        !            25: _sccsid:
        !            26: .asciz "@(#)cbrt.s     5.4     (ucb.elefunt) 6/30/88"
        !            27: 
        !            28: # double cbrt(double arg)
        !            29: # W. Kahan, 10/13/80. revised 1/13/84 for keeping sign symmetry
        !            30: # Re-coded in tahoe assembly language by Z. Alex Liu (7/13/87)
        !            31: # Max error less than 0.667 ULPs _if_ +,-,*,/ were all correctly rounded...
        !            32:        .globl  _cbrt
        !            33:        .globl  _d_cbrt
        !            34:        .globl  _dcbrt_
        !            35:        .text
        !            36:        .align  2
        !            37: _cbrt:
        !            38: _d_cbrt:
        !            39:        .word   0x01fc          # save r2-r8
        !            40:        movl    4(fp),r0        # r0:r1 = x
        !            41:        movl    8(fp),r1
        !            42:        brb     1f
        !            43: _dcbrt_:
        !            44:        .word   0x01fc          # save r2-r8
        !            45:        movl    4(fp),r8
        !            46:        movl    (r8),r0
        !            47:        movl    4(r8),r1        # r0:r1 = x
        !            48: 
        !            49: 1:     andl3   $0x7f800000,r0,r2       # biased exponent of x
        !            50:        beql    return          # dcbrt(0)=0  dcbrt(res)=res. operand
        !            51:        andl3   $0x80000000,r0,r8       # r8 has sign(x)
        !            52:        xorl2   r8,r0           # r0 is abs(x)
        !            53:        movl    r0,r2           # r2 has abs(x)
        !            54:        divl2   $3,r2           # rough dcbrt with bias/3
        !            55:        addl2   B,r2            # restore bias, diminish fraction
        !            56:        ldf     r2              # acc = |q|=|dcbrt| to 5 bits
        !            57:        mulf    r2              # acc = qq
        !            58:        divf    r0              # acc = qq/|x|
        !            59:        mulf    r2              # acc = qqq/|x|
        !            60:        addf    C               # acc = C+qqq/|x|
        !            61:        stf     r3              # r3 = s = C+qqq/|x|
        !            62:        ldf     D               # acc = D
        !            63:        divf    r3              # acc = D/s
        !            64:        addf    E               # acc = E+D/s
        !            65:        addf    r3              # acc = s+E+D/s
        !            66:        stf     r3              # r3 = s+E+D/s
        !            67:        ldf     F               # acc = F
        !            68:        divf    r3              # acc = F/(s+E+D/s)
        !            69:        addf    G               # acc = G+F/(s+E+D/s)
        !            70:        mulf    r2              # acc = q*(G+F/(s+E+D/s)) = new q to 23 bits
        !            71:        stf     r2              # r2 = q*(G+F/(s+E+D/s)) = new q to 23 bits
        !            72:        clrl    r3              # r2:r3 = q as double float
        !            73:        ldd     r2              # acc = q as double float
        !            74:        muld    r2              # acc = qq exactly
        !            75:        std     r4              # r4:r5 = qq exactly
        !            76:        ldd     r0              # acc = |x|
        !            77:        divd    r4              # acc = |x|/(q*q) rounded
        !            78:        std     r0              # r0:r1 = |x|/(q*q) rounded
        !            79:        subd    r2              # acc = |x|/(q*q)-q exactly
        !            80:        std     r6              # r6:r7 = |x|/(q*q)-q exactly
        !            81:        movl    r2,r4
        !            82:        clrl    r5              # r4:r5 = q as double float
        !            83:        addl2   $0x800000,r4    # r4:r5 = 2*q
        !            84:        ldd     r4              # acc = 2*q
        !            85:        addd    r0              # acc = 2*q+|x|/(q*q)
        !            86:        std     r4              # r4:r5 = 2*q+|x|/(q*q)
        !            87:        ldd     r6              # acc = |x|/(q*q)-q
        !            88:        divd    r4              # acc = (|x|/(q*q)-q)/(2*q+|x|/(q*q))
        !            89:        muld    r2              # acc = q*(|x|/(q*q)-q)/(2*q+|x|/(q*q))
        !            90:        addd    r2              # acc = q+q*(|x|/(q*q)-q)/(2*q+|x|/(q*q))
        !            91:        std     r0              # r0:r1 = |result|
        !            92:        orl2    r8,r0           # restore the sign bit
        !            93: return: ret                    # error less than 0.667ULPs?
        !            94: 
        !            95:        .data
        !            96:        .align  2
        !            97: B :    .long   721142941       #(86-0.03306235651)*(2^23)
        !            98:        .align  2
        !            99: C:     .long   0x400af8b0      #.float 0f0.5428571429  # 19/35
        !           100:        .align  2
        !           101: D:     .long   0xc0348ef1      #.float 0f-0.7053061224 # -864/1225
        !           102:        .align  2
        !           103: E:     .long   0x40b50750      #.float 0f1.414285714   # 99/70
        !           104:        .align  2
        !           105: F:     .long   0x40cdb6db      #.float 0f1.607142857   # 45/28
        !           106:        .align  2
        !           107: G:     .long   0x3fb6db6e      #.float 0f0.3571428571  # 5/14

This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.