1# 2# Copyright (c) 1980 Regents of the University of California. 3# All rights reserved. The Berkeley software License Agreement 4# specifies the terms and conditions for redistribution. 5# 6# @(#)cbrt.s 5.1 (Berkeley) 05/08/85 7# 8# 9# double cbrt(arg) 10# double arg 11# no error exits 12#method: range reduction to [1/8,1], poly appox, newtons method 13# J F Jarvis, August 10,1978 14.globl _cbrt 15.text 16.align 1 17_cbrt: 18 .word 0x00c0 19 bispsw $0xe0 20 clrl r3 21 movd 4(ap),r4 22 jgtr range 23 jeql retz 24 mnegd r4,r4 # |arg| in r0,r1 25 movl $0x100,r3 # sign bit of result 26range: 27 extzv $7,$8,r4,r6 28 insv $128,$7,$8,r4 # 0.5<= frac: r4,r5 <1.0 29 clrl r7 30 ediv $3,r6,r6,r7 # r6= expnt/3; r7= expnt%3 31 addb2 $86,r6 32 bisl2 r3,r6 # sign,exponent of result 33 polyf r4,$3,pcoef # initial estimate is Hart&Cheney CBRT 0642 34 # D=4.1 35 muld3 r0,r0,r2 # Newtons method, iteration 1, H&C 6.1.10 36 divd3 r2,r4,r2 37 subd3 r2,r0,r2 38 muld2 third,r2 39 subd2 r2,r0 # D=8.2 40 muld3 r0,r0,r2 # iteration 2 41 divd3 r2,r4,r2 42 subd3 r2,r0,r2 43 muld2 third,r2 44 subd2 r2,r0 45 muld2 hc[r7],r0 # set range 46 insv r6,$7,$9,r0 # set sign,exponent 47 ret 48retz: 49 clrd r0 50 ret 51.data 52.align 2 53third: .double 0d0.33333333333333333333e+0 54hc: 55 .double 0d1.25992104989487316476e+0 56 .double 0d1.58740105196819947475e+0 57 .double 0d1.0e+0 58pcoef: 59 .float 0f0.1467073818e+0 60 .float 0f-0.5173964673e+0 61 .float 0f0.9319858515e+0 62 .float 0f0.4387762363e+0 63