1/*
2 * Copyright (c) 1992, 1993
3 *	The Regents of the University of California.  All rights reserved.
4 *
5 * This software was developed by the Computer Systems Engineering group
6 * at Lawrence Berkeley Laboratory under DARPA contract BG 91-66 and
7 * contributed to Berkeley.
8 *
9 * Redistribution and use in source and binary forms, with or without
10 * modification, are permitted provided that the following conditions
11 * are met:
12 * 1. Redistributions of source code must retain the above copyright
13 *    notice, this list of conditions and the following disclaimer.
14 * 2. Redistributions in binary form must reproduce the above copyright
15 *    notice, this list of conditions and the following disclaimer in the
16 *    documentation and/or other materials provided with the distribution.
17 * 4. Neither the name of the University nor the names of its contributors
18 *    may be used to endorse or promote products derived from this software
19 *    without specific prior written permission.
20 *
21 * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
22 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
23 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
24 * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
25 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
26 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
27 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
28 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
30 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
31 * SUCH DAMAGE.
32 *
33 * from: Header: modf.s,v 1.3 92/06/20 00:00:54 torek Exp
34 */
35
36#if defined(LIBC_SCCS) && !defined(lint)
37	.asciz "@(#)modf.s	8.1 (Berkeley) 6/4/93"
38#if 0
39	.asciz "$NetBSD: modf.S,v 1.2 2000/07/23 07:12:22 eeh Exp $"
40#endif
41#endif /* LIBC_SCCS and not lint */
42#include <machine/asm.h>
43__FBSDID("$FreeBSD: stable/9/lib/libc/sparc64/gen/modf.S 165903 2007-01-09 00:28:16Z imp $");
44
45#include "assym.s"
46
47/*
48 * double modf(double val, double *iptr)
49 *
50 * Returns the fractional part of `val', storing the integer part of
51 * `val' in *iptr.  Both *iptr and the return value have the same sign
52 * as `val'.
53 *
54 * Method:
55 *
56 * We use the fpu's normalization hardware to compute the integer portion
57 * of the double precision argument.  Sun IEEE double precision numbers
58 * have 52 bits of mantissa, 11 bits of exponent, and one bit of sign,
59 * with the sign occupying bit 31 of word 0, and the exponent bits 30:20
60 * of word 0.  Thus, values >= 2^52 are by definition integers.
61 *
62 * If we take a value that is in the range [+0..2^52) and add 2^52, all
63 * of the fractional bits fall out and all of the integer bits are summed
64 * with 2^52.  If we then subtract 2^52, we get those integer bits back.
65 * This must be done with rounding set to `towards 0' or `towards -inf'.
66 * `Toward -inf' fails when the value is 0 (we get -0 back)....
67 *
68 * Note that this method will work anywhere, but is machine dependent in
69 * various aspects.
70 *
71 * Stack usage:
72 *	4@[%fp + SPOFF - 4]	saved %fsr
73 *	4@[%fp + SPOFF - 8]	new %fsr with rounding set to `towards 0'
74 *	8@[%fp + SPOFF - 16]	space for moving between %i and %f registers
75 * Register usage:
76 *	%f0:f1		double val;
77 *	%l0		scratch
78 *	%l1		sign bit (0x80000000)
79 *	%i1		double *iptr;
80 *	%f2:f3		`magic number' 2^52, in fpu registers
81 *	%f4:f5		double v, in fpu registers
82 *	%f6:f7		double temp.
83 */
84
85	.align	8
86.Lmagic:
87	.word	0x43300000	! sign = 0, exponent = 52 + 1023, mantissa = 0
88	.word	0		! (i.e., .double 0r4503599627370496e+00)
89
90.L0:
91	.word	0		! 0.0
92	.word	0
93
94ENTRY(modf)
95	save	%sp, -CCFSZ - 16, %sp
96	PIC_PROLOGUE(%l6, %l7)
97
98	/*
99	 * First, compute v = abs(val)
100	 */
101	fabsd	%f0, %f4		! %f4:f5 = v
102	fcmped	%fcc1, %f0, %f4		! %fcc1 = (val == abs(val))
103	SET(.Lmagic, %l7, %l0)
104	ldd	[%l0], %f2
105
106	/*
107	 * Is %f4:f5 >= %f2:f3 ?  If so, it is all integer bits.
108	 * It is probably less, though.
109	 */
110	fcmped	%f4, %f2
111	fbuge	.Lbig			! if >= (or unordered), go out
112	 nop
113
114	/*
115	 * v < 2^52, so add 2^52, then subtract 2^52, but do it all
116	 * with rounding set towards zero.  We leave any enabled
117	 * traps enabled, but change the rounding mode.  This might
118	 * not be so good.  Oh well....
119	 */
120	st	%fsr, [%fp + SPOFF - 4]	! %l5 = current FSR mode
121	set	FSR_RD_MASK, %l3	! %l3 = rounding direction mask
122	ld	[%fp + SPOFF - 4], %l5
123	set	FSR_RD_RD_Z, %l4
124	andn	%l5, %l3, %l6
125	or	%l6, %l4, %l6		! round towards zero, please
126	and	%l5, %l3, %l5		! save original rounding mode
127	st	%l6, [%fp + SPOFF - 8]
128	ld	[%fp + SPOFF - 8], %fsr
129
130	faddd	%f4, %f2, %f4		! %f4:f5 += 2^52
131	fsubd	%f4, %f2, %f4		! %f4:f5 -= 2^52
132
133	/*
134	 * Restore %fsr, but leave exceptions accrued.
135	 */
136	st	%fsr, [%fp + SPOFF - 4]
137	ld	[%fp + SPOFF - 4], %l6
138	andn	%l6, %l3, %l6		! %l6 = %fsr & ~FSR_RD_MASK;
139	or	%l5, %l6, %l5		! %l5 |= %l6;
140	st	%l5, [%fp + SPOFF - 4]
141	ld	[%fp + SPOFF - 4], %fsr	! restore %fsr, leaving accrued stuff
142
143	/*
144	 * Now insert the original sign in %f4:f5.
145	 * %fcc1 should still have the results of (val == abs(val))
146	 * from above, so we use a conditional move on %fcc1 to:
147	 *
148	 *	%f4 = (val == abs(val)) ? %f4 : -%f4
149	 *
150	 */
151	fnegd	%f4, %f6
152	fmovdnz	%fcc1, %f6, %f4
1531:
154
155	/*
156	 * The value in %f4:f5 is now the integer portion of the original
157	 * argument.  We need to store this in *ival (%i1), subtract it
158	 * from the original value argument (%d0), and return the result.
159	 */
160	std	%f4, [%i1]		! *ival = %f4:f5;
161	fsubd	%f0, %f4, %f0		! %f0:f1 -= %f4:f5;
162	ret
163	restore
164
165.Lbig:
166	/*
167	 * We get here if the original comparison of %f4:f5 (v) to
168	 * %f2:f3 (2^52) came out `greater or unordered'.  In this
169	 * case the integer part is the original value, and the
170	 * fractional part is 0.
171	 */
172	SET(.L0, %l7, %l0)
173	std	%f0, [%i1]		! *ival = val;
174	ldd	[%l0], %f0		! return 0.0;
175	ret
176	restore
177END(modf)
178