Return to sqrt.s CVS log

Up to [CSRG BSD Unix] / 3BSD / lib / libnm

Annotation of 3BSD/lib/libnm/sqrt.s, revision 1.1.1.1

1.1       root        1: # double sqrt(arg)
                      2: # double arg
                      3: # if(arg<0.0) { _errno=EDOM; return(0.0) }
                      4: # J. F. Jarvis August 2, 1978
                      5: .set   EDOM,98
                      6: .text
                      7: .align 1
                      8: .globl _sqrt
                      9: .globl _errno
                     10: _sqrt:
                     11:        .word   0x0c0
                     12:        movd    4(ap),r4
                     13:        jgtr    range
                     14:        jeql    retz
                     15:        movl    $EDOM,_errno
                     16: retz:
                     17:        clrd    r0
                     18:        ret
                     19: #
                     20: range:
                     21:        extzv   $7,$8,r4,r6     # r6 is exponent of arg
                     22:        insv    $128,$7,$8,r4   # r2: 0.5<=fraction(arg)<1.0
                     23:        incl    r6
                     24:        clrl    r7
                     25:        ediv    $2,r6,r6,r7     # r6=(exp+1)/2; r7=(exp+1)%2
                     26:        addb2   $64,r6  # r6 is correct exponent for result
                     27:        polyf   r4,$4,pcoef     # init estimate of sqrt(frac)
                     28:                                                # Hart&Cheney SQRT 0132 D=5.1
                     29:        divd3   r0,r4,r2        # Newtons method, 2 iterations
                     30:        addd2   r2,r0
                     31:        muld2   $0d0.5e+0,r0
                     32:        divd3   r0,r4,r2        # Hart&Cheney 6.1.7
                     33:        addd2   r2,r0   #d=21 at exit
                     34:        muld2   hc[r7],r0       # *sqrt(2) requ for even org exp.
                     35:        insv    r6,$7,$8,r0     # insert correct exp.
                     36:        ret
                     37: .data
                     38: .align 3
                     39: pcoef:
                     40:        .float 0f-0.1214683825e+0
                     41:        .float 0f0.5010420763e+0
                     42:        .float 0f-0.9093210498e+0
                     43:        .float 0f0.1300669049e+1
                     44:        .float 0f0.2290699453e+0
                     45: hc:            .double 0d0.35355339059327376220e+0     # sqrt(2)/4
                     46:                .double 0d0.5e+0

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