modf.S revision 1.4
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.4 1996/11/16 20:53:20 pk 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 <machine/asm.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 %hi(0x80000000), %l1 ! sign bit 103 andn %i0, %l1, %l0 104 st %l0, [%fp - 16] 105 #ifdef PIC 106 PICCY_SET(Lmagic, %l0, %o7) 107 ldd [%l0], %f2 108 #else 109 sethi %hi(Lmagic), %l0 110 ldd [%l0 + %lo(Lmagic)], %f2 111 #endif 112 st %i1, [%fp - 12] 113 ldd [%fp - 16], %f4 ! %f4:f5 = v 114 115 /* 116 * Is %f4:f5 >= %f2:f3 ? If so, it is all integer bits. 117 * It is probably less, though. 118 */ 119 fcmped %f4, %f2 120 nop ! fpop2 delay 121 fbuge Lbig ! if >= (or unordered), go out 122 nop 123 124 /* 125 * v < 2^52, so add 2^52, then subtract 2^52, but do it all 126 * with rounding set towards zero. We leave any enabled 127 * traps enabled, but change the rounding mode. This might 128 * not be so good. Oh well.... 129 */ 130 st %fsr, [%fp - 4] ! %l5 = current FSR mode 131 set FSR_RD, %l3 ! %l3 = rounding direction mask 132 ld [%fp - 4], %l5 133 set FSR_RD_RZ << FSR_RD_SHIFT, %l4 134 andn %l5, %l3, %l6 135 or %l6, %l4, %l6 ! round towards zero, please 136 and %l5, %l3, %l5 ! save original rounding mode 137 st %l6, [%fp - 8] 138 ld [%fp - 8], %fsr 139 140 faddd %f4, %f2, %f4 ! %f4:f5 += 2^52 141 fsubd %f4, %f2, %f4 ! %f4:f5 -= 2^52 142 143 /* 144 * Restore %fsr, but leave exceptions accrued. 145 */ 146 st %fsr, [%fp - 4] 147 ld [%fp - 4], %l6 148 andn %l6, %l3, %l6 ! %l6 = %fsr & ~FSR_RD; 149 or %l5, %l6, %l5 ! %l5 |= %l6; 150 st %l5, [%fp - 4] 151 ld [%fp - 4], %fsr ! restore %fsr, leaving accrued stuff 152 153 /* 154 * Now insert the original sign in %f4:f5. 155 * This is a lot of work, so it is conditional here. 156 */ 157 btst %l1, %i0 158 be 1f 159 nop 160 st %f4, [%fp - 16] 161 ld [%fp - 16], %g1 162 or %l1, %g1, %g1 163 st %g1, [%fp - 16] 164 ld [%fp - 16], %f4 165 1: 166 167 /* 168 * The value in %f4:f5 is now the integer portion of the original 169 * argument. We need to store this in *ival (%i2), subtract it 170 * from the original value argument (%i0:i1), and return the result. 171 */ 172 std %f4, [%i2] ! *ival = %f4:f5; 173 std %i0, [%fp - 16] 174 ldd [%fp - 16], %f0 ! %f0:f1 = val; 175 fsubd %f0, %f4, %f0 ! %f0:f1 -= %f4:f5; 176 ret 177 restore 178 179 Lbig: 180 /* 181 * We get here if the original comparison of %f4:f5 (v) to 182 * %f2:f3 (2^52) came out `greater or unordered'. In this 183 * case the integer part is the original value, and the 184 * fractional part is 0. 185 */ 186 #ifdef PIC 187 PICCY_SET(L0, %l0, %o7) 188 std %f0, [%i2] ! *ival = val; 189 ldd [%l0], %f0 ! return 0.0; 190 #else 191 sethi %hi(L0), %l0 192 std %f0, [%i2] ! *ival = val; 193 ldd [%l0 + %lo(L0)], %f0 ! return 0.0; 194 #endif 195 ret 196 restore 197
Indexes created Mon Sep 29 03:10:08 GMT 2025