| blob | e88c7c15ad7c03e347afb7def4cea273a208a482 |
1 * $NetBSD: decbin.sa,v 1.4 2001/12/09 01:43:13 briggs Exp $
3 * MOTOROLA MICROPROCESSOR & MEMORY TECHNOLOGY GROUP
4 * M68000 Hi-Performance Microprocessor Division
5 * M68040 Software Package
6 *
7 * M68040 Software Package Copyright (c) 1993, 1994 Motorola Inc.
8 * All rights reserved.
9 *
10 * THE SOFTWARE is provided on an "AS IS" basis and without warranty.
11 * To the maximum extent permitted by applicable law,
12 * MOTOROLA DISCLAIMS ALL WARRANTIES WHETHER EXPRESS OR IMPLIED,
13 * INCLUDING IMPLIED WARRANTIES OF MERCHANTABILITY OR FITNESS FOR A
14 * PARTICULAR PURPOSE and any warranty against infringement with
15 * regard to the SOFTWARE (INCLUDING ANY MODIFIED VERSIONS THEREOF)
16 * and any accompanying written materials.
17 *
18 * To the maximum extent permitted by applicable law,
19 * IN NO EVENT SHALL MOTOROLA BE LIABLE FOR ANY DAMAGES WHATSOEVER
20 * (INCLUDING WITHOUT LIMITATION, DAMAGES FOR LOSS OF BUSINESS
21 * PROFITS, BUSINESS INTERRUPTION, LOSS OF BUSINESS INFORMATION, OR
22 * OTHER PECUNIARY LOSS) ARISING OF THE USE OR INABILITY TO USE THE
23 * SOFTWARE. Motorola assumes no responsibility for the maintenance
24 * and support of the SOFTWARE.
25 *
26 * You are hereby granted a copyright license to use, modify, and
27 * distribute the SOFTWARE so long as this entire notice is retained
28 * without alteration in any modified and/or redistributed versions,
29 * and that such modified versions are clearly identified as such.
30 * No licenses are granted by implication, estoppel or otherwise
31 * under any patents or trademarks of Motorola, Inc.
33 *
34 * decbin.sa 3.3 12/19/90
35 *
36 * Description: Converts normalized packed bcd value pointed to by
37 * register A6 to extended-precision value in FP0.
38 *
39 * Input: Normalized packed bcd value in ETEMP(a6).
40 *
41 * Output: Exact floating-point representation of the packed bcd value.
42 *
43 * Saves and Modifies: D2-D5
44 *
45 * Speed: The program decbin takes ??? cycles to execute.
46 *
47 * Object Size:
48 *
49 * External Reference(s): None.
50 *
51 * Algorithm:
52 * Expected is a normal bcd (i.e. non-exceptional; all inf, zero,
53 * and NaN operands are dispatched without entering this routine)
54 * value in 68881/882 format at location ETEMP(A6).
55 *
56 * A1. Convert the bcd exponent to binary by successive adds and muls.
57 * Set the sign according to SE. Subtract 16 to compensate
58 * for the mantissa which is to be interpreted as 17 integer
59 * digits, rather than 1 integer and 16 fraction digits.
60 * Note: this operation can never overflow.
61 *
62 * A2. Convert the bcd mantissa to binary by successive
63 * adds and muls in FP0. Set the sign according to SM.
64 * The mantissa digits will be converted with the decimal point
65 * assumed following the least-significant digit.
66 * Note: this operation can never overflow.
67 *
68 * A3. Count the number of leading/trailing zeros in the
69 * bcd string. If SE is positive, count the leading zeros;
70 * if negative, count the trailing zeros. Set the adjusted
71 * exponent equal to the exponent from A1 and the zero count
72 * added if SM = 1 and subtracted if SM = 0. Scale the
73 * mantissa the equivalent of forcing in the bcd value:
74 *
75 * SM = 0 a non-zero digit in the integer position
76 * SM = 1 a non-zero digit in Mant0, lsd of the fraction
77 *
78 * this will insure that any value, regardless of its
79 * representation (ex. 0.1E2, 1E1, 10E0, 100E-1), is converted
80 * consistently.
81 *
82 * A4. Calculate the factor 10^exp in FP1 using a table of
83 * 10^(2^n) values. To reduce the error in forming factors
84 * greater than 10^27, a directed rounding scheme is used with
85 * tables rounded to RN, RM, and RP, according to the table
86 * in the comments of the pwrten section.
87 *
88 * A5. Form the final binary number by scaling the mantissa by
89 * the exponent factor. This is done by multiplying the
90 * mantissa in FP0 by the factor in FP1 if the adjusted
91 * exponent sign is positive, and dividing FP0 by FP1 if
92 * it is negative.
93 *
94 * Clean up and return. Check if the final mul or div resulted
95 * in an inex2 exception. If so, set inex1 in the fpsr and
96 * check if the inex1 exception is enabled. If so, set d7 upper
97 * word to $0100. This will signal unimp.sa that an enabled inex1
98 * exception occurred. Unimp will fix the stack.
99 *
101 DECBIN IDNT 2,1 Motorola 040 Floating Point Software Package
103 section 8
105 include fpsp.h
107 *
108 * PTENRN, PTENRM, and PTENRP are arrays of powers of 10 rounded
109 * to nearest, minus, and plus, respectively. The tables include
110 * 10**{1,2,4,8,16,32,64,128,256,512,1024,2048,4096}. No rounding
111 * is required until the power is greater than 27, however, all
112 * tables include the first 5 for ease of indexing.
113 *
114 xref PTENRN
115 xref PTENRM
116 xref PTENRP
118 RTABLE dc.b 0,0,0,0
119 dc.b 2,3,2,3
120 dc.b 2,3,3,2
121 dc.b 3,2,2,3
123 xdef decbin
124 xdef calc_e
125 xdef pwrten
126 xdef calc_m
127 xdef norm
128 xdef ap_st_z
129 xdef ap_st_n
130 *
131 FNIBS equ 7
132 FSTRT equ 0
133 *
134 ESTRT equ 4
135 EDIGITS equ 2
136 *
137 * Constants in single precision
138 FZERO dc.l $00000000
139 FONE dc.l $3F800000
140 FTEN dc.l $41200000
142 TEN equ 10
144 *
145 decbin:
146 fmove.l #0,FPCR ;clr real fpcr
147 movem.l d2-d5,-(a7)
148 *
149 * Calculate exponent:
150 * 1. Copy bcd value in memory for use as a working copy.
151 * 2. Calculate absolute value of exponent in d1 by mul and add.
152 * 3. Correct for exponent sign.
153 * 4. Subtract 16 to compensate for interpreting the mant as all integer digits.
154 * (i.e., all digits assumed left of the decimal point.)
155 *
156 * Register usage:
157 *
158 * calc_e:
159 * (*) d0: temp digit storage
160 * (*) d1: accumulator for binary exponent
161 * (*) d2: digit count
162 * (*) d3: offset pointer
163 * ( ) d4: first word of bcd
164 * ( ) a0: pointer to working bcd value
165 * ( ) a6: pointer to original bcd value
166 * (*) FP_SCR1: working copy of original bcd value
167 * (*) L_SCR1: copy of original exponent word
168 *
169 calc_e:
170 move.l #EDIGITS,d2 ;# of nibbles (digits) in fraction part
171 moveq.l #ESTRT,d3 ;counter to pick up digits
172 lea.l FP_SCR1(a6),a0 ;load tmp bcd storage address
173 move.l ETEMP(a6),(a0) ;save input bcd value
174 move.l ETEMP_HI(a6),4(a0) ;save words 2 and 3
175 move.l ETEMP_LO(a6),8(a0) ;and work with these
176 move.l (a0),d4 ;get first word of bcd
177 clr.l d1 ;zero d1 for accumulator
178 e_gd:
179 mulu.l #TEN,d1 ;mul partial product by one digit place
180 bfextu d4{d3:4},d0 ;get the digit and zero extend into d0
181 add.l d0,d1 ;d1 = d1 + d0
182 addq.b #4,d3 ;advance d3 to the next digit
183 dbf.w d2,e_gd ;if we have used all 3 digits, exit loop
184 btst #30,d4 ;get SE
185 beq.b e_pos ;don't negate if pos
186 neg.l d1 ;negate before subtracting
187 e_pos:
188 sub.l #16,d1 ;sub to compensate for shift of mant
189 bge.b e_save ;if still pos, do not neg
190 neg.l d1 ;now negative, make pos and set SE
191 or.l #$40000000,d4 ;set SE in d4,
192 or.l #$40000000,(a0) ;and in working bcd
193 e_save:
194 move.l d1,L_SCR1(a6) ;save exp in memory
195 *
196 *
197 * Calculate mantissa:
198 * 1. Calculate absolute value of mantissa in fp0 by mul and add.
199 * 2. Correct for mantissa sign.
200 * (i.e., all digits assumed left of the decimal point.)
201 *
202 * Register usage:
203 *
204 * calc_m:
205 * (*) d0: temp digit storage
206 * (*) d1: lword counter
207 * (*) d2: digit count
208 * (*) d3: offset pointer
209 * ( ) d4: words 2 and 3 of bcd
210 * ( ) a0: pointer to working bcd value
211 * ( ) a6: pointer to original bcd value
212 * (*) fp0: mantissa accumulator
213 * ( ) FP_SCR1: working copy of original bcd value
214 * ( ) L_SCR1: copy of original exponent word
215 *
216 calc_m:
217 moveq.l #1,d1 ;word counter, init to 1
218 fmove.s FZERO,fp0 ;accumulator
219 *
220 *
221 * Since the packed number has a long word between the first & second parts,
222 * get the integer digit then skip down & get the rest of the
223 * mantissa. We will unroll the loop once.
224 *
225 bfextu (a0){28:4},d0 ;integer part is ls digit in long word
226 fadd.b d0,fp0 ;add digit to sum in fp0
227 *
228 *
229 * Get the rest of the mantissa.
230 *
231 loadlw:
232 move.l (a0,d1.L*4),d4 ;load mantissa lonqword into d4
233 moveq.l #FSTRT,d3 ;counter to pick up digits
234 moveq.l #FNIBS,d2 ;reset number of digits per a0 ptr
235 md2b:
236 fmul.s FTEN,fp0 ;fp0 = fp0 * 10
237 bfextu d4{d3:4},d0 ;get the digit and zero extend
238 fadd.b d0,fp0 ;fp0 = fp0 + digit
239 *
240 *
241 * If all the digits (8) in that long word have been converted (d2=0),
242 * then inc d1 (=2) to point to the next long word and reset d3 to 0
243 * to initialize the digit offset, and set d2 to 7 for the digit count;
244 * else continue with this long word.
245 *
246 addq.b #4,d3 ;advance d3 to the next digit
247 dbf.w d2,md2b ;check for last digit in this lw
248 nextlw:
249 addq.l #1,d1 ;inc lw pointer in mantissa
250 cmp.l #2,d1 ;test for last lw
251 ble loadlw ;if not, get last one
253 *
254 * Check the sign of the mant and make the value in fp0 the same sign.
255 *
256 m_sign:
257 btst #31,(a0) ;test sign of the mantissa
258 beq.b short_ap_st_z ;if clear, go to append/strip zeros
259 fneg.x fp0 ;if set, negate fp0
261 *
262 * Append/strip zeros:
263 *
264 * For adjusted exponents which have an absolute value greater than 27*,
265 * this routine calculates the amount needed to normalize the mantissa
266 * for the adjusted exponent. That number is subtracted from the exp
267 * if the exp was positive, and added if it was negative. The purpose
268 * of this is to reduce the value of the exponent and the possibility
269 * of error in calculation of pwrten.
270 *
271 * 1. Branch on the sign of the adjusted exponent.
272 * 2p.(positive exp)
273 * 2. Check M16 and the digits in lwords 2 and 3 in decending order.
274 * 3. Add one for each zero encountered until a non-zero digit.
275 * 4. Subtract the count from the exp.
276 * 5. Check if the exp has crossed zero in #3 above; make the exp abs
277 * and set SE.
278 * 6. Multiply the mantissa by 10**count.
279 * 2n.(negative exp)
280 * 2. Check the digits in lwords 3 and 2 in decending order.
281 * 3. Add one for each zero encountered until a non-zero digit.
282 * 4. Add the count to the exp.
283 * 5. Check if the exp has crossed zero in #3 above; clear SE.
284 * 6. Divide the mantissa by 10**count.
285 *
286 * *Why 27? If the adjusted exponent is within -28 < expA < 28, than
287 * any adjustment due to append/strip zeros will drive the resultane
288 * exponent towards zero. Since all pwrten constants with a power
289 * of 27 or less are exact, there is no need to use this routine to
290 * attempt to lessen the resultant exponent.
291 *
292 * Register usage:
293 *
294 * ap_st_z:
295 * (*) d0: temp digit storage
296 * (*) d1: zero count
297 * (*) d2: digit count
298 * (*) d3: offset pointer
299 * ( ) d4: first word of bcd
300 * (*) d5: lword counter
301 * ( ) a0: pointer to working bcd value
302 * ( ) FP_SCR1: working copy of original bcd value
303 * ( ) L_SCR1: copy of original exponent word
304 *
305 *
306 * First check the absolute value of the exponent to see if this
307 * routine is necessary. If so, then check the sign of the exponent
308 * and do append (+) or strip (-) zeros accordingly.
309 * This section handles a positive adjusted exponent.
310 *
311 ap_st_z:
312 short_ap_st_z:
313 move.l L_SCR1(a6),d1 ;load expA for range test
314 cmp.l #27,d1 ;test is with 27
315 ble.w pwrten ;if abs(expA) <28, skip ap/st zeros
316 btst #30,(a0) ;check sign of exp
317 bne.b short_ap_st_n ;if neg, go to neg side
318 clr.l d1 ;zero count reg
319 move.l (a0),d4 ;load lword 1 to d4
320 bfextu d4{28:4},d0 ;get M16 in d0
321 bne.b ap_p_fx ;if M16 is non-zero, go fix exp
322 addq.l #1,d1 ;inc zero count
323 moveq.l #1,d5 ;init lword counter
324 move.l (a0,d5.L*4),d4 ;get lword 2 to d4
325 bne.b ap_p_cl ;if lw 2 is zero, skip it
326 addq.l #8,d1 ;and inc count by 8
327 addq.l #1,d5 ;inc lword counter
328 move.l (a0,d5.L*4),d4 ;get lword 3 to d4
329 ap_p_cl:
330 clr.l d3 ;init offset reg
331 moveq.l #7,d2 ;init digit counter
332 ap_p_gd:
333 bfextu d4{d3:4},d0 ;get digit
334 bne.b ap_p_fx ;if non-zero, go to fix exp
335 addq.l #4,d3 ;point to next digit
336 addq.l #1,d1 ;inc digit counter
337 dbf.w d2,ap_p_gd ;get next digit
338 ap_p_fx:
339 move.l d1,d0 ;copy counter to d2
340 move.l L_SCR1(a6),d1 ;get adjusted exp from memory
341 sub.l d0,d1 ;subtract count from exp
342 bge.b ap_p_fm ;if still pos, go to pwrten
343 neg.l d1 ;now its neg; get abs
344 move.l (a0),d4 ;load lword 1 to d4
345 or.l #$40000000,d4 ; and set SE in d4
346 or.l #$40000000,(a0) ; and in memory
347 *
348 * Calculate the mantissa multiplier to compensate for the striping of
349 * zeros from the mantissa.
350 *
351 ap_p_fm:
352 move.l #PTENRN,a1 ;get address of power-of-ten table
353 clr.l d3 ;init table index
354 fmove.s FONE,fp1 ;init fp1 to 1
355 moveq.l #3,d2 ;init d2 to count bits in counter
356 ap_p_el:
357 asr.l #1,d0 ;shift lsb into carry
358 bcc.b ap_p_en ;if 1, mul fp1 by pwrten factor
359 fmul.x (a1,d3),fp1 ;mul by 10**(d3_bit_no)
360 ap_p_en:
361 add.l #12,d3 ;inc d3 to next rtable entry
362 tst.l d0 ;check if d0 is zero
363 bne.b ap_p_el ;if not, get next bit
364 fmul.x fp1,fp0 ;mul mantissa by 10**(no_bits_shifted)
365 bra.b short_pwrten ;go calc pwrten
366 *
367 * This section handles a negative adjusted exponent.
368 *
369 ap_st_n:
370 short_ap_st_n:
371 clr.l d1 ;clr counter
372 moveq.l #2,d5 ;set up d5 to point to lword 3
373 move.l (a0,d5.L*4),d4 ;get lword 3
374 bne.b ap_n_cl ;if not zero, check digits
375 sub.l #1,d5 ;dec d5 to point to lword 2
376 addq.l #8,d1 ;inc counter by 8
377 move.l (a0,d5.L*4),d4 ;get lword 2
378 ap_n_cl:
379 move.l #28,d3 ;point to last digit
380 moveq.l #7,d2 ;init digit counter
381 ap_n_gd:
382 bfextu d4{d3:4},d0 ;get digit
383 bne.b ap_n_fx ;if non-zero, go to exp fix
384 subq.l #4,d3 ;point to previous digit
385 addq.l #1,d1 ;inc digit counter
386 dbf.w d2,ap_n_gd ;get next digit
387 ap_n_fx:
388 move.l d1,d0 ;copy counter to d0
389 move.l L_SCR1(a6),d1 ;get adjusted exp from memory
390 sub.l d0,d1 ;subtract count from exp
391 bgt.b ap_n_fm ;if still pos, go fix mantissa
392 neg.l d1 ;take abs of exp and clr SE
393 move.l (a0),d4 ;load lword 1 to d4
394 and.l #$bfffffff,d4 ; and clr SE in d4
395 and.l #$bfffffff,(a0) ; and in memory
396 *
397 * Calculate the mantissa multiplier to compensate for the appending of
398 * zeros to the mantissa.
399 *
400 ap_n_fm:
401 move.l #PTENRN,a1 ;get address of power-of-ten table
402 clr.l d3 ;init table index
403 fmove.s FONE,fp1 ;init fp1 to 1
404 moveq.l #3,d2 ;init d2 to count bits in counter
405 ap_n_el:
406 asr.l #1,d0 ;shift lsb into carry
407 bcc.b ap_n_en ;if 1, mul fp1 by pwrten factor
408 fmul.x (a1,d3),fp1 ;mul by 10**(d3_bit_no)
409 ap_n_en:
410 add.l #12,d3 ;inc d3 to next rtable entry
411 tst.l d0 ;check if d0 is zero
412 bne.b ap_n_el ;if not, get next bit
413 fdiv.x fp1,fp0 ;div mantissa by 10**(no_bits_shifted)
414 *
415 *
416 * Calculate power-of-ten factor from adjusted and shifted exponent.
417 *
418 * Register usage:
419 *
420 * pwrten:
421 * (*) d0: temp
422 * ( ) d1: exponent
423 * (*) d2: {FPCR[6:5],SM,SE} as index in RTABLE; temp
424 * (*) d3: FPCR work copy
425 * ( ) d4: first word of bcd
426 * (*) a1: RTABLE pointer
427 * calc_p:
428 * (*) d0: temp
429 * ( ) d1: exponent
430 * (*) d3: PWRTxx table index
431 * ( ) a0: pointer to working copy of bcd
432 * (*) a1: PWRTxx pointer
433 * (*) fp1: power-of-ten accumulator
434 *
435 * Pwrten calculates the exponent factor in the selected rounding mode
436 * according to the following table:
437 *
438 * Sign of Mant Sign of Exp Rounding Mode PWRTEN Rounding Mode
439 *
440 * ANY ANY RN RN
441 *
442 * + + RP RP
443 * - + RP RM
444 * + - RP RM
445 * - - RP RP
446 *
447 * + + RM RM
448 * - + RM RP
449 * + - RM RP
450 * - - RM RM
451 *
452 * + + RZ RM
453 * - + RZ RM
454 * + - RZ RP
455 * - - RZ RP
456 *
457 *
458 pwrten:
459 short_pwrten:
460 move.l USER_FPCR(a6),d3 ;get user's FPCR
461 bfextu d3{26:2},d2 ;isolate rounding mode bits
462 move.l (a0),d4 ;reload 1st bcd word to d4
463 asl.l #2,d2 ;format d2 to be
464 bfextu d4{0:2},d0 ; {FPCR[6],FPCR[5],SM,SE}
465 add.l d0,d2 ;in d2 as index into RTABLE
466 lea.l RTABLE,a1 ;load rtable base
467 move.b (a1,d2),d0 ;load new rounding bits from table
468 clr.l d3 ;clear d3 to force no exc and extended
469 bfins d0,d3{26:2} ;stuff new rounding bits in FPCR
470 fmove.l d3,FPCR ;write new FPCR
471 asr.l #1,d0 ;write correct PTENxx table
472 bcc.b not_rp ;to a1
473 lea.l PTENRP,a1 ;it is RP
474 bra.b calc_p ;go to init section
475 not_rp:
476 asr.l #1,d0 ;keep checking
477 bcc.b not_rm
478 lea.l PTENRM,a1 ;it is RM
479 bra.b calc_p ;go to init section
480 not_rm:
481 lea.l PTENRN,a1 ;it is RN
482 calc_p:
483 move.l d1,d0 ;copy exp to d0;use d0
484 bpl.b no_neg ;if exp is negative,
485 neg.l d0 ;invert it
486 or.l #$40000000,(a0) ;and set SE bit
487 no_neg:
488 clr.l d3 ;table index
489 fmove.s FONE,fp1 ;init fp1 to 1
490 e_loop:
491 asr.l #1,d0 ;shift next bit into carry
492 bcc.b e_next ;if zero, skip the mul
493 fmul.x (a1,d3),fp1 ;mul by 10**(d3_bit_no)
494 e_next:
495 add.l #12,d3 ;inc d3 to next rtable entry
496 tst.l d0 ;check if d0 is zero
497 bne.b e_loop ;not zero, continue shifting
498 *
499 *
500 * Check the sign of the adjusted exp and make the value in fp0 the
501 * same sign. If the exp was pos then multiply fp1*fp0;
502 * else divide fp0/fp1.
503 *
504 * Register Usage:
505 * norm:
506 * ( ) a0: pointer to working bcd value
507 * (*) fp0: mantissa accumulator
508 * ( ) fp1: scaling factor - 10**(abs(exp))
509 *
510 norm:
511 btst #30,(a0) ;test the sign of the exponent
512 beq.b mul ;if clear, go to multiply
513 div:
514 fdiv.x fp1,fp0 ;exp is negative, so divide mant by exp
515 bra.b end_dec
516 mul:
517 fmul.x fp1,fp0 ;exp is positive, so multiply by exp
518 *
519 *
520 * Clean up and return with result in fp0.
521 *
522 * If the final mul/div in decbin incurred an inex exception,
523 * it will be inex2, but will be reported as inex1 by get_op.
524 *
525 end_dec:
526 fmove.l FPSR,d0 ;get status register
527 bclr.l #inex2_bit+8,d0 ;test for inex2 and clear it
528 fmove.l d0,FPSR ;return status reg w/o inex2
529 beq.b no_exc ;skip this if no exc
530 or.l #inx1a_mask,USER_FPSR(a6) ;set inex1/ainex
531 no_exc:
532 movem.l (a7)+,d2-d5
533 rts
534 end