xref: /csrg-svn/lib/libm/mc68881/support.s (revision 42210)
1/*-
2 * Copyright (c) 1990 The Regents of the University of California.
3 * All rights reserved.
4 *
5 * This code is derived from software contributed to Berkeley by
6 * the Systems Programming Group of the University of Utah Computer
7 * Science Department.
8 *
9 * %sccs.include.redist.c%
10 *
11 *	@(#)support.s	5.1 (Berkeley) 05/17/90
12 */
13
14	.text
15	.globl	_copysign, _finite, _scalb, _logb, _drem, _pow_p, _atan2__A
16
17| copysign(x,y)
18| returns x with the sign of y.
19_copysign:
20	movl	sp@(4),d0
21	movl	sp@(8),d1
22	tstw	sp@(12)
23	jmi	Lneg
24	bclr	#31,d0
25	rts
26Lneg:
27	bset	#31,d0
28	rts
29
30| finite(x)
31| returns the value TRUE if -INF < x < +INF and returns FALSE otherwise.
32_finite:
33	movw	#0x7FF0,d0
34	movw	sp@(4),d1
35	andw	d0,d1
36	cmpw	d0,d1
37	beq	Lnotfin
38	moveq	#1,d0
39	rts
40Lnotfin:
41	clrl	d0
42	rts
43
44| scalb(x, N)
45| returns  x * (2**N), for integer values N.
46_scalb:
47	fmoved	sp@(4),fp0
48	fbeq	Ldone
49	ftwotoxl	sp@(12),fp1
50	fmoved	fp1,sp@-
51	fmuld	sp@+,fp0
52Ldone:
53	fmoved	fp0,sp@-
54	movel	sp@+,d0
55	movel	sp@+,d1
56	rts
57
58| logb(x)
59| returns the unbiased exponent of x, a signed integer in double precision,
60| except that logb(0) is -INF, logb(INF) is +INF, and logb(NAN) is that NAN.
61_logb:
62	movw	sp@(4),d0
63	movw	#0x7FF0,d1	| exponent bits
64	andw	d1,d0		| mask off all else
65	cmpw	d1,d0		| max exponent?
66	bne	Lfinite		| no, is finite
67	fmoved	sp@(4),fp0	| yes, infinite or NaN
68	fbun	Ldone		| NaN returns NaN
69	fabsx	fp0		| +-inf returns inf
70	jra	Ldone
71Lfinite:
72	fmoved	sp@(4),fp0	| get entire number
73	fbne	Lnonz		| zero?
74	flog2x	fp0		| yes, log(0) a convenient source of -inf
75	jra	Ldone
76Lnonz:
77	fgetexpx	fp0	| get exponent
78	jra	Ldone
79
80| drem(x,y)
81| returns  x REM y  =  x - [x/y]*y , where [x/y] is the integer nearest x/y;
82| in half way case, choose the even one.
83_drem:
84	fmoved	sp@(4),fp0
85	fremd	sp@(12),fp0
86	fmoved	fp0,sp@-
87	movel	sp@+,d0
88	movel	sp@+,d1
89	rts
90
91| pow_p(x,y)
92| return x**y for x with sign=1 and finite y
93_pow_p:
94	flognd	sp@(4),fp0
95	fmuld	sp@(12),fp0
96	fetoxx	fp0
97	fmoved	fp0,sp@-
98	movel	sp@+,d0
99	movel	sp@+,d1
100	rts
101
102| atan2__A(y,x)
103| compute atan2(y,x) where x,y are finite and non-zero
104| called by atan2() after weeding out all the special cases
105_atan2__A:
106	moveq	#0,d0		| sign of result
107	fmoved	sp@(4),fp0	| get y
108	fboge	Lypos		| <0?
109	moveq	#1,d0		| yes, result is neg
110	fnegx	fp0		| make y pos
111Lypos:
112	fmoved	sp@(12),fp1	| get x
113	fboge	Lxpos		| <0?
114	fnegx	fp1		| yes, make x pos
115	fdivx	fp1,fp0		| y/x
116	fatanx	fp0,fp1		| atan(y/x)
117	fmovecr	#0,fp0		| get pi
118	fsubx	fp1,fp0		| pi - atan(y/x)
119	jra	Lsetsign
120Lxpos:
121	fdivx	fp1,fp0		| y/x
122	fatanx	fp0		| atan(y/x)
123Lsetsign:
124	tstl	d0		| should be neg?
125	jeq	Lrpos		| no, all done
126	fnegx	fp0		| yes, negate
127Lrpos:
128	fmoved	fp0,sp@-
129	movel	sp@+,d0
130	movel	sp@+,d1
131	rts
132