xref: /csrg-svn/lib/libm/tahoe/cabs.s (revision 61310)
1*61310Sbostic# Copyright (c) 1987, 1993
2*61310Sbostic#	The Regents of the University of California.  All rights reserved.
332004Szliu#
445306Sbostic# %sccs.include.redist.sh%
534128Sbostic#
6*61310Sbostic#	@(#)cabs.s	8.1 (Berkeley) 06/04/93
734128Sbostic#
832004Szliu	.data
932004Szliu	.align	2
1032004Szliu_sccsid:
11*61310Sbostic.asciz	"@(#)cabs.s	8.1	8.1	(ucb.elefunt)	06/04/93"
1232004Szliu
1332004Szliu# double precision complex absolute value
1432004Szliu# CABS by W. Kahan, 9/7/80.
1532004Szliu# Revised for reserved operands by E. LeBlanc, 8/18/82
1632004Szliu# argument for complex absolute value by reference, *4(fp)
1732004Szliu# argument for cabs and hypot (C fcns) by value, 4(fp)
1832004Szliu# output is in r0:r1
1932004Szliu
2032004Szliu	.text
2132004Szliu	.align	2
2232004Szliu	.globl  _cabs
2332004Szliu	.globl  _hypot
2432004Szliu	.globl	_z_abs
2532004Szliu
2632004Szliu#	entry for c functions cabs and hypot
2732004Szliu_cabs:
2832004Szliu_hypot:
2932004Szliu	.word	0x807c		# save r2-r6, enable floating overflow
3032004Szliu	movl	16(fp),r3
3132004Szliu	movl	12(fp),r2	# r2:3 = y
3232004Szliu	movl	8(fp),r1
3332004Szliu	movl	4(fp),r0	# r0:1 = x
3432004Szliu	brb	1f
3532004Szliu#	entry for Fortran use, call by:   d = abs(z)
3632004Szliu_z_abs:
3732004Szliu	.word	0x807c		# save r2-r6, enable floating overflow
3832004Szliu	movl	4(fp),r4	# indirect addressing is necessary here
3932004Szliu	movl	12(r4),r3	#
4032004Szliu	movl	8(r4),r2	# r2:3 = y
4132004Szliu	movl	4(r4),r1	#
4232004Szliu	movl	(r4),r0		# r0:1 = x
4332004Szliu1:	andl3	$0xff800000,r0,r4	# r4 has signed biased exp of x
4432004Szliu	cmpl	$0x80000000,r4
4532004Szliu	beql	2f		# x is a reserved operand, so return it
4632004Szliu	andl3	$0xff800000,r2,r5	# r5 has signed biased exp of y
4732004Szliu	cmpl	$0x80000000,r5
4832004Szliu	bneq	3f		# y isn't a reserved operand
4932004Szliu	movl	r3,r1
5032004Szliu	movl	r2,r0		# return y if it's reserved
5132004Szliu2:	ret
5232004Szliu
5332004Szliu3:	callf	$4,regs_set	# r0:1 = dsqrt(x^2+y^2)/2^r6
5432004Szliu	addl2	r6,r0		# unscaled cdabs in r0:1
5532004Szliu	jvc	2b		# unless it overflows
5632004Szliu	subl2	$0x800000,r0	# halve r0 to get meaningful overflow
5732004Szliu	ldd	r0
5832004Szliu	addd	r0		# overflow; r0 is half of true abs value
5932004Szliu	ret
6032004Szliu
6132004Szliuregs_set:
6232004Szliu	.word	0x0000
6332004Szliu	andl2	$0x7fffffff,r0	# r0:r1 = dabs(x)
6432004Szliu	andl2	$0x7fffffff,r2	# r2:r3 = dabs(y)
6532004Szliu	cmpl	r0,r2
6632004Szliu	bgeq	4f
6732004Szliu	movl	r1,r5
6832004Szliu	movl	r0,r4
6932004Szliu	movl	r3,r1
7032004Szliu	movl	r2,r0
7132004Szliu	movl	r5,r3
7232004Szliu	movl	r4,r2		# force y's exp <= x's exp
7332004Szliu4:	andl3	$0xff800000,r0,r6	# r6 = exponent(x) + bias(129)
7432004Szliu	beql	5f		# if x = y = 0 then cdabs(x,y) = 0
7532004Szliu	subl2	$0x47800000,r6	# r6 = exponent(x) - 14
7632004Szliu	subl2	r6,r0		# 2^14 <= scaled x < 2^15
7732004Szliu	bitl	$0xff800000,r2
7832004Szliu	beql	5f		# if y = 0 return dabs(x)
7932004Szliu	subl2	r6,r2
8032004Szliu	cmpl	$0x37800000,r2	# if scaled y < 2^-18
8132004Szliu	bgtr	5f		#   return dabs(x)
8232004Szliu	ldd	r0
8332004Szliu	muld	r0
8432004Szliu	std	r0		# r0:1 = scaled x^2
8532004Szliu	ldd	r2
8632004Szliu	muld	r2		# acc = scaled y^2
8732004Szliu	addd	r0
8832004Szliu	std	r0
8932004Szliu	pushl	r1
9032004Szliu	pushl	r0
9132004Szliu	callf	$12,_sqrt	# r0:1 = dsqrt(x^2+y^2)/2^r6
9232004Szliu5:	ret
93