ldexp.S revision 1.1
1 /*- 2 * Copyright (c) 1991, 1993 3 * The Regents of the University of California. All rights reserved. 4 * 5 * This code is derived from software contributed to Berkeley by 6 * Ralph Campbell. 7 * 8 * Redistribution and use in source and binary forms, with or without 9 * modification, are permitted provided that the following conditions 10 * are met: 11 * 1. Redistributions of source code must retain the above copyright 12 * notice, this list of conditions and the following disclaimer. 13 * 2. Redistributions in binary form must reproduce the above copyright 14 * notice, this list of conditions and the following disclaimer in the 15 * documentation and/or other materials provided with the distribution. 16 * 3. All advertising materials mentioning features or use of this software 17 * must display the following acknowledgement: 18 * This product includes software developed by the University of 19 * California, Berkeley and its contributors. 20 * 4. Neither the name of the University nor the names of its contributors 21 * may be used to endorse or promote products derived from this software 22 * without specific prior written permission. 23 * 24 * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND 25 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 26 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 27 * ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE 28 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 29 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS 30 * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 31 * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 32 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY 33 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF 34 * SUCH DAMAGE. 35 */ 36 37 #include <machine/machAsmDefs.h> 38 39 #if defined(LIBC_SCCS) && !defined(lint) 40 ASMSTR("from: @(#)ldexp.s 8.1 (Berkeley) 6/4/93") 41 ASMSTR("$Id: ldexp.S,v 1.1 1994/05/24 07:12:33 glass Exp $") 42 #endif /* LIBC_SCCS and not lint */ 43 44 #define DEXP_INF 0x7ff 45 #define DEXP_BIAS 1023 46 #define DEXP_MIN -1022 47 #define DEXP_MAX 1023 48 #define DFRAC_BITS 52 49 #define DIMPL_ONE 0x00100000 50 #define DLEAD_ZEROS 31 - 20 51 #define STICKYBIT 1 52 #define GUARDBIT 0x80000000 53 #define DSIGNAL_NAN 0x00040000 54 #define DQUIET_NAN0 0x0007ffff 55 #define DQUIET_NAN1 0xffffffff 56 57 /* 58 * double ldexp(x, N) 59 * double x; int N; 60 * 61 * Return x * (2**N), for integer values N. 62 */ 63 LEAF(ldexp) 64 mfc1 v1, $f13 # get MSW of x 65 mfc1 t3, $f12 # get LSW of x 66 sll t1, v1, 1 # get x exponent 67 srl t1, t1, 32 - 11 68 beq t1, DEXP_INF, 9f # is it a NAN or infinity? 69 beq t1, zero, 1f # zero or denormalized number? 70 addu t1, t1, a2 # scale exponent 71 sll v0, a2, 20 # position N for addition 72 bge t1, DEXP_INF, 8f # overflow? 73 addu v0, v0, v1 # multiply by (2**N) 74 ble t1, zero, 4f # underflow? 75 mtc1 v0, $f1 # save MSW of result 76 mtc1 t3, $f0 # save LSW of result 77 j ra 78 1: 79 sll t2, v1, 32 - 20 # get x fraction 80 srl t2, t2, 32 - 20 81 srl t0, v1, 31 # get x sign 82 bne t2, zero, 1f 83 beq t3, zero, 9f # result is zero 84 1: 85 /* 86 * Find out how many leading zero bits are in t2,t3 and put in t9. 87 */ 88 move v0, t2 89 move t9, zero 90 bne t2, zero, 1f 91 move v0, t3 92 addu t9, 32 93 1: 94 srl t4, v0, 16 95 bne t4, zero, 1f 96 addu t9, 16 97 sll v0, 16 98 1: 99 srl t4, v0, 24 100 bne t4, zero, 1f 101 addu t9, 8 102 sll v0, 8 103 1: 104 srl t4, v0, 28 105 bne t4, zero, 1f 106 addu t9, 4 107 sll v0, 4 108 1: 109 srl t4, v0, 30 110 bne t4, zero, 1f 111 addu t9, 2 112 sll v0, 2 113 1: 114 srl t4, v0, 31 115 bne t4, zero, 1f 116 addu t9, 1 117 /* 118 * Now shift t2,t3 the correct number of bits. 119 */ 120 1: 121 subu t9, t9, DLEAD_ZEROS # dont count normal leading zeros 122 li t1, DEXP_MIN + DEXP_BIAS 123 subu t1, t1, t9 # adjust exponent 124 addu t1, t1, a2 # scale exponent 125 li v0, 32 126 blt t9, v0, 1f 127 subu t9, t9, v0 # shift fraction left >= 32 bits 128 sll t2, t3, t9 129 move t3, zero 130 b 2f 131 1: 132 subu v0, v0, t9 # shift fraction left < 32 bits 133 sll t2, t2, t9 134 srl t4, t3, v0 135 or t2, t2, t4 136 sll t3, t3, t9 137 2: 138 bge t1, DEXP_INF, 8f # overflow? 139 ble t1, zero, 4f # underflow? 140 sll t2, t2, 32 - 20 # clear implied one bit 141 srl t2, t2, 32 - 20 142 3: 143 sll t1, t1, 31 - 11 # reposition exponent 144 sll t0, t0, 31 # reposition sign 145 or t0, t0, t1 # put result back together 146 or t0, t0, t2 147 mtc1 t0, $f1 # save MSW of result 148 mtc1 t3, $f0 # save LSW of result 149 j ra 150 4: 151 li v0, 0x80000000 152 ble t1, -52, 7f # is result too small for denorm? 153 sll t2, v1, 31 - 20 # clear exponent, extract fraction 154 or t2, t2, v0 # set implied one bit 155 blt t1, -30, 2f # will all bits in t3 be shifted out? 156 srl t2, t2, 31 - 20 # shift fraction back to normal position 157 subu t1, t1, 1 158 sll t4, t2, t1 # shift right t2,t3 based on exponent 159 srl t8, t3, t1 # save bits shifted out 160 negu t1 161 srl t3, t3, t1 162 or t3, t3, t4 163 srl t2, t2, t1 164 bge t8, zero, 1f # does result need to be rounded? 165 addu t3, t3, 1 # round result 166 sltu t4, t3, 1 167 sll t8, t8, 1 168 addu t2, t2, t4 169 bne t8, zero, 1f # round result to nearest 170 and t3, t3, ~1 171 1: 172 mtc1 t3, $f0 # save denormalized result (LSW) 173 mtc1 t2, $f1 # save denormalized result (MSW) 174 bge v1, zero, 1f # should result be negative? 175 neg.d $f0, $f0 # negate result 176 1: 177 j ra 178 2: 179 mtc1 zero, $f1 # exponent and upper fraction 180 addu t1, t1, 20 # compute amount to shift right by 181 sll t8, t2, t1 # save bits shifted out 182 negu t1 183 srl t3, t2, t1 184 bge t8, zero, 1f # does result need to be rounded? 185 addu t3, t3, 1 # round result 186 sltu t4, t3, 1 187 sll t8, t8, 1 188 mtc1 t4, $f1 # exponent and upper fraction 189 bne t8, zero, 1f # round result to nearest 190 and t3, t3, ~1 191 1: 192 mtc1 t3, $f0 193 bge v1, zero, 1f # is result negative? 194 neg.d $f0, $f0 # negate result 195 1: 196 j ra 197 7: 198 mtc1 zero, $f0 # result is zero 199 mtc1 zero, $f1 200 beq t0, zero, 1f # is result positive? 201 neg.d $f0, $f0 # negate result 202 1: 203 j ra 204 8: 205 li t1, 0x7ff00000 # result is infinity (MSW) 206 mtc1 t1, $f1 207 mtc1 zero, $f0 # result is infinity (LSW) 208 bge v1, zero, 1f # should result be negative infinity? 209 neg.d $f0, $f0 # result is negative infinity 210 1: 211 add.d $f0, $f0 # cause overflow faults if enabled 212 j ra 213 9: 214 mov.d $f0, $f12 # yes, result is just x 215 j ra 216 END(ldexp) 217
Indexes created Thu Oct 02 14:10:14 GMT 2025