|
|
1.1 ! root 1: ////////// ! 2: / libc/crt/i386/daddsub.s ! 3: / i386 C runtime library. ! 4: / IEEE software floating point support. ! 5: ////////// ! 6: ! 7: ////////// ! 8: / double _dadd(double d) ! 9: / Return d + %edx:eax in %edx:%eax. ! 10: / ! 11: / double _dsub(double d) ! 12: / Return d - %edx:eax in %edx:%eax. ! 13: / ! 14: / double _drsub(double d) ! 15: / Return %edx:eax - d in %edx:%eax. ! 16: / ! 17: / Addition is easy: align the mantissa of the smaller operand with the larger, ! 18: / possibly rounding up; add the aligned mantissas; shift right 1 bit ! 19: / if necessary to normalize the result mantissa. ! 20: / ! 21: / Subtraction is a bit trickier: normalizing the result is nontrivial, ! 22: / and we must be careful not to lose precision. A simple example: ! 23: / A = 0x10 0000 0000 0001 53 bits ! 24: / B = 0x8 0000 0000 0000 . 8 53 bits ! 25: / A-B = 0x8 0000 0000 0000 . 8 53 bits ! 26: / If we justify the mantissa of B with A before subtracting, we round it to ! 27: / B' = 0x8 0000 0000 0001 52 bits ! 28: / and the result of the subtraction is ! 29: / A-B' = 0x8 0000 0000 0000 52 bits ! 30: / so we have lost a bit of precision. ! 31: / Shifting the larger operand left one bit eliminates this problem: ! 32: / A' = 0x10 0000 0000 0001 . 0 54 bits (1 after '.') ! 33: / A'-B = 0x8 0000 0000 0000 . 8 53 bits ! 34: / which is the correct result. Shifting once is enough, because subtraction ! 35: / can clear at most the first result bit when operands are not aligned. ! 36: ////////// ! 37: ! 38: d .equ 8 ! 39: EXPMASK .equ 0x7FF00000 ! 40: MANMASK .equ 0x000FFFFF ! 41: SGNMASK .equ 0x80000000 ! 42: HIDDEN .equ 0x00100000 ! 43: DMBITS .equ 53 ! 44: ! 45: .globl _dadd ! 46: .globl _dsub ! 47: .globl _drsub ! 48: ! 49: _drsub: ! 50: xchgl %edx, 8(%esp) ! 51: xchgl %eax, 4(%esp) / exchange arg order ! 52: / jmp _dsub / and fall through to subtract ! 53: ! 54: _dsub: ! 55: xorl %edx, $SGNMASK / flip rhs sign ! 56: / jmp _dadd / and fall through to add ! 57: ! 58: _dadd: ! 59: push %ebp ! 60: mov %ebp, %esp ! 61: push %esi ! 62: push %edi ! 63: push %ebx ! 64: push %ecx ! 65: ! 66: movl %esi, d+4(%ebp) ! 67: movl %edi, d(%ebp) / d to ESI:EDI ! 68: / now done with EBP as index register ! 69: ! 70: / Check for special cases +-0.0, +-infinity, NaN on each side. ! 71: / Ignore denormals. ! 72: movl %ebx, %esi ! 73: andl %ebx, $EXPMASK ! 74: jz ?done / lhs 0.0, return rhs; ignore denormal ! 75: cmpl %ebx, $EXPMASK ! 76: jz ?retlhs / lhs +-infinity or NaN, return it ! 77: movl %ecx, %edx ! 78: andl %ecx, $EXPMASK ! 79: jz ?retlhs / rhs 0.0, return lhs; ignore denormal ! 80: cmpl %ecx, $EXPMASK ! 81: jz ?done / rhs +-infinity or NaN, return it ! 82: ! 83: / Force arg with larger exponent to EDX:EAX. ! 84: shrl %ebx, $20 / lhs biased exp in EBX ! 85: shrl %ecx, $20 / rhs biased exp in ECX ! 86: cmpl %ecx, %ebx ! 87: jge ?expdiff ! 88: xchgl %esi, %edx ! 89: xchgl %edi, %eax / force |EDX:EAX| >= |ESI:EDI| (maybe) ! 90: xchgl %ebx, %ecx / and exchange exponents accordingly ! 91: ! 92: / Compute the exponent difference in preparation for shifting ESI:EDI. ! 93: / We know that exponent(EDX:EAX) >= exponent(ESI:EDI). ! 94: ?expdiff: ! 95: movl %ebp, %ecx / save EDX:EAX exponent in EBP ! 96: subl %ecx, %ebx / nonnegative exponent difference to ECX ! 97: cmpl %ecx, $DMBITS+1 ! 98: jg ?done / ESI:EDI insignificant, return EDX:EAX as is ! 99: ! 100: movl %ebx, %edx ! 101: xorl %ebx, %esi / sign bit 0 iff arg signs match ! 102: pushfl / save as add/sub flag ! 103: movl %ebx, %edx ! 104: andl %ebx, $SGNMASK / save result sign bit in EBX ! 105: andl %edx, $MANMASK ! 106: orl %edx, $HIDDEN / extract mantissa, restore hidden bit ! 107: andl %esi, $MANMASK ! 108: orl %esi, $HIDDEN / extract mantissa, restore hidden bit ! 109: ! 110: / Align mantissas by shifting ESI:EDI to the right position. ! 111: jecxz ?addsub / skip shift if exponents match ! 112: popfl / recover add/sub flag ! 113: pushfl ! 114: jns ?align / addition, just adjust ! 115: ! 116: / For subtraction with differing exponents, shift args left one bit ! 117: / to avoid loss of precision; see comment at top. ! 118: shld %esi, %edi, $1 ! 119: shll %edi, $1 / shift ESI:EDI left 1 bit ! 120: shld %edx, %eax, $1 ! 121: shll %eax, $1 / shift EDX:EAX left 1 bit ! 122: decl %ebp / adjust result exponent ! 123: ! 124: ?align: ! 125: cmpl %ecx, $32 ! 126: jl ?adj ! 127: movl %edi, %esi ! 128: subl %esi, %esi / shift ESI:EDI right by 32 bits ! 129: ?adj: ! 130: shrd %edi, %esi, %cl ! 131: pushfl ! 132: shrl %esi, %cl / shift ESI:EDI right by CL mod 32 ! 133: popfl / carry indicates rounding ! 134: adcl %edi, $0 / round up ! 135: adcl %esi, $0 ! 136: ! 137: ?addsub: ! 138: popfl / restore add/sub flag ! 139: js ?sub ! 140: addl %eax, %edi / add lo mantissas ! 141: adcl %edx, %esi / and hi mantissas ! 142: cmpl %edx, $HIDDEN<<1 / check for carry past hidden bit ! 143: jl ?pack / done if none ! 144: ! 145: / Shift result mantissa right one bit. ! 146: ?rshift: ! 147: incl %ebp / bump result exponent ! 148: cmpl %ebp, $EXPMASK ! 149: je ?inf / exponent overflow, return +-infinity ! 150: shrd %eax, %edx, $1 ! 151: pushfl ! 152: shrl %edx, $1 / shift result right one bit ! 153: popfl ! 154: adcl %eax, $0 ! 155: adcl %edx, $0 / round up, carry to hidden bit impossible ! 156: ! 157: / Pack result mantissa from EDX:EAX, exponent from EBP, sign from EBX. ! 158: ?pack: ! 159: andl %edx, $MANMASK / mask off hidden bit ! 160: orl %edx, %ebx / pack with sign ! 161: shll %ebp, $20 / position exponent ! 162: orl %edx, %ebp / pack with exponent ! 163: ! 164: / Return the value from EDX:EAX. ! 165: ?done: ! 166: pop %ecx ! 167: pop %ebx ! 168: pop %edi ! 169: pop %esi ! 170: pop %ebp ! 171: ret ! 172: ! 173: / Return the value from ESI:EDI. ! 174: ?retlhs: ! 175: xchgl %esi, %edx ! 176: xchgl %edi, %eax ! 177: jmp ?done ! 178: ! 179: / Return +-infinity. ! 180: ?inf: ! 181: subl %eax, %eax ! 182: subl %edx, %edx / zero mantissa ! 183: jmp ?pack ! 184: ! 185: / Subtract. ! 186: / We want the result mantissa to be nonnegative. ! 187: / The exponent check above assures that exp(EDX:EAX) >= exp(ESI:EDI), ! 188: / but when the exponents are equal we must compare mantissas. ! 189: ?sub: ! 190: jcxz ?argtest / compare mantissas if exponents equal ! 191: ! 192: ?dosub: ! 193: subl %eax, %edi / subtract lo mantissas ! 194: sbbl %edx, %esi / and hi mantissas ! 195: ! 196: / Normalize the difference. ! 197: / The code would be simpler if it used iterated left shift by 1 ! 198: / to find hidden bit, but the reverse bit scan used here is faster. ! 199: / Bit 21 can be set if the EDX:EAX operand was shifted left above. ! 200: bsrl %ecx, %edx / hi bit position 0-21 to ECX ! 201: jz ?lotest / hi mantissa is 0 ! 202: cmpl %ecx, $21 ! 203: jz ?rshift / right shift one bit if hi bit was 21 ! 204: ! 205: ?normalize: ! 206: negl %ecx ! 207: addl %ecx, $20 / shift count to CL ! 208: ! 209: / Shift left by CL bits to normalize, adjusting exponent accordingly. ! 210: ?lshift: ! 211: shld %edx, %eax, %cl ! 212: shll %eax, %cl / shift EDX:EAX left CL bits ! 213: subl %ebp, %ecx / adjust result exponent ! 214: jg ?pack / pack and go home ! 215: subl %eax, %eax / exponent underflow, return 0.0 ! 216: subl %edx, %edx ! 217: jmp ?done ! 218: ! 219: / Hi mantissa EDX is 0, check lo mantissa EAX. ! 220: ?lotest: ! 221: bsrl %ecx, %eax / hi bit position 0-31 to ECX ! 222: jz ?done / mantissa difference 0, return 0.0 ! 223: cmpl %ecx, $20 ! 224: jle ?losmall / hi bit of lo mantissa is 0 to 20 ! 225: / else hi bit is bit 21 to 31 ! 226: negl %ecx / -21 to -31 ! 227: addl %ecx, $52 / 31 to 21 shift count in CL ! 228: jmp ?lshift / shift EDX:EAX as above ! 229: ! 230: / EDX is 0, hi bit of EAX is bit 0 to 20; move EAX up to EDX. ! 231: ?losmall: ! 232: xchgl %edx, %eax / shift left 32 bits; EAX becomes 0 ! 233: subl %ebp, $32 / adjust exponent ! 234: jmp ?normalize / and shift as above ! 235: ! 236: / Compare mantissas, exchange args to force |EDX:EAX| >= |ESI:EDI|. ! 237: ?argtest: ! 238: cmpl %edx, %esi / compare hi mantissas ! 239: ja ?dosub ! 240: jb ?flip ! 241: cmpl %eax, %edi / compare lo mantissas ! 242: jae ?dosub ! 243: ! 244: / |ESI:EDI| > |EDX:EAX|, flip args and change result sign. ! 245: ?flip: ! 246: xchgl %esi, %edx ! 247: xchgl %edi, %eax / force |EDX:EAX| >= |ESI:EDI| (for sure) ! 248: xorl %ebx, $SGNMASK / and change the result sign ! 249: jmp ?dosub ! 250: ! 251: / end of libc/crt/i386/daddsub.s
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.