![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dmul.s CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [MW Coherent from dump] / coherent / b / lib / libc / crt / i386|
Return to dmul.s CVS log | Up to [MW Coherent from dump] / coherent / b / lib / libc / crt / i386 |
1.1 root 1: //////////
2: / libc/crt/i386/dmul.s
3: / i386 C runtime library.
4: / IEEE software floating point support.
5: //////////
6:
7: //////////
8: / double _dmul(double d)
9: / Return d * %edx:eax in %edx:%eax.
10: /
11: / In 32-bit dwords:
12: / A = hiA:loA
13: / B = hiB:loB
14: / A*B = hi(hiA*hiB):lo(hiA*hiB)
15: / + hi(hiA*loB):lo(hiA*loB)
16: / + hi(loA*hiB):lo(loA*hiB)
17: / + lo(loA*loB):lo(loA*loB)
18: / Multiplying two n bit quantities produces a 2*n-1 or 2*n bit result.
19: / Thus, multiplying two 53 bit quantities produces a 105 or 106 bit result.
20: / We want a normalized 53-bit result, so the lo-order dword above is irrelevant.
21: //////////
22:
23: d .equ 8
24: BIAS .equ 1023
25: EXPMASK .equ 0x7FF00000
26: MANMASK .equ 0x000FFFFF
27: SGNMASK .equ 0x80000000
28: HIDDEN .equ 0x00100000
29:
30: .globl _dmul
31:
32: _dmul:
33: push %ebp
34: movl %ebp, %esp
35: push %esi
36: push %edi
37: push %ebx
38: push %ecx
39:
40: movl %esi, d+4(%ebp)
41: movl %edi, d(%ebp) / d to ESI:EDI, call it A
42: / now done with EBP as index register
43:
44: / Check for special cases +-0.0 and +-infinity on each side.
45: movl %ebx, %esi
46: andl %ebx, $EXPMASK
47: jz ?retlhs / A is 0.0, return it; ignore denormal
48: cmpl %ebx, $EXPMASK
49: jz ?retlhs / A is +-infinity or NaN, return it
50: movl %ecx, %edx
51: andl %ecx, $EXPMASK
52: jz ?done / B is 0.0, return it; ignore denormal
53: cmpl %ecx, $EXPMASK
54: jz ?done / B is +-infinity or NaN, return it
55:
56: / Compute result sign.
57: movl %ebp, %edx
58: xorl %ebp, %esi
59: andl %ebp, $SGNMASK
60: push %ebp / save result sign bit
61:
62: / Compute result exponent.
63: shrl %ebx, $20 / A biased exp in EBX
64: shrl %ecx, $20 / B biased exp in ECX
65: subl %ecx, $BIAS / unbiased
66: addl %ecx, %ebx / form biased result exponent
67: jl ?zero / underflow, return 0.0
68: cmpl %ecx, $EXPMASK>>20
69: jge ?inf / overflow, return +-infinity
70: movl %ebp, %ecx / save result exponent in EBP
71:
72: / Perform the 64*64 bit multiply.
73: andl %esi, $MANMASK
74: orl %esi, $HIDDEN / extract A mantissa, restore hidden bit
75: andl %edx, $MANMASK
76: orl %edx, $HIDDEN / extract B mantissa, restore hidden bit
77: movl %ecx, %edx
78: movl %ebx, %eax / B mantissa to ECX:EBX
79:
80: mul %edi / loA*loB
81: push %edx / save hi(loA*loB), ignore lo dword
82: movl %eax, %ebx / loB
83: mul %esi / hiA*loB
84: xchgl %edi, %eax / lo(hiA*loB) to EDI, loA to EAX
85: movl %ebx, %edx / hi(hiA*loB) to EBX
86: mul %ecx / loA*hiB
87: xchgl %esi, %eax / lo(loA*hiB) to ESI, hiA to EAX
88: xchgl %ecx, %edx / hi(loA*hiB) to ECX, hiB to EDX
89: mul %edx / hiA*hiB, leave in EDX:EAX
90:
91: / Add partial products to obtain 96 bit product in EDX:EAX:ESI.
92: / This is really a 105-32=73 or 106-32=74 bit value, higher bits are 0.
93: addl %esi, %edi
94: adcl %eax, %ebx
95: adcl %edx, $0
96: pop %edi / hi(loA*loB)
97: addl %esi, %edi
98: adcl %eax, %ecx
99: adcl %edx, $0
100:
101: / Normalize the result mantissa.
102: / The product presently contains 73 or 74 bits.
103: / The normalized result contains 53+32 = 85 bits.
104: movb %cl, $12 / shift count to CL
105: testl %edx, $0x200 / examine product bit 74
106: jz ?L0 / 73 bit product
107: decb %cl / 74 bit product, decrement shift count
108: incl %ebp / and bump the exponent
109:
110: ?L0:
111: / Shift EDX:EAX:ESI left by CL bits.
112: shld %edx, %eax, %cl
113: shld %eax, %esi, %cl
114: shll %esi, %cl
115:
116: / Round up result if appropriate.
117: shll %esi, $1 / high bit of lo 64 bits to CF
118: ?addcarry:
119: adcl %eax, $0
120: adcl %edx, $0
121: testl %edx, $HIDDEN<<1
122: jnz ?carry
123:
124: / Pack result mantissa from EDX:EAX, exponent from EBP, sign from stack.
125: ?pack:
126: andl %edx, $MANMASK / mask off hidden bit
127: orl %ebp, %ebp
128: jle ?zero / exponent underflow, return 0.0
129: cmp %ebp, $EXPMASK>>20
130: jge ?inf / exponent overflow, return infinity
131: shll %ebp, $20 / position exponent
132: orl %edx, %ebp / pack mantissa and exponent
133: pop %ecx
134: orl %edx, %ecx / pack with sign
135:
136: / Return the value from EDX:EAX.
137: ?done:
138: pop %ecx
139: pop %ebx
140: pop %edi
141: pop %esi
142: pop %ebp
143: ret
144:
145: / Rounding generated carry past hidden bit, shift one bit right.
146: ?carry:
147: incl %ebp / bump the exponent
148: shrd %eax, %edx, $1
149: pushfl
150: shrl %edx, $1
151: popfl
152: jmp ?addcarry
153:
154: / Return the value from ESI:EDI.
155: ?retlhs:
156: xchgl %esi, %edx
157: xchgl %edi, %eax
158: jmp ?done
159:
160: ?inf:
161: pop %edx / pop result sign bit
162: orl %edx, $EXPMASK / max exp, zero mantissa for infinity
163: jmp ?L2
164:
165: ?zero:
166: pop %edx / pop sign bit and ignore
167: subl %edx, %edx / return 0.0
168: ?L2:
169: subl %eax, %eax
170: jmp ?done
171:
172: / end of libc/crt/i386/dmul.s
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.