|
|
1.1 root 1: #
2: # atof, strtod: convert ascii to floating
3: #
4: # C usage:
5: #
6: # double atof (s)
7: # char *s;
8: #
9: # double strtod (s, p)
10: # char *s, **p;
11: #
12: # For strtod, if p is non-null, *p gets set to point at the
13: # character after the input string.
14: #
15: # Register usage:
16: #
17: # r0-1: value being developed
18: # r2: first section: pointer to the next character
19: # second section: binary exponent
20: # r3: flags
21: # r4: first section: the current character
22: # second section: scratch
23: # r5: the decimal exponent
24: # r6-7: scratch
25: # r8: a copy of p (strtod) or zero (atof)
26: #
27: # Flag definitions
28: #
29: .set msign,0 # mantissa has negative sign
30: .set esign,1 # exponent has negative sign
31: .set decpt,2 # decimal point encountered
32:
33: .align 2
34: two31: .word 0x5000 # 2 ** 31
35: .word 0 # (=2147483648)
36: .word 0 # in floating-point
37: .word 0 # (so atof doesn't have to convert it)
38: #
39: # Entry points
40: #
41: .text
42: .align 2
43: .globl _strtod
44: _strtod:.word 0x01c0 # Save r8 through r6 (we use r0-r8)
45: movl 8(ap),r8 # Pick up value of p
46: jbr init
47:
48: .align 2
49: .globl _atof
50: _atof: .word 0x01c0 # Save r8 through r6 (we use r0-r8)
51: clrl r8 # This is atof, not strtod
52: #
53: # Initialization
54: #
55: init:
56: clrl r3 # All flags start out false
57: movl 4(ap),r2 # Address the first character
58: clrl r5 # Clear starting exponent
59: #
60: # Skip leading white space
61: #
62: sk0: movzbl (r2)+,r4 # Fetch the next (first) character
63: cmpb $' ,r4 # Is it blank?
64: jeql sk0 # ...yes
65: cmpb r4,$8 # 8 is lowest of white-space group
66: jlss sk1 # Jump if char too low to be white space
67: cmpb r4,$13 # 13 is highest of white-space group
68: jleq sk0 # Jump if character is white space
69: sk1:
70: #
71: # Check for a sign
72: #
73: cmpb $'+,r4 # Positive sign?
74: jeql cs1 # ... yes
75: cmpb $'-,r4 # Negative sign?
76: jneq cs2 # ... no
77: bisb2 $1<msign,r3 # Indicate a negative mantissa
78: cs1: movzbl (r2)+,r4 # Skip the character
79: cs2:
80: #
81: # Accumulate digits, keeping track of the exponent
82: #
83: clrq r0 # Clear the accumulator
84: ad0: cmpb r4,$'0 # Do we have a digit?
85: jlss ad4 # ... no, too small
86: cmpb r4,$'9
87: jgtr ad4 # ... no, too large
88: #
89: # We got a digit. Accumulate it
90: #
91: cmpl r1,$214748364 # Would this digit cause overflow?
92: jgeq ad1 # ... yes
93: #
94: # Multiply (r0,r1) by 10. This is done by developing
95: # (r0,r1)*2 in (r6,r7), shifting (r0,r1) left three bits,
96: # and adding the two quadwords.
97: #
98: ashq $1,r0,r6 # (r6,r7)=(r0,r1)*2
99: ashq $3,r0,r0 # (r0,r1)=(r0,r1)*8
100: addl2 r6,r0 # Add low halves
101: adwc r7,r1 # Add high halves
102: #
103: # Add in the digit
104: #
105: subl2 $'0,r4 # Get the digit value
106: addl2 r4,r0 # Add it into the accumulator
107: adwc $0,r1 # Possible carry into high half
108: jbr ad2 # Join common code
109: #
110: # Here when the digit won't fit in the accumulator
111: #
112: ad1: incl r5 # Ignore the digit, bump exponent
113: #
114: # If we have seen a decimal point, decrease the exponent by 1
115: #
116: ad2: jbc $decpt,r3,ad3 # Jump if decimal point not seen
117: decl r5 # Decrease exponent
118: ad3:
119: #
120: # Fetch the next character, back for more
121: #
122: movzbl (r2)+,r4 # Fetch
123: jbr ad0 # Try again
124: #
125: # Not a digit. Could it be a decimal point?
126: #
127: ad4: cmpb r4,$'. # If it's not a decimal point, either it's
128: jneq ad5 # the end of the number or the start of
129: # the exponent.
130: jbcs $decpt,r3,ad3 # If it IS a decimal point, we record that
131: # we've seen one, and keep collecting
132: # digits if it is the first one.
133: #
134: # Check for an exponent
135: #
136: ad5: clrl r6 # Initialize the exponent accumulator
137:
138: cmpb r4,$'e # We allow both lower case e
139: jeql ex1 # ... and ...
140: cmpb r4,$'E # upper-case E
141: jneq ex7
142: #
143: # Does the exponent have a sign?
144: #
145: ex1: movzbl (r2)+,r4 # Get next character
146: cmpb r4,$'+ # Positive sign?
147: jeql ex2 # ... yes ...
148: cmpb r4,$'- # Negative sign?
149: jneq ex3 # ... no ...
150: bisb2 $1<esign,r3 # Indicate exponent is negative
151: ex2: movzbl (r2)+,r4 # Grab the next character
152: #
153: # Accumulate exponent digits in r6
154: #
155: ex3: cmpb r4,$'0 # A digit is within the range
156: jlss ex4 # '0' through
157: cmpb r4,$'9 # '9',
158: jgtr ex4 # inclusive.
159: cmpl r6,$214748364 # Exponent outrageously large already?
160: jgeq ex2 # ... yes
161: moval (r6)[r6],r6 # r6 *= 5
162: movaw -'0(r4)[r6],r6 # r6 = r6 * 2 + r4 - '0'
163: jbr ex2 # Go 'round again
164: ex4:
165: #
166: # Now get the final exponent and force it within a reasonable
167: # range so our scaling loops don't take forever for values
168: # that will ultimately cause overflow or underflow anyway.
169: # A tight check on over/underflow will be done by ldexp.
170: #
171: jbc $esign,r3,ex5 # Jump if exponent not negative
172: mnegl r6,r6 # If sign, negate exponent
173: ex5: addl2 r6,r5 # Add given exponent to calculated exponent
174: cmpl r5,$-100 # Absurdly small?
175: jgtr ex6 # ... no
176: movl $-100,r5 # ... yes, force within limit
177: ex6: cmpl r5,$100 # Absurdly large?
178: jlss ex7 # ... no
179: movl $100,r5 # ... yes, force within bounds
180: ex7:
181: #
182: # If this is strtod and p is non-null, store a value in *p
183: #
184: tstl r8 # Is there a place to store it?
185: jeql ex8 # ... no, skip store
186: decl r2
187: movl r2,(r8) # ... yes, store value in *p
188: ex8:
189: #
190: # Our number has now been reduced to a mantissa and an exponent.
191: # The mantissa is a 63-bit positive binary integer in r0,r1,
192: # and the exponent is a signed power of 10 in r5. The msign
193: # bit in r3 will be on if the mantissa should ultimately be
194: # considered negative.
195: #
196: # We now have to convert it to a standard format floating point
197: # number. This will be done by accumulating a binary exponent
198: # in r2, as we progressively get r5 closer to zero.
199: #
200: # Don't bother scaling if the mantissa is zero
201: #
202: movq r0,r0 # Mantissa zero?
203: jeql exit # ... yes
204:
205: clrl r2 # Initialize binary exponent
206: tstl r5 # Which way to scale?
207: jleq sd0 # Scale down if decimal exponent <= 0
208: #
209: # Scale up by "multiplying" r0,r1 by 10 as many times as necessary,
210: # as follows:
211: #
212: # Step 1: Shift r0,r1 right as necessary to ensure that no
213: # overflow can occur when multiplying.
214: #
215: su0: cmpl r1,$429496729 # Compare high word to (2**31)/5
216: jlss su1 # Jump out if guaranteed safe
217: ashq $-1,r0,r0 # Else shift right one bit
218: incl r2 # bump exponent to compensate
219: jbr su0 # and go back to test again.
220: #
221: # Step 2: Multiply r0,r1 by 5, by appropriate shifting and
222: # double-precision addition
223: #
224: su1: ashq $2,r0,r6 # (r6,r7) := (r0,r1) * 4
225: addl2 r6,r0 # Add low-order halves
226: adwc r7,r1 # and high-order halves
227: #
228: # Step 3: Increment the binary exponent to take care of the final
229: # factor of 2, and go back if we still need to scale more.
230: #
231: incl r2 # Increment the exponent
232: sobgtr r5,su0 # and back for more (maybe)
233:
234: jbr cm0 # Merge to build final value
235:
236: #
237: # Scale down. We must "divide" r0,r1 by 10 as many times
238: # as needed, as follows:
239: #
240: # Step 0: Right now, the condition codes reflect the state
241: # of r5. If it's zero, we are done.
242: #
243: sd0: jeql cm0 # If finished, build final number
244: #
245: # Step 1: Shift r0,r1 left until the high-order bit (not counting
246: # the sign bit) is nonzero, so that the division will preserve
247: # as much precision as possible.
248: #
249: tstl r1 # Is the entire high-order half zero?
250: jneq sd2 # ...no, go shift one bit at a time
251: ashq $30,r0,r0 # ...yes, shift left 30,
252: subl2 $30,r2 # decrement the exponent to compensate,
253: # and now it's known to be safe to shift
254: # at least once more.
255: sd1: ashq $1,r0,r0 # Shift (r0,r1) left one, and
256: decl r2 # decrement the exponent to compensate
257: sd2: jbc $30,r1,sd1 # If the high-order bit is off, go shift
258: #
259: # Step 2: Divide the high-order part of (r0,r1) by 5,
260: # giving a quotient in r1 and a remainder in r7.
261: #
262: sd3: movl r1,r6 # Copy the high-order part
263: clrl r7 # Zero-extend to 64 bits
264: ediv $5,r6,r1,r7 # Divide (cannot overflow)
265: #
266: # Step 3: Divide the low-order part of (r0,r1) by 5,
267: # using the remainder from step 2 for rounding.
268: # Note that the result of this computation is unsigned,
269: # so we have to allow for the fact that an ordinary division
270: # by 5 could overflow. We make allowance by dividing by 10,
271: # multiplying the quotient by 2, and using the remainder
272: # to adjust the modified quotient.
273: #
274: addl3 $2,r0,r6 # Dividend is low part of (r0,r1) plus
275: adwc $0,r7 # 2 for rounding plus
276: # (2**32) * previous remainder
277: ediv $10,r6,r0,r6 # r0 := quotient, r6 := remainder.
278: addl2 r0,r0 # Make r0 result of dividing by 5
279: cmpl r6,$5 # If remainder is 5 or greater,
280: jlss sd4 # increment the adjustted quotient.
281: incl r0
282: #
283: # Step 4: Increment the decimal exponent, decrement the binary
284: # exponent (to make the division by 5 into a division by 10),
285: # and back for another iteration.
286: #
287: sd4: decl r2 # Binary exponent
288: aoblss $0,r5,sd2
289: #
290: # We now have the following:
291: #
292: # r0: low-order half of a 64-bit integer
293: # r1: high-order half of the same 64-bit integer
294: # r2: a binary exponent
295: #
296: # Our final result is the integer represented by (r0,r1)
297: # multiplied by 2 to the power contained in r2.
298: # We will transform (r0,r1) into a floating-point value,
299: # set the sign appropriately, and let ldexp do the
300: # rest of the work.
301: #
302: # Step 1: if the high-order bit (excluding the sign) of
303: # the high-order half (r1) is 1, then we have 63 bits of
304: # fraction, too many to convert easily. However, we also
305: # know we won't need them all, so we will just throw the
306: # low-order bit away (and adjust the exponent appropriately).
307: #
308: cm0: jbc $30,r1,cm1 # jump if no adjustment needed
309: ashq $-1,r0,r0 # lose the low-order bit
310: incl r2 # increase the exponent to compensate
311: #
312: # Step 2: split the 62-bit number in (r0,r1) into two
313: # 31-bit positive quantities
314: #
315: cm1: ashq $1,r0,r0 # put the high-order bits in r1
316: # and a 0 in the bottom of r0
317: rotl $-1,r0,r0 # right-justify the bits in r0
318: # moving the 0 from the ashq
319: # into the sign bit.
320: #
321: # Step 3: convert both halves to floating point
322: #
323: cvtld r0,r6 # low-order part in r6-r7
324: cvtld r1,r0 # high-order part in r0-r1
325: #
326: # Step 4: multiply the high order part by 2**31 and combine them
327: #
328: muld2 two31,r0 # multiply
329: addd2 r6,r0 # combine
330: #
331: # Step 5: if appropriate, negate the floating value
332: #
333: jbc $msign,r3,cm2 # Jump if mantissa not signed
334: mnegd r0,r0 # If negative, make it so
335: #
336: # Step 6: call ldexp to complete the job
337: #
338: cm2: pushl r2 # Put exponent in parameter list
339: movd r0,-(sp) # and also mantissa
340: calls $3,_ldexp # go combine them
341:
342: exit: ret
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.