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