Annotation of 43BSD/usr.lib/libm/VAX/cabs.s, revision 1.1

1.1     ! root        1: # 
        !             2: # Copyright (c) 1985 Regents of the University of California.
        !             3: # 
        !             4: # Use and reproduction of this software are granted  in  accordance  with
        !             5: # the terms and conditions specified in  the  Berkeley  Software  License
        !             6: # Agreement (in particular, this entails acknowledgement of the programs'
        !             7: # source, and inclusion of this notice) with the additional understanding
        !             8: # that  all  recipients  should regard themselves as participants  in  an
        !             9: # ongoing  research  project and hence should  feel  obligated  to report
        !            10: # their  experiences (good or bad) with these elementary function  codes,
        !            11: # using "sendbug 4bsd-bugs@BERKELEY", to the authors.
        !            12: #
        !            13: 
        !            14: # @(#)cabs.s   1.2 (Berkeley) 8/21/85
        !            15: 
        !            16: # double precision complex absolute value
        !            17: # CABS by W. Kahan, 9/7/80.
        !            18: # Revised for reserved operands by E. LeBlanc, 8/18/82
        !            19: # argument for complex absolute value by reference, *4(ap)
        !            20: # argument for cabs and hypot (C fcns) by value, 4(ap)
        !            21: # output is in r0:r1 (error less than 0.86 ulps)
        !            22: 
        !            23:        .text
        !            24:        .align  1
        !            25:        .globl  _cabs
        !            26:        .globl  _hypot
        !            27:        .globl  _z_abs
        !            28:        .globl  libm$cdabs_r6
        !            29:        .globl  libm$dsqrt_r5
        !            30: 
        !            31: #      entry for c functions cabs and hypot
        !            32: _cabs:
        !            33: _hypot:
        !            34:        .word   0x807c          # save r2-r6, enable floating overflow
        !            35:        movq    4(ap),r0        # r0:1 = x
        !            36:        movq    12(ap),r2       # r2:3 = y
        !            37:        jmp     cabs2
        !            38: #      entry for Fortran use, call by:   d = abs(z)
        !            39: _z_abs:
        !            40:        .word   0x807c          # save r2-r6, enable floating overflow
        !            41:        movl    4(ap),r2        # indirect addressing is necessary here
        !            42:        movq    (r2)+,r0        # r0:1 = x
        !            43:        movq    (r2),r2         # r2:3 = y
        !            44: 
        !            45: cabs2:
        !            46:        bicw3   $0x7f,r0,r4     # r4 has signed biased exp of x
        !            47:        cmpw    $0x8000,r4
        !            48:        jeql    return          # x is a reserved operand, so return it
        !            49:        bicw3   $0x7f,r2,r5     # r5 has signed biased exp of y
        !            50:        cmpw    $0x8000,r5
        !            51:        jneq    cont            # y isn't a reserved operand
        !            52:        movq    r2,r0           # return y if it's reserved
        !            53:        ret
        !            54: 
        !            55: cont:
        !            56:        bsbb    regs_set        # r0:1 = dsqrt(x^2+y^2)/2^r6
        !            57:        addw2   r6,r0           # unscaled cdabs in r0:1
        !            58:        jvc     return          # unless it overflows
        !            59:        subw2   $0x80,r0        # halve r0 to get meaningful overflow
        !            60:        addd2   r0,r0           # overflow; r0 is half of true abs value
        !            61: return:
        !            62:        ret
        !            63: 
        !            64: libm$cdabs_r6:                 # ENTRY POINT for cdsqrt
        !            65:                                # calculates a scaled (factor in r6)
        !            66:                                # complex absolute value
        !            67: 
        !            68:        movq    (r4)+,r0        # r0:r1 = x via indirect addressing
        !            69:        movq    (r4),r2         # r2:r3 = y via indirect addressing
        !            70: 
        !            71:        bicw3   $0x7f,r0,r5     # r5 has signed biased exp of x
        !            72:        cmpw    $0x8000,r5
        !            73:        jeql    cdreserved      # x is a reserved operand
        !            74:        bicw3   $0x7f,r2,r5     # r5 has signed biased exp of y
        !            75:        cmpw    $0x8000,r5
        !            76:        jneq    regs_set        # y isn't a reserved operand either?
        !            77: 
        !            78: cdreserved:
        !            79:        movl    *4(ap),r4       # r4 -> (u,v), if x or y is reserved
        !            80:        movq    r0,(r4)+        # copy u and v as is and return
        !            81:        movq    r2,(r4)         # (again addressing is indirect)
        !            82:        ret
        !            83: 
        !            84: regs_set:
        !            85:        bicw2   $0x8000,r0      # r0:r1 = dabs(x)
        !            86:        bicw2   $0x8000,r2      # r2:r3 = dabs(y)
        !            87:        cmpw    r0,r2
        !            88:        jgeq    ordered
        !            89:        movq    r0,r4
        !            90:        movq    r2,r0
        !            91:        movq    r4,r2           # force y's exp <= x's exp
        !            92: ordered:
        !            93:        bicw3   $0x7f,r0,r6     # r6 = exponent(x) + bias(129)
        !            94:        jeql    retsb           # if x = y = 0 then cdabs(x,y) = 0
        !            95:        subw2   $0x4780,r6      # r6 = exponent(x) - 14
        !            96:        subw2   r6,r0           # 2^14 <= scaled x < 2^15
        !            97:        bitw    $0xff80,r2
        !            98:        jeql    retsb           # if y = 0 return dabs(x)
        !            99:        subw2   r6,r2
        !           100:        cmpw    $0x3780,r2      # if scaled y < 2^-18
        !           101:        jgtr    retsb           #   return dabs(x)
        !           102:        emodd   r0,$0,r0,r4,r0  # r4 + r0:1 = scaled x^2
        !           103:        emodd   r2,$0,r2,r5,r2  # r5 + r2:3 = scaled y^2
        !           104:        addd2   r2,r0
        !           105:        addl2   r5,r4
        !           106:        cvtld   r4,r2
        !           107:        addd2   r2,r0           # r0:1 = scaled x^2 + y^2
        !           108:        jmp     libm$dsqrt_r5   # r0:1 = dsqrt(x^2+y^2)/2^r6
        !           109: retsb:
        !           110:        rsb                     # error < 0.86 ulp

unix.superglobalmegacorp.com

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