|
|
1.1 root 1: /*
2: * Copyright (c) 1999 Apple Computer, Inc. All rights reserved.
3: *
4: * @APPLE_LICENSE_HEADER_START@
5: *
6: * Portions Copyright (c) 1999 Apple Computer, Inc. All Rights
7: * Reserved. This file contains Original Code and/or Modifications of
8: * Original Code as defined in and that are subject to the Apple Public
9: * Source License Version 1.1 (the "License"). You may not use this file
10: * except in compliance with the License. Please obtain a copy of the
11: * License at http://www.apple.com/publicsource and read it before using
12: * this file.
13: *
14: * The Original Code and all software distributed under the License are
15: * distributed on an "AS IS" basis, WITHOUT WARRANTY OF ANY KIND, EITHER
16: * EXPRESS OR IMPLIED, AND APPLE HEREBY DISCLAIMS ALL SUCH WARRANTIES,
17: * INCLUDING WITHOUT LIMITATION, ANY WARRANTIES OF MERCHANTABILITY,
18: * FITNESS FOR A PARTICULAR PURPOSE OR NON-INFRINGEMENT. Please see the
19: * License for the specific language governing rights and limitations
20: * under the License.
21: *
22: * @APPLE_LICENSE_HEADER_END@
23: */
24:
25: .file "divmul.s"
26:
27: .ident "@(#)kern-fp:divmul.s 1.1"
28:
29: // ***********************************************************************
30: //
31: // m u l d i v . m o d
32: // ===================
33: //
34: // ===============================================================
35: // intel corporation proprietary information
36: // this software is supplied under the terms of a license
37: // agreement or non-disclosure agreement with intel corporation
38: // and may not be copied nor disclosed except in accordance with
39: // the terms of that agreement.
40: // ===============================================================
41: //
42: // function:
43: // preforms floating-point divide of unpacked
44: // non-zero, valid numbers.
45: // performs floating-point multiply of unpacked
46: // non-zero, valid numbers.
47: //
48: // public procedures:
49: // mulx mult
50: // divx divid
51: //
52: // ****************************************************************************
53: //
54: //...september 16, 1983...
55: //
56: // .file *a_mdm*
57: //$nolist
58: #include "fp_e80387.h"
59: //$list
60: .text //a_med segment er public
61: //
62: // extrn sticky_right_shift,one_left_normalize
63: // extrn clear_5w,test_5w,test_3w,set_5w
64: // extrn get_precision,set_i_error
65: // extrn left_shift_result_cl
66: //
67: .globl accel_divx
68: .globl divx
69: .globl divid
70: .globl log_divx
71: .globl mulx
72: .globl mult
73: //
74: ALIGN
75: quotient_length:
76: .byte 28,36,57,68 // incremented and changed to
77: log_quotient_length:
78: .byte 28,36,60,68 // a byte table on 12/02/82.
79: // for unknown reasons, divx
80: // doesnt work with 53-bit
81: // precision from log function
82: low_quotient_byte:
83: .byte offset_result+6,offset_result+4
84: .byte offset_result+2,offset_result+1
85: //$eject
86: // *********************************************************************
87: // accel_divx:
88: // **********
89: // function:
90: // fractional divide. result_frac <-- frac1/frac2.
91: //
92: // inputs:
93: // frac2 is assumed to be normalized and non-zero.
94: //
95: // outputs:
96: // the sticky bit is set for result_frac
97: // the remainder for ((frac1)/2)mod(frac2) is left in frac1
98: //
99: // data accessed:
100: // - dword_frac1 offset_operand2
101: // - dword_frac2 result_dword_frac
102: // - lsb_result offset_operand1
103: //
104: // data changed:
105: // - dword_frac1 result_dword_frac
106: // - lsb_result
107: //
108: // procedures called:
109: // sticky_right_shift get_precision
110: //
111: // *********************************************************************
112: ALIGN
113: divisor_dwords:
114: .byte 1,2,2,2
115: quotient_dwords:
116: .byte 1,2,2,3
117:
118: #define r_apprx_j extra_dword_reg(%ebp)
119: #define q_apprx_j extra_dword_reg+4(%ebp)
120: #define carry extra_dword_reg+8(%ebp)
121:
122: ALIGN
123: accel_divx: //proc
124: clr_quo:
125: mov $offset_result,%edi
126: call clear_6w // clears eax as well.
127: mov $offset_operand1,%edi // shift frac1 (here, the dividend)
128: movb $1,%cl // right by 1 bit
129: call sticky_right_shift
130: call get_precision
131: cmpb prec64,%dl
132: je prec_to_base
133: movb dword_frac1+3(%ebp),%cl
134: movb %cl,lsb_result
135: prec_to_base:
136: movzbl %dl,%ebx
137: movl $quotient_dwords,%edi //(edi) has offset quotient_dwords table
138: movzbl %cs:(%ebx,%edi),%ecx // q holds the number of dwords
139: movb %cl,q // to be produced for the quotient.
140: lea dword_frac1+frac64(%ebp),%esi
141: push %ss:4(%esi)
142:
143: testb $1, is16bit(%ebp)
144: jz prec_32
145: movzwl %sp, %edx
146: jmp prec_32c
147: ALIGN
148: prec_32:
149: mov %esp,%edx // edx holds the offset relative to ss of the
150: // interim partial remainders highest dword.
151: prec_32c:
152: dec %ecx
153: jnz push_dvdnd
154: inc %ecx
155: FALLSTHRU
156: push_dvdnd:
157: push %ss:(%esi)
158: sub $4,%esi
159: LOOP(push_dvdnd)
160: set_up_push_dvsr:
161: mov %ecx,%eax // ecx is 0 here.
162: movl $divisor_dwords,%edi
163: movb %cs:(%ebx,%edi),%cl // bit_ct holds the number of
164: movb %cl,bit_ct // dwords in the divisor.
165: FALLSTHRU
166: push_zeroes:
167: push %eax
168: LOOP(push_zeroes)
169: mov %edx,%ebx // ebx holds the offset relative to ss of the
170: // interim partial remainders highest word.
171: mov %eax,%edi // edi is index to words of quotient.
172: movb q,%cl
173: main_loop:
174: push %ecx
175: mov %ss:(%ebx),%edx
176: cmp dword_frac2+frac32(%ebp),%edx // is dvsr(1) = prem(j)?
177: jne do_divide
178: movl $0x0ffffffff,%eax // eax <-- apprx_q(j) = 2**32 - 1
179: movl %ss:-4(%ebx),%edx // edx <-- prem(j+1)
180: jmp get_r_apprx_j
181: ALIGN
182: do_divide:
183: movl %ss:-4(%ebx),%eax // (edx,eax) <-- prem(j) * 2**32 + prem(j+1)
184: divl dword_frac2+frac32(%ebp)
185: mov %eax,q_apprx_j
186: mov %edx,r_apprx_j
187: jmp test_q_apprx_j
188: ALIGN
189: dec_q_apprx_j:
190: mov q_apprx_j,%eax
191: dec %eax
192: mov r_apprx_j,%edx
193: get_r_apprx_j:
194: mov %eax,q_apprx_j
195: add dword_frac2+frac32(%ebp),%edx // r_apprx_j <-- edx + dvsr(1)
196: jc adjst_prem
197: mov %edx,r_apprx_j
198: test_q_apprx_j:
199: mov q_apprx_j,%eax
200: mull dword_frac2+frac64(%ebp) // (edx,eax) <-- dvsr(2)*q_apprx_j
201: cmp r_apprx_j,%edx
202: jb adjst_prem
203: ja dec_q_apprx_j
204: cmp %ss:-8(%ebx),%eax // is eax > prem(j+2)?
205: ja dec_q_apprx_j
206: adjst_prem:
207: xor %eax,%eax
208: mov %eax,carry
209: movb bit_ct,%al
210: mov %eax,%ecx
211: shl $2,%eax
212: neg %eax
213: mov %eax,%esi
214: FALLSTHRU
215: prem_loop:
216: mov sign2(%ebp,%esi),%eax
217: mull q_apprx_j
218: sub %eax,%ss:(%ebx,%esi)
219: jnc sbtrct_carry
220: inc %edx
221: sbtrct_carry:
222: mov carry,%eax
223: sub %eax,%ss:(%ebx,%esi)
224: jnc next_carry
225: inc %edx
226: next_carry:
227: mov %edx,carry
228: add $4,%esi
229: LOOP(prem_loop)
230: //
231: sub %edx,%ss:(%ebx,%esi) // here, esi = 0.
232: jnc next_j
233: decl q_apprx_j
234: movb bit_ct,%cl
235: mov %ecx,%eax
236: shl $2,%eax
237: neg %eax
238: mov %eax,%esi
239: clc
240: FALLSTHRU
241: fix_prem_loop:
242: mov sign2(%ebp,%esi),%eax
243: adc %eax,%ss:(%ebx,%esi)
244: add $4,%esi
245: LOOP(fix_prem_loop)
246: //
247: adc %ecx,%ss:(%ebx,%esi)
248: next_j:
249: mov q_apprx_j,%eax
250: mov %eax,result_dword_frac+frac32(%ebp,%edi)
251: pop %ecx
252: dec %ecx
253: jz get_sticky_bit
254: sub $4,%ebx
255: sub $4,%edi
256: jmp main_loop
257: ALIGN
258: //
259: get_sticky_bit:
260: mov %ecx,%eax
261: mov %eax,%esi
262: movb bit_ct,%cl
263: shlb $1,%cl
264: FALLSTHRU
265: sticky_loop:
266: // orw %ax,%ss:(%ebx+%esi-2) // bug in asm386 hardcoded
267: .byte 0x066,0x36,0x00b,0x44,0x33,0x0fe
268:
269: sub $2,%esi
270: LOOP(sticky_loop)
271: orb %al,%ah
272: orb %ah,lsb_result
273: mov %ecx,%esi
274: movb bit_ct,%cl
275: //
276: FALLSTHRU
277: stor_rmndr:
278: // movl %eax,%ss:(%ebx+%esi-4) // bug in asm386 hard coded
279: .byte 0x36,0x8b,0x44,0x33,0x0fc
280:
281: mov %eax,dword_frac1+frac32(%ebp,%esi)
282: sub $4,%esi
283: LOOP(stor_rmndr)
284: //
285: movb q,%cl
286: addb bit_ct,%cl
287: cmpb $2,%cl
288: jne get_stack_dwords
289: inc %ecx
290: get_stack_dwords:
291: shl $2,%ecx
292: add %ecx,%esp // restore stack
293: //
294: ret
295: //accel_divx endp
296: // *********************************************************************
297: // divx:
298: // ****
299: // function:
300: // fractional divide. result_frac <-- frac1/frac2.
301: //
302: // inputs:
303: // frac2 is assumed to be normalized and non-zero.
304: //
305: // outputs:
306: // the sticky bit is set for result_frac
307: // the remainder is left in frac1
308: //
309: // data accessed:
310: // - word_frac1 offset_operand2
311: // - word_frac2 result_word_frac
312: // - lsb_result
313: //
314: // data changed:
315: // - word_frac1 result_word_frac
316: // - lsb_result
317: //
318: // procedures called:
319: // sticky_right_shift get_precision
320: //
321: // *********************************************************************
322: ALIGN
323: divx: //proc
324: push $quotient_length // set normal entry ptr
325: jmp clear_quotient
326: ALIGN
327: log_divx:
328: push $log_quotient_length // set log entry ptr
329: clear_quotient:
330: mov $offset_result,%edi
331: call clear_6w
332: mov $offset_operand2,%edi //shift frac2 (divisor)
333: movb $1,%cl //right by 1 bit
334: call sticky_right_shift
335: call get_precision
336: movzbl %dl, %ebx
337: pop %edi //retrieve offset quotient length table
338: movb %cs:(%ebx,%edi),%al // bit_ct holds
339: movb %al,bit_ct // the quotient bits
340: movzbl %cs:low_quotient_byte(%ebx),%eax //(edi) points to low-
341: //order byte of quotient
342: mov %eax,%edi
343: movl dword_frac1(%ebp),%eax //load dividend (frac1)
344: movl dword_frac1+frac64(%ebp),%ebx // in frac1
345: movl dword_frac1+frac32(%ebp),%edx
346: .long 0x02d8ac0f //shrd eax, ebx, 2
347: .long 0x02d3ac0f //shrd ebx, edx, 2
348: shrl $2,%edx //shift dvdnd rt 2 bits
349: jmp subtract_divisor
350: ALIGN
351: frac_divide_loop:
352: shll $1,%eax //shift dvdnd (partial
353: rcll $1,%ebx //remainder) lft one bit
354: rcll $1,%edx
355: jc quotient_bit_0 //jump if cy from shift
356: orb $0x20,frac80(%ebp,%edi) // *shift in* a 1-bit
357: subtract_divisor:
358: subl dword_frac2(%ebp),%eax //subtract divisor from
359: sbbl dword_frac2+frac64(%ebp),%ebx //partial remainder
360: sbbl dword_frac2+frac32(%ebp),%edx
361: jmp shift_quotient_left
362: ALIGN
363: quotient_bit_0:
364: call add_divisor // add divisor to prem
365: shift_quotient_left:
366: call shift_result_left
367: decb bit_ct // decrement bit count
368: jnz frac_divide_loop // next quotient bit
369: andl %edx,%edx // branch if remainder
370: js adjust_remainder // is negative
371: orb $0x20,frac80(%ebp,%edi) // *shift in* last 1-bit
372: jmp store_remainder
373: ALIGN
374: adjust_remainder:
375: call add_divisor //add divisor to prem
376: store_remainder:
377: movl %eax,dword_frac1(%ebp) //store partial rmndr
378: movl %ebx,dword_frac1+frac64(%ebp) // in frac1
379: movl %edx,dword_frac1+frac32(%ebp)
380: orl %ebx,%eax //set sticky bits if
381: orl %edx,%eax //partial rmndr non-zero
382: .long 0x10c2a40f // shld edx, eax, 16
383: orw %dx, %ax
384: orb %al,%ah
385: movb %ah,lsb_result
386: ret
387: ALIGN
388: add_divisor:
389: addl dword_frac2(%ebp),%eax //store partial rmndr
390: adcl dword_frac2+frac64(%ebp),%ebx // in frac1
391: adcl dword_frac2+frac32(%ebp),%edx
392: ret
393: ALIGN
394: shift_result_left:
395: shll $1,result_dword_frac(%ebp)
396: rcll $1,result_dword_frac+frac64(%ebp)
397: rcll $1,result_dword_frac+frac32(%ebp)
398: ret
399: //divx endp
400: //$eject
401: // *************************************************************************
402: // divid:
403: // ******
404: // function:
405: // floating-point divide.
406: //
407: // inputs:
408: // assumes operands are unpacked, non-zero, and valid
409: //
410: // outputs:
411: // calculates unpacked result and returns with al
412: // set to true if underflow is possible, false if
413: // overflow is possible. the quotient is left in
414: // the result, and the remainder is left in frac1.
415: //
416: // data accessed:
417: // - expon1 expon2
418: // - msb_frac2 offset_result
419: // - result_sign result_expon
420: // - result_word_frac
421: //
422: // data changed:
423: // - result_sign result_expon
424: // - result_word_frac
425: //
426: // procedures called:
427: // divx one_left_normalize
428: // set_i_error get_precision
429: // left_shift_result_cl
430: //
431: // *************************************************************************
432: ALIGN
433: divid: //proc
434: movl dword_expon1,%eax // stack underflow possible
435: subl dword_expon2,%eax // flag (sign bit)
436: addl $exponent_bias,%eax //form biased exponent
437: movl %eax,dword_result_expon
438: movb sign1(%ebp),%al // sign = '+' if sign1 = sign2
439: xorb sign2(%ebp),%al
440: movb %al,result_sign(%ebp)
441: fractional_divide:
442: call accel_divx
443: // The following five lines are not to be used when accel_divx
444: // does the division instead of divx.
445: // call get_precision
446: // cmpb prec53,%dl /if double precision, shift
447: // jne norm_quotient / quotient left 3 bits
448: // movb $3,%cl
449: // call left_shift_result_cl
450: norm_quotient:
451: mov $offset_result,%edi //normalize by 1 left
452: jmp one_left_normalize // shift, if unnormalized
453: ALIGN
454: //divid endp
455: //$eject
456: // **********************************************************************
457: // mult:
458: // ****
459: // function:
460: // floating point multiply.
461: //
462: // inputs:
463: // assumes operands are unpacked, valid, non-zero.
464: //
465: // outputs:
466: // calculates unpacked result and returns with al set
467: // to true if underflow is possible, false if overflow
468: // is possible.
469: //
470: // data accessed:
471: // - sign1 expon1
472: // - sign2 expon2
473: // - offset_result result_sign
474: // - result_expon
475: //
476: // data changed:
477: // - result
478: //
479: // procedures called:
480: // mulx norm_quotient
481: //
482: // ****************************************************************************
483: ALIGN
484: mult: //proc
485: movl dword_expon1,%eax //form doubly-biased exponent
486: addl dword_expon2,%eax //if high bit set, underflow
487: subl $exponent_bias-1,%eax //form singly-biased exponent
488: movl %eax,dword_result_expon
489: movb sign1(%ebp),%al // result sign = sign1 xor sign2
490: xorb sign2(%ebp),%al
491: movb %al,result_sign(%ebp)
492: call mulx
493: movl $offset_result, %edi //normalize by 1 left
494: jmp one_left_normalize // shift if unnormalized
495: //mult endp
496: //$eject
497: // **********************************************************************
498: // mulx:
499: // *****
500: // function:
501: // fractional multiply. result_frac <-- frac1 * frac2.
502: //
503: // inputs:
504: // assumes the operands are unpacked, valid, non-zero.
505: //
506: // outputs:
507: // product in result_frac (sticky indicator left in low bit)
508: //
509: // data accessed:
510: // - word_frac1 offset_operand1
511: // - word_frac2 offset_operand2
512: // - extra_word_reg lsb_result
513: // - offset_result
514: //
515: // data changed:
516: // - extra_word_reg offset_result
517: //
518: // procedures called:
519: // clear_5w set_5w
520: // test_3w test_5w
521: //
522: // ************************************************************************
523: ALIGN
524: mulx: //proc
525: push %ds // save a?msr
526: movl $extra_dword_reg,%edi // ss:bp+di => extra_word_reg
527: push %ss
528: pop %ds
529: lea (%ebp,%edi),%ebx // ds:bx also => extra_word_reg
530: call clear_6w // clear extra_word_reg
531: movl $offset_result,%edi // clear result_frac
532: call set_6w
533: movl $frac32,%ecx // load s.p. offset
534: movl $offset_operand2,%edi //di => to multiplier
535: call test_4w //if low 3 words <> zero,
536: jnz examine_frac1 // branch if non single
537: add %ecx,%edi //frac2 is s. p.,
538: add %ecx,%ebx //so adjust pointers
539: examine_frac1:
540: movl $offset_operand1,%esi //si points to multiplicand
541: movl frac64(%ebp,%esi),%eax //if low 3 words = zero,
542: orl (%ebp,%esi),%eax //then single precision
543: jnz do_frac_multiply
544: add %ecx,%esi //frac1 is s. p.,
545: add %ecx,%ebx //so adjust pointers
546: do_frac_multiply:
547: push %edi // save frac2 offset
548: movl (%ebp,%edi),%edi // load multiplier
549: xorl %ecx,%ecx //clear cx
550: movl (%ebp,%esi),%eax // multiply first word
551: mull %edi
552: addl %eax,(%ebx) //add to partial product
553: adcl %edx,%ecx // cx initially 0
554: cmpl $offset_operand1, %esi
555: je mult_scnd
556: movl %ecx, 4(%ebx) //multiplicand is s. p. so
557: jmp end_of_mul_loop //go to multipliers next word
558: mult_scnd:
559: movl 4(%ebp,%esi),%eax // multiply second word
560: mull %edi
561: addl %ecx,%eax
562: adcl $0,%edx
563: xorl %ecx,%ecx
564: addl %eax,4(%ebx) //add to partial product
565: adcl %edx,%ecx
566: movl 8(%ebp,%esi), %eax
567: mull %edi
568: addl %ecx, %eax
569: adcl $0, %edx
570: addl %eax, 8(%ebx)
571: adcl $0, %edx
572: movl %edx, 12(%ebx)
573: end_of_mul_loop:
574: pop %edi // reload frac2 offset
575: addl $4, %ebx //adjust pointers for next iteration
576: addl $4, %edi
577: cmp $(offset_operand2+12),%edi
578: jne do_frac_multiply
579: movl $extra_dword_reg,%edi //set sticky bit if any extra
580: call test_6w // reg words are nonzero
581: jz frac_mult_done
582: orb $0x01,lsb_result
583: frac_mult_done:
584: pop %ds // restore a?msr
585: ret
586: //mulx endp
587: //
588: //a_med ends
589: //
590: // end
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.