modf.S revision 1.1
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 * 3. All advertising materials mentioning features or use of this software 18 * must display the following acknowledgement: 19 * This product includes software developed by the University of 20 * California, Berkeley and its contributors. 21 * 4. Neither the name of the University nor the names of its contributors 22 * may be used to endorse or promote products derived from this software 23 * without specific prior written permission. 24 * 25 * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND 26 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 27 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 28 * ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE 29 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 30 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS 31 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 32 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 33 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY 34 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF 35 * SUCH DAMAGE. 36 * 37 * from: Header: modf.s,v 1.3 92/06/20 00:00:54 torek Exp 38 * $Id: modf.S,v 1.1 1993/10/07 00:21:37 cgd Exp $ 39 */ 40 41 #if defined(LIBC_SCCS) && !defined(lint) 42 .asciz "@(#)modf.s 8.1 (Berkeley) 6/4/93" 43 #endif /* LIBC_SCCS and not lint */ 44 45 #include "DEFS.h" 46 #include <machine/fsr.h> 47 48 /* 49 * double modf(double val, double *iptr) 50 * 51 * Returns the fractional part of `val', storing the integer part of 52 * `val' in *iptr. Both *iptr and the return value have the same sign 53 * as `val'. 54 * 55 * Method: 56 * 57 * We use the fpu's normalization hardware to compute the integer portion 58 * of the double precision argument. Sun IEEE double precision numbers 59 * have 52 bits of mantissa, 11 bits of exponent, and one bit of sign, 60 * with the sign occupying bit 31 of word 0, and the exponent bits 30:20 61 * of word 0. Thus, values >= 2^52 are by definition integers. 62 * 63 * If we take a value that is in the range [+0..2^52) and add 2^52, all 64 * of the fractional bits fall out and all of the integer bits are summed 65 * with 2^52. If we then subtract 2^52, we get those integer bits back. 66 * This must be done with rounding set to `towards 0' or `towards -inf'. 67 * `Toward -inf' fails when the value is 0 (we get -0 back).... 68 * 69 * Note that this method will work anywhere, but is machine dependent in 70 * various aspects. 71 * 72 * Stack usage: 73 * 4@[%fp - 4] saved %fsr 74 * 4@[%fp - 8] new %fsr with rounding set to `towards 0' 75 * 8@[%fp - 16] space for moving between %i and %f registers 76 * Register usage: 77 * %i0%i1 double val; 78 * %l0 scratch 79 * %l1 sign bit (0x80000000) 80 * %i2 double *iptr; 81 * %f2:f3 `magic number' 2^52, in fpu registers 82 * %f4:f5 double v, in fpu registers 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 94 ENTRY(modf) 95 save %sp, -64-16, %sp 96 97 /* 98 * First, compute v = abs(val) by clearing sign bit, 99 * and then set up the fpu registers. This would be 100 * much easier if we could do alu operations on fpu registers! 101 */ 102 sethi 0x80000000, %l1 ! sign bit 103 andn %i0, %l1, %l0 104 st %l0, [%fp - 16] 105 sethi %hi(Lmagic), %l0 106 ldd [%l0 + %lo(Lmagic)], %f2 107 st %i1, [%fp - 12] 108 ldd [%fp - 16], %f4 ! %f4:f5 = v 109 110 /* 111 * Is %f4:f5 >= %f2:f3 ? If so, it is all integer bits. 112 * It is probably less, though. 113 */ 114 fcmped %f4, %f2 115 nop ! fpop2 delay 116 fbuge Lbig ! if >= (or unordered), go out 117 nop 118 119 /* 120 * v < 2^52, so add 2^52, then subtract 2^52, but do it all 121 * with rounding set towards zero. We leave any enabled 122 * traps enabled, but change the rounding mode. This might 123 * not be so good. Oh well.... 124 */ 125 st %fsr, [%fp - 4] ! %l5 = current FSR mode 126 set FSR_RD, %l3 ! %l3 = rounding direction mask 127 ld [%fp - 4], %l5 128 set FSR_RD_RZ << FSR_RD_SHIFT, %l4 129 andn %l5, %l3, %l6 130 or %l6, %l4, %l6 ! round towards zero, please 131 and %l5, %l3, %l5 ! save original rounding mode 132 st %l6, [%fp - 8] 133 ld [%fp - 8], %fsr 134 135 faddd %f4, %f2, %f4 ! %f4:f5 += 2^52 136 fsubd %f4, %f2, %f4 ! %f4:f5 -= 2^52 137 138 /* 139 * Restore %fsr, but leave exceptions accrued. 140 */ 141 st %fsr, [%fp - 4] 142 ld [%fp - 4], %l6 143 andn %l6, %l3, %l6 ! %l6 = %fsr & ~FSR_RD; 144 or %l5, %l6, %l5 ! %l5 |= %l6; 145 st %l5, [%fp - 4] 146 ld [%fp - 4], %fsr ! restore %fsr, leaving accrued stuff 147 148 /* 149 * Now insert the original sign in %f4:f5. 150 * This is a lot of work, so it is conditional here. 151 */ 152 btst %l1, %i0 153 be 1f 154 nop 155 st %f4, [%fp - 16] 156 ld [%fp - 16], %g1 157 or %l1, %g1, %g1 158 st %g1, [%fp - 16] 159 ld [%fp - 16], %f4 160 1: 161 162 /* 163 * The value in %f4:f5 is now the integer portion of the original 164 * argument. We need to store this in *ival (%i2), subtract it 165 * from the original value argument (%i0:i1), and return the result. 166 */ 167 std %f4, [%i2] ! *ival = %f4:f5; 168 std %i0, [%fp - 16] 169 ldd [%fp - 16], %f0 ! %f0:f1 = val; 170 fsubd %f0, %f4, %f0 ! %f0:f1 -= %f4:f5; 171 ret 172 restore 173 174 Lbig: 175 /* 176 * We get here if the original comparison of %f4:f5 (v) to 177 * %f2:f3 (2^52) came out `greater or unordered'. In this 178 * case the integer part is the original value, and the 179 * fractional part is 0. 180 */ 181 sethi %hi(L0), %l0 182 std %f0, [%i2] ! *ival = val; 183 ldd [%l0 + %lo(L0)], %f0 ! return 0.0; 184 ret 185 restore 186
Indexes created Thu Sep 25 16:09:42 GMT 2025