![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to lb1sf68.asm CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [GCC 2.x] / gcc / config / m68k|
Return to lb1sf68.asm CVS log | Up to [GCC 2.x] / gcc / config / m68k |
1.1 ! root 1: /* libgcc1 routines for 68000 w/o floating-point hardware. */ ! 2: /* Copyright (C) 1994 Free Software Foundation, Inc. ! 3: ! 4: This file is free software; you can redistribute it and/or modify it ! 5: under the terms of the GNU General Public License as published by the ! 6: Free Software Foundation; either version 2, or (at your option) any ! 7: later version. ! 8: ! 9: In addition to the permissions in the GNU General Public License, the ! 10: Free Software Foundation gives you unlimited permission to link the ! 11: compiled version of this file with other programs, and to distribute ! 12: those programs without any restriction coming from the use of this ! 13: file. (The General Public License restrictions do apply in other ! 14: respects; for example, they cover modification of the file, and ! 15: distribution when not linked into another program.) ! 16: ! 17: This file is distributed in the hope that it will be useful, but ! 18: WITHOUT ANY WARRANTY; without even the implied warranty of ! 19: MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU ! 20: General Public License for more details. ! 21: ! 22: You should have received a copy of the GNU General Public License ! 23: along with this program; see the file COPYING. If not, write to ! 24: the Free Software Foundation, 59 Temple Place - Suite 330, ! 25: Boston, MA 02111-1307, USA. */ ! 26: ! 27: /* As a special exception, if you link this library with files ! 28: compiled with GCC to produce an executable, this does not cause ! 29: the resulting executable to be covered by the GNU General Public License. ! 30: This exception does not however invalidate any other reasons why ! 31: the executable file might be covered by the GNU General Public License. */ ! 32: ! 33: /* Use this one for any 680x0; assumes no floating point hardware. ! 34: The trailing " '" appearing on some lines is for ANSI preprocessors. Yuk. ! 35: Some of this code comes from MINIX, via the folks at ericsson. ! 36: D. V. Henkel-Wallace ([email protected]) Fete Bastille, 1992 ! 37: */ ! 38: ! 39: /* These are predefined by new versions of GNU cpp. */ ! 40: ! 41: #ifndef __USER_LABEL_PREFIX__ ! 42: #define __USER_LABEL_PREFIX__ _ ! 43: #endif ! 44: ! 45: #ifndef __REGISTER_PREFIX__ ! 46: #define __REGISTER_PREFIX__ ! 47: #endif ! 48: ! 49: #ifndef __IMMEDIATE_PREFIX__ ! 50: #define __IMMEDIATE_PREFIX__ # ! 51: #endif ! 52: ! 53: /* ANSI concatenation macros. */ ! 54: ! 55: #define CONCAT1(a, b) CONCAT2(a, b) ! 56: #define CONCAT2(a, b) a ## b ! 57: ! 58: /* Use the right prefix for global labels. */ ! 59: ! 60: #define SYM(x) CONCAT1 (__USER_LABEL_PREFIX__, x) ! 61: ! 62: /* Use the right prefix for registers. */ ! 63: ! 64: #define REG(x) CONCAT1 (__REGISTER_PREFIX__, x) ! 65: ! 66: /* Use the right prefix for immediate values. */ ! 67: ! 68: #define IMM(x) CONCAT1 (__IMMEDIATE_PREFIX__, x) ! 69: ! 70: #define d0 REG (d0) ! 71: #define d1 REG (d1) ! 72: #define d2 REG (d2) ! 73: #define d3 REG (d3) ! 74: #define d4 REG (d4) ! 75: #define d5 REG (d5) ! 76: #define d6 REG (d6) ! 77: #define d7 REG (d7) ! 78: #define a0 REG (a0) ! 79: #define a1 REG (a1) ! 80: #define a2 REG (a2) ! 81: #define a3 REG (a3) ! 82: #define a4 REG (a4) ! 83: #define a5 REG (a5) ! 84: #define a6 REG (a6) ! 85: #define fp REG (fp) ! 86: #define sp REG (sp) ! 87: ! 88: #ifdef L_floatex ! 89: ! 90: | This is an attempt at a decent floating point (single, double and ! 91: | extended double) code for the GNU C compiler. It should be easy to ! 92: | adapt to other compilers (but beware of the local labels!). ! 93: ! 94: | Starting date: 21 October, 1990 ! 95: ! 96: | It is convenient to introduce the notation (s,e,f) for a floating point ! 97: | number, where s=sign, e=exponent, f=fraction. We will call a floating ! 98: | point number fpn to abbreviate, independently of the precision. ! 99: | Let MAX_EXP be in each case the maximum exponent (255 for floats, 1023 ! 100: | for doubles and 16383 for long doubles). We then have the following ! 101: | different cases: ! 102: | 1. Normalized fpns have 0 < e < MAX_EXP. They correspond to ! 103: | (-1)^s x 1.f x 2^(e-bias-1). ! 104: | 2. Denormalized fpns have e=0. They correspond to numbers of the form ! 105: | (-1)^s x 0.f x 2^(-bias). ! 106: | 3. +/-INFINITY have e=MAX_EXP, f=0. ! 107: | 4. Quiet NaN (Not a Number) have all bits set. ! 108: | 5. Signaling NaN (Not a Number) have s=0, e=MAX_EXP, f=1. ! 109: ! 110: |============================================================================= ! 111: | exceptions ! 112: |============================================================================= ! 113: ! 114: | This is the floating point condition code register (_fpCCR): ! 115: | ! 116: | struct { ! 117: | short _exception_bits; ! 118: | short _trap_enable_bits; ! 119: | short _sticky_bits; ! 120: | short _rounding_mode; ! 121: | short _format; ! 122: | short _last_operation; ! 123: | union { ! 124: | float sf; ! 125: | double df; ! 126: | } _operand1; ! 127: | union { ! 128: | float sf; ! 129: | double df; ! 130: | } _operand2; ! 131: | } _fpCCR; ! 132: ! 133: .data ! 134: .even ! 135: ! 136: .globl SYM (_fpCCR) ! 137: ! 138: SYM (_fpCCR): ! 139: __exception_bits: ! 140: .word 0 ! 141: __trap_enable_bits: ! 142: .word 0 ! 143: __sticky_bits: ! 144: .word 0 ! 145: __rounding_mode: ! 146: .word ROUND_TO_NEAREST ! 147: __format: ! 148: .word NIL ! 149: __last_operation: ! 150: .word NOOP ! 151: __operand1: ! 152: .long 0 ! 153: .long 0 ! 154: __operand2: ! 155: .long 0 ! 156: .long 0 ! 157: ! 158: | Offsets: ! 159: EBITS = __exception_bits - SYM (_fpCCR) ! 160: TRAPE = __trap_enable_bits - SYM (_fpCCR) ! 161: STICK = __sticky_bits - SYM (_fpCCR) ! 162: ROUND = __rounding_mode - SYM (_fpCCR) ! 163: FORMT = __format - SYM (_fpCCR) ! 164: LASTO = __last_operation - SYM (_fpCCR) ! 165: OPER1 = __operand1 - SYM (_fpCCR) ! 166: OPER2 = __operand2 - SYM (_fpCCR) ! 167: ! 168: | The following exception types are supported: ! 169: INEXACT_RESULT = 0x0001 ! 170: UNDERFLOW = 0x0002 ! 171: OVERFLOW = 0x0004 ! 172: DIVIDE_BY_ZERO = 0x0008 ! 173: INVALID_OPERATION = 0x0010 ! 174: ! 175: | The allowed rounding modes are: ! 176: UNKNOWN = -1 ! 177: ROUND_TO_NEAREST = 0 | round result to nearest representable value ! 178: ROUND_TO_ZERO = 1 | round result towards zero ! 179: ROUND_TO_PLUS = 2 | round result towards plus infinity ! 180: ROUND_TO_MINUS = 3 | round result towards minus infinity ! 181: ! 182: | The allowed values of format are: ! 183: NIL = 0 ! 184: SINGLE_FLOAT = 1 ! 185: DOUBLE_FLOAT = 2 ! 186: LONG_FLOAT = 3 ! 187: ! 188: | The allowed values for the last operation are: ! 189: NOOP = 0 ! 190: ADD = 1 ! 191: MULTIPLY = 2 ! 192: DIVIDE = 3 ! 193: NEGATE = 4 ! 194: COMPARE = 5 ! 195: EXTENDSFDF = 6 ! 196: TRUNCDFSF = 7 ! 197: ! 198: |============================================================================= ! 199: | __clear_sticky_bits ! 200: |============================================================================= ! 201: ! 202: | The sticky bits are normally not cleared (thus the name), whereas the ! 203: | exception type and exception value reflect the last computation. ! 204: | This routine is provided to clear them (you can also write to _fpCCR, ! 205: | since it is globally visible). ! 206: ! 207: .globl SYM (__clear_sticky_bit) ! 208: ! 209: .text ! 210: .even ! 211: ! 212: | void __clear_sticky_bits(void); ! 213: SYM (__clear_sticky_bit): ! 214: lea SYM (_fpCCR),a0 ! 215: movew IMM (0),a0@(STICK) ! 216: rts ! 217: ! 218: |============================================================================= ! 219: | $_exception_handler ! 220: |============================================================================= ! 221: ! 222: .globl $_exception_handler ! 223: ! 224: .text ! 225: .even ! 226: ! 227: | This is the common exit point if an exception occurs. ! 228: | NOTE: it is NOT callable from C! ! 229: | It expects the exception type in d7, the format (SINGLE_FLOAT, ! 230: | DOUBLE_FLOAT or LONG_FLOAT) in d6, and the last operation code in d5. ! 231: | It sets the corresponding exception and sticky bits, and the format. ! 232: | Depending on the format if fills the corresponding slots for the ! 233: | operands which produced the exception (all this information is provided ! 234: | so if you write your own exception handlers you have enough information ! 235: | to deal with the problem). ! 236: | Then checks to see if the corresponding exception is trap-enabled, ! 237: | in which case it pushes the address of _fpCCR and traps through ! 238: | trap FPTRAP (15 for the moment). ! 239: ! 240: FPTRAP = 15 ! 241: ! 242: $_exception_handler: ! 243: lea SYM (_fpCCR),a0 ! 244: movew d7,a0@(EBITS) | set __exception_bits ! 245: orw d7,a0@(STICK) | and __sticky_bits ! 246: movew d6,a0@(FORMT) | and __format ! 247: movew d5,a0@(LASTO) | and __last_operation ! 248: ! 249: | Now put the operands in place: ! 250: cmpw IMM (SINGLE_FLOAT),d6 ! 251: beq 1f ! 252: movel a6@(8),a0@(OPER1) ! 253: movel a6@(12),a0@(OPER1+4) ! 254: movel a6@(16),a0@(OPER2) ! 255: movel a6@(20),a0@(OPER2+4) ! 256: bra 2f ! 257: 1: movel a6@(8),a0@(OPER1) ! 258: movel a6@(12),a0@(OPER2) ! 259: 2: ! 260: | And check whether the exception is trap-enabled: ! 261: andw a0@(TRAPE),d7 | is exception trap-enabled? ! 262: beq 1f | no, exit ! 263: pea SYM (_fpCCR) | yes, push address of _fpCCR ! 264: trap IMM (FPTRAP) | and trap ! 265: 1: moveml sp@+,d2-d7 | restore data registers ! 266: unlk a6 | and return ! 267: rts ! 268: #endif /* L_floatex */ ! 269: ! 270: #ifdef L_mulsi3 ! 271: .text ! 272: .proc ! 273: .globl SYM (__mulsi3) ! 274: SYM (__mulsi3): ! 275: movew sp@(4), d0 /* x0 -> d0 */ ! 276: muluw sp@(10), d0 /* x0*y1 */ ! 277: movew sp@(6), d1 /* x1 -> d1 */ ! 278: muluw sp@(8), d1 /* x1*y0 */ ! 279: addw d1, d0 ! 280: swap d0 ! 281: clrw d0 ! 282: movew sp@(6), d1 /* x1 -> d1 */ ! 283: muluw sp@(10), d1 /* x1*y1 */ ! 284: addl d1, d0 ! 285: ! 286: rts ! 287: #endif /* L_mulsi3 */ ! 288: ! 289: #ifdef L_udivsi3 ! 290: .text ! 291: .proc ! 292: .globl SYM (__udivsi3) ! 293: SYM (__udivsi3): ! 294: movel d2, sp@- ! 295: movel sp@(12), d1 /* d1 = divisor */ ! 296: movel sp@(8), d0 /* d0 = dividend */ ! 297: ! 298: cmpl IMM (0x10000), d1 /* divisor >= 2 ^ 16 ? */ ! 299: jcc L3 /* then try next algorithm */ ! 300: movel d0, d2 ! 301: clrw d2 ! 302: swap d2 ! 303: divu d1, d2 /* high quotient in lower word */ ! 304: movew d2, d0 /* save high quotient */ ! 305: swap d0 ! 306: movew sp@(10), d2 /* get low dividend + high rest */ ! 307: divu d1, d2 /* low quotient */ ! 308: movew d2, d0 ! 309: jra L6 ! 310: ! 311: L3: movel d1, d2 /* use d2 as divisor backup */ ! 312: L4: lsrl IMM (1), d1 /* shift divisor */ ! 313: lsrl IMM (1), d0 /* shift dividend */ ! 314: cmpl IMM (0x10000), d1 /* still divisor >= 2 ^ 16 ? */ ! 315: jcc L4 ! 316: divu d1, d0 /* now we have 16 bit divisor */ ! 317: andl IMM (0xffff), d0 /* mask out divisor, ignore remainder */ ! 318: ! 319: /* Multiply the 16 bit tentative quotient with the 32 bit divisor. Because of ! 320: the operand ranges, this might give a 33 bit product. If this product is ! 321: greater than the dividend, the tentative quotient was too large. */ ! 322: movel d2, d1 ! 323: mulu d0, d1 /* low part, 32 bits */ ! 324: swap d2 ! 325: mulu d0, d2 /* high part, at most 17 bits */ ! 326: swap d2 /* align high part with low part */ ! 327: btst IMM (0), d2 /* high part 17 bits? */ ! 328: jne L5 /* if 17 bits, quotient was too large */ ! 329: addl d2, d1 /* add parts */ ! 330: jcs L5 /* if sum is 33 bits, quotient was too large */ ! 331: cmpl sp@(8), d1 /* compare the sum with the dividend */ ! 332: jls L6 /* if sum > dividend, quotient was too large */ ! 333: L5: subql IMM (1), d0 /* adjust quotient */ ! 334: ! 335: L6: movel sp@+, d2 ! 336: rts ! 337: #endif /* L_udivsi3 */ ! 338: ! 339: #ifdef L_divsi3 ! 340: .text ! 341: .proc ! 342: .globl SYM (__divsi3) ! 343: SYM (__divsi3): ! 344: movel d2, sp@- ! 345: ! 346: moveb IMM (1), d2 /* sign of result stored in d2 (=1 or =-1) */ ! 347: movel sp@(12), d1 /* d1 = divisor */ ! 348: jpl L1 ! 349: negl d1 ! 350: negb d2 /* change sign because divisor <0 */ ! 351: L1: movel sp@(8), d0 /* d0 = dividend */ ! 352: jpl L2 ! 353: negl d0 ! 354: negb d2 ! 355: ! 356: L2: movel d1, sp@- ! 357: movel d0, sp@- ! 358: jbsr SYM (__udivsi3) /* divide abs(dividend) by abs(divisor) */ ! 359: addql IMM (8), sp ! 360: ! 361: tstb d2 ! 362: jpl L3 ! 363: negl d0 ! 364: ! 365: L3: movel sp@+, d2 ! 366: rts ! 367: #endif /* L_divsi3 */ ! 368: ! 369: #ifdef L_umodsi3 ! 370: .text ! 371: .proc ! 372: .globl SYM (__umodsi3) ! 373: SYM (__umodsi3): ! 374: movel sp@(8), d1 /* d1 = divisor */ ! 375: movel sp@(4), d0 /* d0 = dividend */ ! 376: movel d1, sp@- ! 377: movel d0, sp@- ! 378: jbsr SYM (__udivsi3) ! 379: addql IMM (8), sp ! 380: movel sp@(8), d1 /* d1 = divisor */ ! 381: movel d1, sp@- ! 382: movel d0, sp@- ! 383: jbsr SYM (__mulsi3) /* d0 = (a/b)*b */ ! 384: addql IMM (8), sp ! 385: movel sp@(4), d1 /* d1 = dividend */ ! 386: subl d0, d1 /* d1 = a - (a/b)*b */ ! 387: movel d1, d0 ! 388: rts ! 389: #endif /* L_umodsi3 */ ! 390: ! 391: #ifdef L_modsi3 ! 392: .text ! 393: .proc ! 394: .globl SYM (__modsi3) ! 395: SYM (__modsi3): ! 396: movel sp@(8), d1 /* d1 = divisor */ ! 397: movel sp@(4), d0 /* d0 = dividend */ ! 398: movel d1, sp@- ! 399: movel d0, sp@- ! 400: jbsr SYM (__divsi3) ! 401: addql IMM (8), sp ! 402: movel sp@(8), d1 /* d1 = divisor */ ! 403: movel d1, sp@- ! 404: movel d0, sp@- ! 405: jbsr SYM (__mulsi3) /* d0 = (a/b)*b */ ! 406: addql IMM (8), sp ! 407: movel sp@(4), d1 /* d1 = dividend */ ! 408: subl d0, d1 /* d1 = a - (a/b)*b */ ! 409: movel d1, d0 ! 410: rts ! 411: #endif /* L_modsi3 */ ! 412: ! 413: ! 414: #ifdef L_double ! 415: ! 416: .globl SYM (_fpCCR) ! 417: .globl $_exception_handler ! 418: ! 419: QUIET_NaN = 0xffffffff ! 420: ! 421: D_MAX_EXP = 0x07ff ! 422: D_BIAS = 1022 ! 423: DBL_MAX_EXP = D_MAX_EXP - D_BIAS ! 424: DBL_MIN_EXP = 1 - D_BIAS ! 425: DBL_MANT_DIG = 53 ! 426: ! 427: INEXACT_RESULT = 0x0001 ! 428: UNDERFLOW = 0x0002 ! 429: OVERFLOW = 0x0004 ! 430: DIVIDE_BY_ZERO = 0x0008 ! 431: INVALID_OPERATION = 0x0010 ! 432: ! 433: DOUBLE_FLOAT = 2 ! 434: ! 435: NOOP = 0 ! 436: ADD = 1 ! 437: MULTIPLY = 2 ! 438: DIVIDE = 3 ! 439: NEGATE = 4 ! 440: COMPARE = 5 ! 441: EXTENDSFDF = 6 ! 442: TRUNCDFSF = 7 ! 443: ! 444: UNKNOWN = -1 ! 445: ROUND_TO_NEAREST = 0 | round result to nearest representable value ! 446: ROUND_TO_ZERO = 1 | round result towards zero ! 447: ROUND_TO_PLUS = 2 | round result towards plus infinity ! 448: ROUND_TO_MINUS = 3 | round result towards minus infinity ! 449: ! 450: | Entry points: ! 451: ! 452: .globl SYM (__adddf3) ! 453: .globl SYM (__subdf3) ! 454: .globl SYM (__muldf3) ! 455: .globl SYM (__divdf3) ! 456: .globl SYM (__negdf2) ! 457: .globl SYM (__cmpdf2) ! 458: ! 459: .text ! 460: .even ! 461: ! 462: | These are common routines to return and signal exceptions. ! 463: ! 464: Ld$den: ! 465: | Return and signal a denormalized number ! 466: orl d7,d0 ! 467: movew IMM (UNDERFLOW),d7 ! 468: orw IMM (INEXACT_RESULT),d7 ! 469: movew IMM (DOUBLE_FLOAT),d6 ! 470: jmp $_exception_handler ! 471: ! 472: Ld$infty: ! 473: Ld$overflow: ! 474: | Return a properly signed INFINITY and set the exception flags ! 475: movel IMM (0x7ff00000),d0 ! 476: movel IMM (0),d1 ! 477: orl d7,d0 ! 478: movew IMM (OVERFLOW),d7 ! 479: orw IMM (INEXACT_RESULT),d7 ! 480: movew IMM (DOUBLE_FLOAT),d6 ! 481: jmp $_exception_handler ! 482: ! 483: Ld$underflow: ! 484: | Return 0 and set the exception flags ! 485: movel IMM (0),d0 ! 486: movel d0,d1 ! 487: movew IMM (UNDERFLOW),d7 ! 488: orw IMM (INEXACT_RESULT),d7 ! 489: movew IMM (DOUBLE_FLOAT),d6 ! 490: jmp $_exception_handler ! 491: ! 492: Ld$inop: ! 493: | Return a quiet NaN and set the exception flags ! 494: movel IMM (QUIET_NaN),d0 ! 495: movel d0,d1 ! 496: movew IMM (INVALID_OPERATION),d7 ! 497: orw IMM (INEXACT_RESULT),d7 ! 498: movew IMM (DOUBLE_FLOAT),d6 ! 499: jmp $_exception_handler ! 500: ! 501: Ld$div$0: ! 502: | Return a properly signed INFINITY and set the exception flags ! 503: movel IMM (0x7ff00000),d0 ! 504: movel IMM (0),d1 ! 505: orl d7,d0 ! 506: movew IMM (DIVIDE_BY_ZERO),d7 ! 507: orw IMM (INEXACT_RESULT),d7 ! 508: movew IMM (DOUBLE_FLOAT),d6 ! 509: jmp $_exception_handler ! 510: ! 511: |============================================================================= ! 512: |============================================================================= ! 513: | double precision routines ! 514: |============================================================================= ! 515: |============================================================================= ! 516: ! 517: | A double precision floating point number (double) has the format: ! 518: | ! 519: | struct _double { ! 520: | unsigned int sign : 1; /* sign bit */ ! 521: | unsigned int exponent : 11; /* exponent, shifted by 126 */ ! 522: | unsigned int fraction : 52; /* fraction */ ! 523: | } double; ! 524: | ! 525: | Thus sizeof(double) = 8 (64 bits). ! 526: | ! 527: | All the routines are callable from C programs, and return the result ! 528: | in the register pair d0-d1. They also preserve all registers except ! 529: | d0-d1 and a0-a1. ! 530: ! 531: |============================================================================= ! 532: | __subdf3 ! 533: |============================================================================= ! 534: ! 535: | double __subdf3(double, double); ! 536: SYM (__subdf3): ! 537: bchg IMM (31),sp@(12) | change sign of second operand ! 538: | and fall through, so we always add ! 539: |============================================================================= ! 540: | __adddf3 ! 541: |============================================================================= ! 542: ! 543: | double __adddf3(double, double); ! 544: SYM (__adddf3): ! 545: link a6,IMM (0) | everything will be done in registers ! 546: moveml d2-d7,sp@- | save all data registers and a2 (but d0-d1) ! 547: movel a6@(8),d0 | get first operand ! 548: movel a6@(12),d1 | ! 549: movel a6@(16),d2 | get second operand ! 550: movel a6@(20),d3 | ! 551: ! 552: movel d0,d7 | get d0's sign bit in d7 ' ! 553: addl d1,d1 | check and clear sign bit of a, and gain one ! 554: addxl d0,d0 | bit of extra precision ! 555: beq Ladddf$b | if zero return second operand ! 556: ! 557: movel d2,d6 | save sign in d6 ! 558: addl d3,d3 | get rid of sign bit and gain one bit of ! 559: addxl d2,d2 | extra precision ! 560: beq Ladddf$a | if zero return first operand ! 561: ! 562: andl IMM (0x80000000),d7 | isolate a's sign bit ' ! 563: swap d6 | and also b's sign bit ' ! 564: andw IMM (0x8000),d6 | ! 565: orw d6,d7 | and combine them into d7, so that a's sign ' ! 566: | bit is in the high word and b's is in the ' ! 567: | low word, so d6 is free to be used ! 568: movel d7,a0 | now save d7 into a0, so d7 is free to ! 569: | be used also ! 570: ! 571: | Get the exponents and check for denormalized and/or infinity. ! 572: ! 573: movel IMM (0x001fffff),d6 | mask for the fraction ! 574: movel IMM (0x00200000),d7 | mask to put hidden bit back ! 575: ! 576: movel d0,d4 | ! 577: andl d6,d0 | get fraction in d0 ! 578: notl d6 | make d6 into mask for the exponent ! 579: andl d6,d4 | get exponent in d4 ! 580: beq Ladddf$a$den | branch if a is denormalized ! 581: cmpl d6,d4 | check for INFINITY or NaN ! 582: beq Ladddf$nf | ! 583: orl d7,d0 | and put hidden bit back ! 584: Ladddf$1: ! 585: swap d4 | shift right exponent so that it starts ! 586: lsrw IMM (5),d4 | in bit 0 and not bit 20 ! 587: | Now we have a's exponent in d4 and fraction in d0-d1 ' ! 588: movel d2,d5 | save b to get exponent ! 589: andl d6,d5 | get exponent in d5 ! 590: beq Ladddf$b$den | branch if b is denormalized ! 591: cmpl d6,d5 | check for INFINITY or NaN ! 592: beq Ladddf$nf ! 593: notl d6 | make d6 into mask for the fraction again ! 594: andl d6,d2 | and get fraction in d2 ! 595: orl d7,d2 | and put hidden bit back ! 596: Ladddf$2: ! 597: swap d5 | shift right exponent so that it starts ! 598: lsrw IMM (5),d5 | in bit 0 and not bit 20 ! 599: ! 600: | Now we have b's exponent in d5 and fraction in d2-d3. ' ! 601: ! 602: | The situation now is as follows: the signs are combined in a0, the ! 603: | numbers are in d0-d1 (a) and d2-d3 (b), and the exponents in d4 (a) ! 604: | and d5 (b). To do the rounding correctly we need to keep all the ! 605: | bits until the end, so we need to use d0-d1-d2-d3 for the first number ! 606: | and d4-d5-d6-d7 for the second. To do this we store (temporarily) the ! 607: | exponents in a2-a3. ! 608: ! 609: moveml a2-a3,sp@- | save the address registers ! 610: ! 611: movel d4,a2 | save the exponents ! 612: movel d5,a3 | ! 613: ! 614: movel IMM (0),d7 | and move the numbers around ! 615: movel d7,d6 | ! 616: movel d3,d5 | ! 617: movel d2,d4 | ! 618: movel d7,d3 | ! 619: movel d7,d2 | ! 620: ! 621: | Here we shift the numbers until the exponents are the same, and put ! 622: | the largest exponent in a2. ! 623: exg d4,a2 | get exponents back ! 624: exg d5,a3 | ! 625: cmpw d4,d5 | compare the exponents ! 626: beq Ladddf$3 | if equal don't shift ' ! 627: bhi 9f | branch if second exponent is higher ! 628: ! 629: | Here we have a's exponent larger than b's, so we have to shift b. We do ! 630: | this by using as counter d2: ! 631: 1: movew d4,d2 | move largest exponent to d2 ! 632: subw d5,d2 | and subtract second exponent ! 633: exg d4,a2 | get back the longs we saved ! 634: exg d5,a3 | ! 635: | if difference is too large we don't shift (actually, we can just exit) ' ! 636: cmpw IMM (DBL_MANT_DIG+2),d2 ! 637: bge Ladddf$b$small ! 638: cmpw IMM (32),d2 | if difference >= 32, shift by longs ! 639: bge 5f ! 640: 2: cmpw IMM (16),d2 | if difference >= 16, shift by words ! 641: bge 6f ! 642: bra 3f | enter dbra loop ! 643: ! 644: 4: lsrl IMM (1),d4 ! 645: roxrl IMM (1),d5 ! 646: roxrl IMM (1),d6 ! 647: roxrl IMM (1),d7 ! 648: 3: dbra d2,4b ! 649: movel IMM (0),d2 ! 650: movel d2,d3 ! 651: bra Ladddf$4 ! 652: 5: ! 653: movel d6,d7 ! 654: movel d5,d6 ! 655: movel d4,d5 ! 656: movel IMM (0),d4 ! 657: subw IMM (32),d2 ! 658: bra 2b ! 659: 6: ! 660: movew d6,d7 ! 661: swap d7 ! 662: movew d5,d6 ! 663: swap d6 ! 664: movew d4,d5 ! 665: swap d5 ! 666: movew IMM (0),d4 ! 667: swap d4 ! 668: subw IMM (16),d2 ! 669: bra 3b ! 670: ! 671: 9: exg d4,d5 ! 672: movew d4,d6 ! 673: subw d5,d6 | keep d5 (largest exponent) in d4 ! 674: exg d4,a2 ! 675: exg d5,a3 ! 676: | if difference is too large we don't shift (actually, we can just exit) ' ! 677: cmpw IMM (DBL_MANT_DIG+2),d6 ! 678: bge Ladddf$a$small ! 679: cmpw IMM (32),d6 | if difference >= 32, shift by longs ! 680: bge 5f ! 681: 2: cmpw IMM (16),d6 | if difference >= 16, shift by words ! 682: bge 6f ! 683: bra 3f | enter dbra loop ! 684: ! 685: 4: lsrl IMM (1),d0 ! 686: roxrl IMM (1),d1 ! 687: roxrl IMM (1),d2 ! 688: roxrl IMM (1),d3 ! 689: 3: dbra d6,4b ! 690: movel IMM (0),d7 ! 691: movel d7,d6 ! 692: bra Ladddf$4 ! 693: 5: ! 694: movel d2,d3 ! 695: movel d1,d2 ! 696: movel d0,d1 ! 697: movel IMM (0),d0 ! 698: subw IMM (32),d6 ! 699: bra 2b ! 700: 6: ! 701: movew d2,d3 ! 702: swap d3 ! 703: movew d1,d2 ! 704: swap d2 ! 705: movew d0,d1 ! 706: swap d1 ! 707: movew IMM (0),d0 ! 708: swap d0 ! 709: subw IMM (16),d6 ! 710: bra 3b ! 711: Ladddf$3: ! 712: exg d4,a2 ! 713: exg d5,a3 ! 714: Ladddf$4: ! 715: | Now we have the numbers in d0--d3 and d4--d7, the exponent in a2, and ! 716: | the signs in a4. ! 717: ! 718: | Here we have to decide whether to add or subtract the numbers: ! 719: exg d7,a0 | get the signs ! 720: exg d6,a3 | a3 is free to be used ! 721: movel d7,d6 | ! 722: movew IMM (0),d7 | get a's sign in d7 ' ! 723: swap d6 | ! 724: movew IMM (0),d6 | and b's sign in d6 ' ! 725: eorl d7,d6 | compare the signs ! 726: bmi Lsubdf$0 | if the signs are different we have ! 727: | to subtract ! 728: exg d7,a0 | else we add the numbers ! 729: exg d6,a3 | ! 730: addl d7,d3 | ! 731: addxl d6,d2 | ! 732: addxl d5,d1 | ! 733: addxl d4,d0 | ! 734: ! 735: movel a2,d4 | return exponent to d4 ! 736: movel a0,d7 | ! 737: andl IMM (0x80000000),d7 | d7 now has the sign ! 738: ! 739: moveml sp@+,a2-a3 ! 740: ! 741: | Before rounding normalize so bit #DBL_MANT_DIG is set (we will consider ! 742: | the case of denormalized numbers in the rounding routine itself). ! 743: | As in the addition (not in the subtraction!) we could have set ! 744: | one more bit we check this: ! 745: btst IMM (DBL_MANT_DIG+1),d0 ! 746: beq 1f ! 747: lsrl IMM (1),d0 ! 748: roxrl IMM (1),d1 ! 749: roxrl IMM (1),d2 ! 750: roxrl IMM (1),d3 ! 751: addw IMM (1),d4 ! 752: 1: ! 753: lea Ladddf$5,a0 | to return from rounding routine ! 754: lea SYM (_fpCCR),a1 | check the rounding mode ! 755: movew a1@(6),d6 | rounding mode in d6 ! 756: beq Lround$to$nearest ! 757: cmpw IMM (ROUND_TO_PLUS),d6 ! 758: bhi Lround$to$minus ! 759: blt Lround$to$zero ! 760: bra Lround$to$plus ! 761: Ladddf$5: ! 762: | Put back the exponent and check for overflow ! 763: cmpw IMM (0x7ff),d4 | is the exponent big? ! 764: bge 1f ! 765: bclr IMM (DBL_MANT_DIG-1),d0 ! 766: lslw IMM (4),d4 | put exponent back into position ! 767: swap d0 | ! 768: orw d4,d0 | ! 769: swap d0 | ! 770: bra Ladddf$ret ! 771: 1: ! 772: movew IMM (ADD),d5 ! 773: bra Ld$overflow ! 774: ! 775: Lsubdf$0: ! 776: | Here we do the subtraction. ! 777: exg d7,a0 | put sign back in a0 ! 778: exg d6,a3 | ! 779: subl d7,d3 | ! 780: subxl d6,d2 | ! 781: subxl d5,d1 | ! 782: subxl d4,d0 | ! 783: beq Ladddf$ret$1 | if zero just exit ! 784: bpl 1f | if positive skip the following ! 785: exg d7,a0 | ! 786: bchg IMM (31),d7 | change sign bit in d7 ! 787: exg d7,a0 | ! 788: negl d3 | ! 789: negxl d2 | ! 790: negxl d1 | and negate result ! 791: negxl d0 | ! 792: 1: ! 793: movel a2,d4 | return exponent to d4 ! 794: movel a0,d7 ! 795: andl IMM (0x80000000),d7 | isolate sign bit ! 796: moveml sp@+,a2-a3 | ! 797: ! 798: | Before rounding normalize so bit #DBL_MANT_DIG is set (we will consider ! 799: | the case of denormalized numbers in the rounding routine itself). ! 800: | As in the addition (not in the subtraction!) we could have set ! 801: | one more bit we check this: ! 802: btst IMM (DBL_MANT_DIG+1),d0 ! 803: beq 1f ! 804: lsrl IMM (1),d0 ! 805: roxrl IMM (1),d1 ! 806: roxrl IMM (1),d2 ! 807: roxrl IMM (1),d3 ! 808: addw IMM (1),d4 ! 809: 1: ! 810: lea Lsubdf$1,a0 | to return from rounding routine ! 811: lea SYM (_fpCCR),a1 | check the rounding mode ! 812: movew a1@(6),d6 | rounding mode in d6 ! 813: beq Lround$to$nearest ! 814: cmpw IMM (ROUND_TO_PLUS),d6 ! 815: bhi Lround$to$minus ! 816: blt Lround$to$zero ! 817: bra Lround$to$plus ! 818: Lsubdf$1: ! 819: | Put back the exponent and sign (we don't have overflow). ' ! 820: bclr IMM (DBL_MANT_DIG-1),d0 ! 821: lslw IMM (4),d4 | put exponent back into position ! 822: swap d0 | ! 823: orw d4,d0 | ! 824: swap d0 | ! 825: bra Ladddf$ret ! 826: ! 827: | If one of the numbers was too small (difference of exponents >= ! 828: | DBL_MANT_DIG+1) we return the other (and now we don't have to ' ! 829: | check for finiteness or zero). ! 830: Ladddf$a$small: ! 831: moveml sp@+,a2-a3 ! 832: movel a6@(16),d0 ! 833: movel a6@(20),d1 ! 834: lea SYM (_fpCCR),a0 ! 835: movew IMM (0),a0@ ! 836: moveml sp@+,d2-d7 | restore data registers ! 837: unlk a6 | and return ! 838: rts ! 839: ! 840: Ladddf$b$small: ! 841: moveml sp@+,a2-a3 ! 842: movel a6@(8),d0 ! 843: movel a6@(12),d1 ! 844: lea SYM (_fpCCR),a0 ! 845: movew IMM (0),a0@ ! 846: moveml sp@+,d2-d7 | restore data registers ! 847: unlk a6 | and return ! 848: rts ! 849: ! 850: Ladddf$a$den: ! 851: movel d7,d4 | d7 contains 0x00200000 ! 852: bra Ladddf$1 ! 853: ! 854: Ladddf$b$den: ! 855: movel d7,d5 | d7 contains 0x00200000 ! 856: notl d6 ! 857: bra Ladddf$2 ! 858: ! 859: Ladddf$b: ! 860: | Return b (if a is zero) ! 861: movel d2,d0 ! 862: movel d3,d1 ! 863: bra 1f ! 864: Ladddf$a: ! 865: movel a6@(8),d0 ! 866: movel a6@(12),d1 ! 867: 1: ! 868: movew IMM (ADD),d5 ! 869: | Check for NaN and +/-INFINITY. ! 870: movel d0,d7 | ! 871: andl IMM (0x80000000),d7 | ! 872: bclr IMM (31),d0 | ! 873: cmpl IMM (0x7ff00000),d0 | ! 874: bge 2f | ! 875: movel d0,d0 | check for zero, since we don't ' ! 876: bne Ladddf$ret | want to return -0 by mistake ! 877: bclr IMM (31),d7 | ! 878: bra Ladddf$ret | ! 879: 2: ! 880: andl IMM (0x000fffff),d0 | check for NaN (nonzero fraction) ! 881: orl d1,d0 | ! 882: bne Ld$inop | ! 883: bra Ld$infty | ! 884: ! 885: Ladddf$ret$1: ! 886: moveml sp@+,a2-a3 | restore regs and exit ! 887: ! 888: Ladddf$ret: ! 889: | Normal exit. ! 890: lea SYM (_fpCCR),a0 ! 891: movew IMM (0),a0@ ! 892: orl d7,d0 | put sign bit back ! 893: moveml sp@+,d2-d7 ! 894: unlk a6 ! 895: rts ! 896: ! 897: Ladddf$ret$den: ! 898: | Return a denormalized number. ! 899: lsrl IMM (1),d0 | shift right once more ! 900: roxrl IMM (1),d1 | ! 901: bra Ladddf$ret ! 902: ! 903: Ladddf$nf: ! 904: movew IMM (ADD),d5 ! 905: | This could be faster but it is not worth the effort, since it is not ! 906: | executed very often. We sacrifice speed for clarity here. ! 907: movel a6@(8),d0 | get the numbers back (remember that we ! 908: movel a6@(12),d1 | did some processing already) ! 909: movel a6@(16),d2 | ! 910: movel a6@(20),d3 | ! 911: movel IMM (0x7ff00000),d4 | useful constant (INFINITY) ! 912: movel d0,d7 | save sign bits ! 913: movel d2,d6 | ! 914: bclr IMM (31),d0 | clear sign bits ! 915: bclr IMM (31),d2 | ! 916: | We know that one of them is either NaN of +/-INFINITY ! 917: | Check for NaN (if either one is NaN return NaN) ! 918: cmpl d4,d0 | check first a (d0) ! 919: bhi Ld$inop | if d0 > 0x7ff00000 or equal and ! 920: bne 2f ! 921: tstl d1 | d1 > 0, a is NaN ! 922: bne Ld$inop | ! 923: 2: cmpl d4,d2 | check now b (d1) ! 924: bhi Ld$inop | ! 925: bne 3f ! 926: tstl d3 | ! 927: bne Ld$inop | ! 928: 3: ! 929: | Now comes the check for +/-INFINITY. We know that both are (maybe not ! 930: | finite) numbers, but we have to check if both are infinite whether we ! 931: | are adding or subtracting them. ! 932: eorl d7,d6 | to check sign bits ! 933: bmi 1f ! 934: andl IMM (0x80000000),d7 | get (common) sign bit ! 935: bra Ld$infty ! 936: 1: ! 937: | We know one (or both) are infinite, so we test for equality between the ! 938: | two numbers (if they are equal they have to be infinite both, so we ! 939: | return NaN). ! 940: cmpl d2,d0 | are both infinite? ! 941: bne 1f | if d0 <> d2 they are not equal ! 942: cmpl d3,d1 | if d0 == d2 test d3 and d1 ! 943: beq Ld$inop | if equal return NaN ! 944: 1: ! 945: andl IMM (0x80000000),d7 | get a's sign bit ' ! 946: cmpl d4,d0 | test now for infinity ! 947: beq Ld$infty | if a is INFINITY return with this sign ! 948: bchg IMM (31),d7 | else we know b is INFINITY and has ! 949: bra Ld$infty | the opposite sign ! 950: ! 951: |============================================================================= ! 952: | __muldf3 ! 953: |============================================================================= ! 954: ! 955: | double __muldf3(double, double); ! 956: SYM (__muldf3): ! 957: link a6,IMM (0) ! 958: moveml d2-d7,sp@- ! 959: movel a6@(8),d0 | get a into d0-d1 ! 960: movel a6@(12),d1 | ! 961: movel a6@(16),d2 | and b into d2-d3 ! 962: movel a6@(20),d3 | ! 963: movel d0,d7 | d7 will hold the sign of the product ! 964: eorl d2,d7 | ! 965: andl IMM (0x80000000),d7 | ! 966: movel d7,a0 | save sign bit into a0 ! 967: movel IMM (0x7ff00000),d7 | useful constant (+INFINITY) ! 968: movel d7,d6 | another (mask for fraction) ! 969: notl d6 | ! 970: bclr IMM (31),d0 | get rid of a's sign bit ' ! 971: movel d0,d4 | ! 972: orl d1,d4 | ! 973: beq Lmuldf$a$0 | branch if a is zero ! 974: movel d0,d4 | ! 975: bclr IMM (31),d2 | get rid of b's sign bit ' ! 976: movel d2,d5 | ! 977: orl d3,d5 | ! 978: beq Lmuldf$b$0 | branch if b is zero ! 979: movel d2,d5 | ! 980: cmpl d7,d0 | is a big? ! 981: bhi Lmuldf$inop | if a is NaN return NaN ! 982: beq Lmuldf$a$nf | we still have to check d1 and b ... ! 983: cmpl d7,d2 | now compare b with INFINITY ! 984: bhi Lmuldf$inop | is b NaN? ! 985: beq Lmuldf$b$nf | we still have to check d3 ... ! 986: | Here we have both numbers finite and nonzero (and with no sign bit). ! 987: | Now we get the exponents into d4 and d5. ! 988: andl d7,d4 | isolate exponent in d4 ! 989: beq Lmuldf$a$den | if exponent zero, have denormalized ! 990: andl d6,d0 | isolate fraction ! 991: orl IMM (0x00100000),d0 | and put hidden bit back ! 992: swap d4 | I like exponents in the first byte ! 993: lsrw IMM (4),d4 | ! 994: Lmuldf$1: ! 995: andl d7,d5 | ! 996: beq Lmuldf$b$den | ! 997: andl d6,d2 | ! 998: orl IMM (0x00100000),d2 | and put hidden bit back ! 999: swap d5 | ! 1000: lsrw IMM (4),d5 | ! 1001: Lmuldf$2: | ! 1002: addw d5,d4 | add exponents ! 1003: subw IMM (D_BIAS+1),d4 | and subtract bias (plus one) ! 1004: ! 1005: | We are now ready to do the multiplication. The situation is as follows: ! 1006: | both a and b have bit 52 ( bit 20 of d0 and d2) set (even if they were ! 1007: | denormalized to start with!), which means that in the product bit 104 ! 1008: | (which will correspond to bit 8 of the fourth long) is set. ! 1009: ! 1010: | Here we have to do the product. ! 1011: | To do it we have to juggle the registers back and forth, as there are not ! 1012: | enough to keep everything in them. So we use the address registers to keep ! 1013: | some intermediate data. ! 1014: ! 1015: moveml a2-a3,sp@- | save a2 and a3 for temporary use ! 1016: movel IMM (0),a2 | a2 is a null register ! 1017: movel d4,a3 | and a3 will preserve the exponent ! 1018: ! 1019: | First, shift d2-d3 so bit 20 becomes bit 31: ! 1020: rorl IMM (5),d2 | rotate d2 5 places right ! 1021: swap d2 | and swap it ! 1022: rorl IMM (5),d3 | do the same thing with d3 ! 1023: swap d3 | ! 1024: movew d3,d6 | get the rightmost 11 bits of d3 ! 1025: andw IMM (0x07ff),d6 | ! 1026: orw d6,d2 | and put them into d2 ! 1027: andw IMM (0xf800),d3 | clear those bits in d3 ! 1028: ! 1029: movel d2,d6 | move b into d6-d7 ! 1030: movel d3,d7 | move a into d4-d5 ! 1031: movel d0,d4 | and clear d0-d1-d2-d3 (to put result) ! 1032: movel d1,d5 | ! 1033: movel IMM (0),d3 | ! 1034: movel d3,d2 | ! 1035: movel d3,d1 | ! 1036: movel d3,d0 | ! 1037: ! 1038: | We use a1 as counter: ! 1039: movel IMM (DBL_MANT_DIG-1),a1 ! 1040: exg d7,a1 ! 1041: ! 1042: 1: exg d7,a1 | put counter back in a1 ! 1043: addl d3,d3 | shift sum once left ! 1044: addxl d2,d2 | ! 1045: addxl d1,d1 | ! 1046: addxl d0,d0 | ! 1047: addl d7,d7 | ! 1048: addxl d6,d6 | ! 1049: bcc 2f | if bit clear skip the following ! 1050: exg d7,a2 | ! 1051: addl d5,d3 | else add a to the sum ! 1052: addxl d4,d2 | ! 1053: addxl d7,d1 | ! 1054: addxl d7,d0 | ! 1055: exg d7,a2 | ! 1056: 2: exg d7,a1 | put counter in d7 ! 1057: dbf d7,1b | decrement and branch ! 1058: ! 1059: movel a3,d4 | restore exponent ! 1060: moveml sp@+,a2-a3 ! 1061: ! 1062: | Now we have the product in d0-d1-d2-d3, with bit 8 of d0 set. The ! 1063: | first thing to do now is to normalize it so bit 8 becomes bit ! 1064: | DBL_MANT_DIG-32 (to do the rounding); later we will shift right. ! 1065: swap d0 ! 1066: swap d1 ! 1067: movew d1,d0 ! 1068: swap d2 ! 1069: movew d2,d1 ! 1070: swap d3 ! 1071: movew d3,d2 ! 1072: movew IMM (0),d3 ! 1073: lsrl IMM (1),d0 ! 1074: roxrl IMM (1),d1 ! 1075: roxrl IMM (1),d2 ! 1076: roxrl IMM (1),d3 ! 1077: lsrl IMM (1),d0 ! 1078: roxrl IMM (1),d1 ! 1079: roxrl IMM (1),d2 ! 1080: roxrl IMM (1),d3 ! 1081: lsrl IMM (1),d0 ! 1082: roxrl IMM (1),d1 ! 1083: roxrl IMM (1),d2 ! 1084: roxrl IMM (1),d3 ! 1085: ! 1086: | Now round, check for over- and underflow, and exit. ! 1087: movel a0,d7 | get sign bit back into d7 ! 1088: movew IMM (MULTIPLY),d5 ! 1089: ! 1090: btst IMM (DBL_MANT_DIG+1-32),d0 ! 1091: beq Lround$exit ! 1092: lsrl IMM (1),d0 ! 1093: roxrl IMM (1),d1 ! 1094: addw IMM (1),d4 ! 1095: bra Lround$exit ! 1096: ! 1097: Lmuldf$inop: ! 1098: movew IMM (MULTIPLY),d5 ! 1099: bra Ld$inop ! 1100: ! 1101: Lmuldf$b$nf: ! 1102: movew IMM (MULTIPLY),d5 ! 1103: movel a0,d7 | get sign bit back into d7 ! 1104: tstl d3 | we know d2 == 0x7ff00000, so check d3 ! 1105: bne Ld$inop | if d3 <> 0 b is NaN ! 1106: bra Ld$overflow | else we have overflow (since a is finite) ! 1107: ! 1108: Lmuldf$a$nf: ! 1109: movew IMM (MULTIPLY),d5 ! 1110: movel a0,d7 | get sign bit back into d7 ! 1111: tstl d1 | we know d0 == 0x7ff00000, so check d1 ! 1112: bne Ld$inop | if d1 <> 0 a is NaN ! 1113: bra Ld$overflow | else signal overflow ! 1114: ! 1115: | If either number is zero return zero, unless the other is +/-INFINITY or ! 1116: | NaN, in which case we return NaN. ! 1117: Lmuldf$b$0: ! 1118: movew IMM (MULTIPLY),d5 ! 1119: exg d2,d0 | put b (==0) into d0-d1 ! 1120: exg d3,d1 | and a (with sign bit cleared) into d2-d3 ! 1121: bra 1f ! 1122: Lmuldf$a$0: ! 1123: movel a6@(16),d2 | put b into d2-d3 again ! 1124: movel a6@(20),d3 | ! 1125: bclr IMM (31),d2 | clear sign bit ! 1126: 1: cmpl IMM (0x7ff00000),d2 | check for non-finiteness ! 1127: bge Ld$inop | in case NaN or +/-INFINITY return NaN ! 1128: lea SYM (_fpCCR),a0 ! 1129: movew IMM (0),a0@ ! 1130: moveml sp@+,d2-d7 ! 1131: unlk a6 ! 1132: rts ! 1133: ! 1134: | If a number is denormalized we put an exponent of 1 but do not put the ! 1135: | hidden bit back into the fraction; instead we shift left until bit 21 ! 1136: | (the hidden bit) is set, adjusting the exponent accordingly. We do this ! 1137: | to ensure that the product of the fractions is close to 1. ! 1138: Lmuldf$a$den: ! 1139: movel IMM (1),d4 ! 1140: andl d6,d0 ! 1141: 1: addl d1,d1 | shift a left until bit 20 is set ! 1142: addxl d0,d0 | ! 1143: subw IMM (1),d4 | and adjust exponent ! 1144: btst IMM (20),d0 | ! 1145: bne Lmuldf$1 | ! 1146: bra 1b ! 1147: ! 1148: Lmuldf$b$den: ! 1149: movel IMM (1),d5 ! 1150: andl d6,d2 ! 1151: 1: addl d3,d3 | shift b left until bit 20 is set ! 1152: addxl d2,d2 | ! 1153: subw IMM (1),d5 | and adjust exponent ! 1154: btst IMM (20),d2 | ! 1155: bne Lmuldf$2 | ! 1156: bra 1b ! 1157: ! 1158: ! 1159: |============================================================================= ! 1160: | __divdf3 ! 1161: |============================================================================= ! 1162: ! 1163: | double __divdf3(double, double); ! 1164: SYM (__divdf3): ! 1165: link a6,IMM (0) ! 1166: moveml d2-d7,sp@- ! 1167: movel a6@(8),d0 | get a into d0-d1 ! 1168: movel a6@(12),d1 | ! 1169: movel a6@(16),d2 | and b into d2-d3 ! 1170: movel a6@(20),d3 | ! 1171: movel d0,d7 | d7 will hold the sign of the result ! 1172: eorl d2,d7 | ! 1173: andl IMM (0x80000000),d7 ! 1174: movel d7,a0 | save sign into a0 ! 1175: movel IMM (0x7ff00000),d7 | useful constant (+INFINITY) ! 1176: movel d7,d6 | another (mask for fraction) ! 1177: notl d6 | ! 1178: bclr IMM (31),d0 | get rid of a's sign bit ' ! 1179: movel d0,d4 | ! 1180: orl d1,d4 | ! 1181: beq Ldivdf$a$0 | branch if a is zero ! 1182: movel d0,d4 | ! 1183: bclr IMM (31),d2 | get rid of b's sign bit ' ! 1184: movel d2,d5 | ! 1185: orl d3,d5 | ! 1186: beq Ldivdf$b$0 | branch if b is zero ! 1187: movel d2,d5 ! 1188: cmpl d7,d0 | is a big? ! 1189: bhi Ldivdf$inop | if a is NaN return NaN ! 1190: beq Ldivdf$a$nf | if d0 == 0x7ff00000 we check d1 ! 1191: cmpl d7,d2 | now compare b with INFINITY ! 1192: bhi Ldivdf$inop | if b is NaN return NaN ! 1193: beq Ldivdf$b$nf | if d2 == 0x7ff00000 we check d3 ! 1194: | Here we have both numbers finite and nonzero (and with no sign bit). ! 1195: | Now we get the exponents into d4 and d5 and normalize the numbers to ! 1196: | ensure that the ratio of the fractions is around 1. We do this by ! 1197: | making sure that both numbers have bit #DBL_MANT_DIG-32-1 (hidden bit) ! 1198: | set, even if they were denormalized to start with. ! 1199: | Thus, the result will satisfy: 2 > result > 1/2. ! 1200: andl d7,d4 | and isolate exponent in d4 ! 1201: beq Ldivdf$a$den | if exponent is zero we have a denormalized ! 1202: andl d6,d0 | and isolate fraction ! 1203: orl IMM (0x00100000),d0 | and put hidden bit back ! 1204: swap d4 | I like exponents in the first byte ! 1205: lsrw IMM (4),d4 | ! 1206: Ldivdf$1: | ! 1207: andl d7,d5 | ! 1208: beq Ldivdf$b$den | ! 1209: andl d6,d2 | ! 1210: orl IMM (0x00100000),d2 ! 1211: swap d5 | ! 1212: lsrw IMM (4),d5 | ! 1213: Ldivdf$2: | ! 1214: subw d5,d4 | subtract exponents ! 1215: addw IMM (D_BIAS),d4 | and add bias ! 1216: ! 1217: | We are now ready to do the division. We have prepared things in such a way ! 1218: | that the ratio of the fractions will be less than 2 but greater than 1/2. ! 1219: | At this point the registers in use are: ! 1220: | d0-d1 hold a (first operand, bit DBL_MANT_DIG-32=0, bit ! 1221: | DBL_MANT_DIG-1-32=1) ! 1222: | d2-d3 hold b (second operand, bit DBL_MANT_DIG-32=1) ! 1223: | d4 holds the difference of the exponents, corrected by the bias ! 1224: | a0 holds the sign of the ratio ! 1225: ! 1226: | To do the rounding correctly we need to keep information about the ! 1227: | nonsignificant bits. One way to do this would be to do the division ! 1228: | using four registers; another is to use two registers (as originally ! 1229: | I did), but use a sticky bit to preserve information about the ! 1230: | fractional part. Note that we can keep that info in a1, which is not ! 1231: | used. ! 1232: movel IMM (0),d6 | d6-d7 will hold the result ! 1233: movel d6,d7 | ! 1234: movel IMM (0),a1 | and a1 will hold the sticky bit ! 1235: ! 1236: movel IMM (DBL_MANT_DIG-32+1),d5 ! 1237: ! 1238: 1: cmpl d0,d2 | is a < b? ! 1239: bhi 3f | if b > a skip the following ! 1240: beq 4f | if d0==d2 check d1 and d3 ! 1241: 2: subl d3,d1 | ! 1242: subxl d2,d0 | a <-- a - b ! 1243: bset d5,d6 | set the corresponding bit in d6 ! 1244: 3: addl d1,d1 | shift a by 1 ! 1245: addxl d0,d0 | ! 1246: dbra d5,1b | and branch back ! 1247: bra 5f ! 1248: 4: cmpl d1,d3 | here d0==d2, so check d1 and d3 ! 1249: bhi 3b | if d1 > d2 skip the subtraction ! 1250: bra 2b | else go do it ! 1251: 5: ! 1252: | Here we have to start setting the bits in the second long. ! 1253: movel IMM (31),d5 | again d5 is counter ! 1254: ! 1255: 1: cmpl d0,d2 | is a < b? ! 1256: bhi 3f | if b > a skip the following ! 1257: beq 4f | if d0==d2 check d1 and d3 ! 1258: 2: subl d3,d1 | ! 1259: subxl d2,d0 | a <-- a - b ! 1260: bset d5,d7 | set the corresponding bit in d7 ! 1261: 3: addl d1,d1 | shift a by 1 ! 1262: addxl d0,d0 | ! 1263: dbra d5,1b | and branch back ! 1264: bra 5f ! 1265: 4: cmpl d1,d3 | here d0==d2, so check d1 and d3 ! 1266: bhi 3b | if d1 > d2 skip the subtraction ! 1267: bra 2b | else go do it ! 1268: 5: ! 1269: | Now go ahead checking until we hit a one, which we store in d2. ! 1270: movel IMM (DBL_MANT_DIG),d5 ! 1271: 1: cmpl d2,d0 | is a < b? ! 1272: bhi 4f | if b < a, exit ! 1273: beq 3f | if d0==d2 check d1 and d3 ! 1274: 2: addl d1,d1 | shift a by 1 ! 1275: addxl d0,d0 | ! 1276: dbra d5,1b | and branch back ! 1277: movel IMM (0),d2 | here no sticky bit was found ! 1278: movel d2,d3 ! 1279: bra 5f ! 1280: 3: cmpl d1,d3 | here d0==d2, so check d1 and d3 ! 1281: bhi 2b | if d1 > d2 go back ! 1282: 4: ! 1283: | Here put the sticky bit in d2-d3 (in the position which actually corresponds ! 1284: | to it; if you don't do this the algorithm loses in some cases). ' ! 1285: movel IMM (0),d2 ! 1286: movel d2,d3 ! 1287: subw IMM (DBL_MANT_DIG),d5 ! 1288: addw IMM (63),d5 ! 1289: cmpw IMM (31),d5 ! 1290: bhi 2f ! 1291: 1: bset d5,d3 ! 1292: bra 5f ! 1293: subw IMM (32),d5 ! 1294: 2: bset d5,d2 ! 1295: 5: ! 1296: | Finally we are finished! Move the longs in the address registers to ! 1297: | their final destination: ! 1298: movel d6,d0 ! 1299: movel d7,d1 ! 1300: movel IMM (0),d3 ! 1301: ! 1302: | Here we have finished the division, with the result in d0-d1-d2-d3, with ! 1303: | 2^21 <= d6 < 2^23. Thus bit 23 is not set, but bit 22 could be set. ! 1304: | If it is not, then definitely bit 21 is set. Normalize so bit 22 is ! 1305: | not set: ! 1306: btst IMM (DBL_MANT_DIG-32+1),d0 ! 1307: beq 1f ! 1308: lsrl IMM (1),d0 ! 1309: roxrl IMM (1),d1 ! 1310: roxrl IMM (1),d2 ! 1311: roxrl IMM (1),d3 ! 1312: addw IMM (1),d4 ! 1313: 1: ! 1314: | Now round, check for over- and underflow, and exit. ! 1315: movel a0,d7 | restore sign bit to d7 ! 1316: movew IMM (DIVIDE),d5 ! 1317: bra Lround$exit ! 1318: ! 1319: Ldivdf$inop: ! 1320: movew IMM (DIVIDE),d5 ! 1321: bra Ld$inop ! 1322: ! 1323: Ldivdf$a$0: ! 1324: | If a is zero check to see whether b is zero also. In that case return ! 1325: | NaN; then check if b is NaN, and return NaN also in that case. Else ! 1326: | return zero. ! 1327: movew IMM (DIVIDE),d5 ! 1328: bclr IMM (31),d2 | ! 1329: movel d2,d4 | ! 1330: orl d3,d4 | ! 1331: beq Ld$inop | if b is also zero return NaN ! 1332: cmpl IMM (0x7ff00000),d2 | check for NaN ! 1333: bhi Ld$inop | ! 1334: blt 1f | ! 1335: tstl d3 | ! 1336: bne Ld$inop | ! 1337: 1: movel IMM (0),d0 | else return zero ! 1338: movel d0,d1 | ! 1339: lea SYM (_fpCCR),a0 | clear exception flags ! 1340: movew IMM (0),a0@ | ! 1341: moveml sp@+,d2-d7 | ! 1342: unlk a6 | ! 1343: rts | ! 1344: ! 1345: Ldivdf$b$0: ! 1346: movew IMM (DIVIDE),d5 ! 1347: | If we got here a is not zero. Check if a is NaN; in that case return NaN, ! 1348: | else return +/-INFINITY. Remember that a is in d0 with the sign bit ! 1349: | cleared already. ! 1350: movel a0,d7 | put a's sign bit back in d7 ' ! 1351: cmpl IMM (0x7ff00000),d0 | compare d0 with INFINITY ! 1352: bhi Ld$inop | if larger it is NaN ! 1353: tstl d1 | ! 1354: bne Ld$inop | ! 1355: bra Ld$div$0 | else signal DIVIDE_BY_ZERO ! 1356: ! 1357: Ldivdf$b$nf: ! 1358: movew IMM (DIVIDE),d5 ! 1359: | If d2 == 0x7ff00000 we have to check d3. ! 1360: tstl d3 | ! 1361: bne Ld$inop | if d3 <> 0, b is NaN ! 1362: bra Ld$underflow | else b is +/-INFINITY, so signal underflow ! 1363: ! 1364: Ldivdf$a$nf: ! 1365: movew IMM (DIVIDE),d5 ! 1366: | If d0 == 0x7ff00000 we have to check d1. ! 1367: tstl d1 | ! 1368: bne Ld$inop | if d1 <> 0, a is NaN ! 1369: | If a is INFINITY we have to check b ! 1370: cmpl d7,d2 | compare b with INFINITY ! 1371: bge Ld$inop | if b is NaN or INFINITY return NaN ! 1372: tstl d3 | ! 1373: bne Ld$inop | ! 1374: bra Ld$overflow | else return overflow ! 1375: ! 1376: | If a number is denormalized we put an exponent of 1 but do not put the ! 1377: | bit back into the fraction. ! 1378: Ldivdf$a$den: ! 1379: movel IMM (1),d4 ! 1380: andl d6,d0 ! 1381: 1: addl d1,d1 | shift a left until bit 20 is set ! 1382: addxl d0,d0 ! 1383: subw IMM (1),d4 | and adjust exponent ! 1384: btst IMM (DBL_MANT_DIG-32-1),d0 ! 1385: bne Ldivdf$1 ! 1386: bra 1b ! 1387: ! 1388: Ldivdf$b$den: ! 1389: movel IMM (1),d5 ! 1390: andl d6,d2 ! 1391: 1: addl d3,d3 | shift b left until bit 20 is set ! 1392: addxl d2,d2 ! 1393: subw IMM (1),d5 | and adjust exponent ! 1394: btst IMM (DBL_MANT_DIG-32-1),d2 ! 1395: bne Ldivdf$2 ! 1396: bra 1b ! 1397: ! 1398: Lround$exit: ! 1399: | This is a common exit point for __muldf3 and __divdf3. When they enter ! 1400: | this point the sign of the result is in d7, the result in d0-d1, normalized ! 1401: | so that 2^21 <= d0 < 2^22, and the exponent is in the lower byte of d4. ! 1402: ! 1403: | First check for underlow in the exponent: ! 1404: cmpw IMM (-DBL_MANT_DIG-1),d4 ! 1405: blt Ld$underflow ! 1406: | It could happen that the exponent is less than 1, in which case the ! 1407: | number is denormalized. In this case we shift right and adjust the ! 1408: | exponent until it becomes 1 or the fraction is zero (in the latter case ! 1409: | we signal underflow and return zero). ! 1410: movel d7,a0 | ! 1411: movel IMM (0),d6 | use d6-d7 to collect bits flushed right ! 1412: movel d6,d7 | use d6-d7 to collect bits flushed right ! 1413: cmpw IMM (1),d4 | if the exponent is less than 1 we ! 1414: bge 2f | have to shift right (denormalize) ! 1415: 1: addw IMM (1),d4 | adjust the exponent ! 1416: lsrl IMM (1),d0 | shift right once ! 1417: roxrl IMM (1),d1 | ! 1418: roxrl IMM (1),d2 | ! 1419: roxrl IMM (1),d3 | ! 1420: roxrl IMM (1),d6 | ! 1421: roxrl IMM (1),d7 | ! 1422: cmpw IMM (1),d4 | is the exponent 1 already? ! 1423: beq 2f | if not loop back ! 1424: bra 1b | ! 1425: bra Ld$underflow | safety check, shouldn't execute ' ! 1426: 2: orl d6,d2 | this is a trick so we don't lose ' ! 1427: orl d7,d3 | the bits which were flushed right ! 1428: movel a0,d7 | get back sign bit into d7 ! 1429: | Now call the rounding routine (which takes care of denormalized numbers): ! 1430: lea Lround$0,a0 | to return from rounding routine ! 1431: lea SYM (_fpCCR),a1 | check the rounding mode ! 1432: movew a1@(6),d6 | rounding mode in d6 ! 1433: beq Lround$to$nearest ! 1434: cmpw IMM (ROUND_TO_PLUS),d6 ! 1435: bhi Lround$to$minus ! 1436: blt Lround$to$zero ! 1437: bra Lround$to$plus ! 1438: Lround$0: ! 1439: | Here we have a correctly rounded result (either normalized or denormalized). ! 1440: ! 1441: | Here we should have either a normalized number or a denormalized one, and ! 1442: | the exponent is necessarily larger or equal to 1 (so we don't have to ' ! 1443: | check again for underflow!). We have to check for overflow or for a ! 1444: | denormalized number (which also signals underflow). ! 1445: | Check for overflow (i.e., exponent >= 0x7ff). ! 1446: cmpw IMM (0x07ff),d4 ! 1447: bge Ld$overflow ! 1448: | Now check for a denormalized number (exponent==0): ! 1449: movew d4,d4 ! 1450: beq Ld$den ! 1451: 1: ! 1452: | Put back the exponents and sign and return. ! 1453: lslw IMM (4),d4 | exponent back to fourth byte ! 1454: bclr IMM (DBL_MANT_DIG-32-1),d0 ! 1455: swap d0 | and put back exponent ! 1456: orw d4,d0 | ! 1457: swap d0 | ! 1458: orl d7,d0 | and sign also ! 1459: ! 1460: lea SYM (_fpCCR),a0 ! 1461: movew IMM (0),a0@ ! 1462: moveml sp@+,d2-d7 ! 1463: unlk a6 ! 1464: rts ! 1465: ! 1466: |============================================================================= ! 1467: | __negdf2 ! 1468: |============================================================================= ! 1469: ! 1470: | double __negdf2(double, double); ! 1471: SYM (__negdf2): ! 1472: link a6,IMM (0) ! 1473: moveml d2-d7,sp@- ! 1474: movew IMM (NEGATE),d5 ! 1475: movel a6@(8),d0 | get number to negate in d0-d1 ! 1476: movel a6@(12),d1 | ! 1477: bchg IMM (31),d0 | negate ! 1478: movel d0,d2 | make a positive copy (for the tests) ! 1479: bclr IMM (31),d2 | ! 1480: movel d2,d4 | check for zero ! 1481: orl d1,d4 | ! 1482: beq 2f | if zero (either sign) return +zero ! 1483: cmpl IMM (0x7ff00000),d2 | compare to +INFINITY ! 1484: blt 1f | if finite, return ! 1485: bhi Ld$inop | if larger (fraction not zero) is NaN ! 1486: tstl d1 | if d2 == 0x7ff00000 check d1 ! 1487: bne Ld$inop | ! 1488: movel d0,d7 | else get sign and return INFINITY ! 1489: andl IMM (0x80000000),d7 ! 1490: bra Ld$infty ! 1491: 1: lea SYM (_fpCCR),a0 ! 1492: movew IMM (0),a0@ ! 1493: moveml sp@+,d2-d7 ! 1494: unlk a6 ! 1495: rts ! 1496: 2: bclr IMM (31),d0 ! 1497: bra 1b ! 1498: ! 1499: |============================================================================= ! 1500: | __cmpdf2 ! 1501: |============================================================================= ! 1502: ! 1503: GREATER = 1 ! 1504: LESS = -1 ! 1505: EQUAL = 0 ! 1506: ! 1507: | int __cmpdf2(double, double); ! 1508: SYM (__cmpdf2): ! 1509: link a6,IMM (0) ! 1510: moveml d2-d7,sp@- | save registers ! 1511: movew IMM (COMPARE),d5 ! 1512: movel a6@(8),d0 | get first operand ! 1513: movel a6@(12),d1 | ! 1514: movel a6@(16),d2 | get second operand ! 1515: movel a6@(20),d3 | ! 1516: | First check if a and/or b are (+/-) zero and in that case clear ! 1517: | the sign bit. ! 1518: movel d0,d6 | copy signs into d6 (a) and d7(b) ! 1519: bclr IMM (31),d0 | and clear signs in d0 and d2 ! 1520: movel d2,d7 | ! 1521: bclr IMM (31),d2 | ! 1522: cmpl IMM (0x7fff0000),d0 | check for a == NaN ! 1523: bhi Ld$inop | if d0 > 0x7ff00000, a is NaN ! 1524: beq Lcmpdf$a$nf | if equal can be INFINITY, so check d1 ! 1525: movel d0,d4 | copy into d4 to test for zero ! 1526: orl d1,d4 | ! 1527: beq Lcmpdf$a$0 | ! 1528: Lcmpdf$0: ! 1529: cmpl IMM (0x7fff0000),d2 | check for b == NaN ! 1530: bhi Ld$inop | if d2 > 0x7ff00000, b is NaN ! 1531: beq Lcmpdf$b$nf | if equal can be INFINITY, so check d3 ! 1532: movel d2,d4 | ! 1533: orl d3,d4 | ! 1534: beq Lcmpdf$b$0 | ! 1535: Lcmpdf$1: ! 1536: | Check the signs ! 1537: eorl d6,d7 ! 1538: bpl 1f ! 1539: | If the signs are not equal check if a >= 0 ! 1540: tstl d6 ! 1541: bpl Lcmpdf$a$gt$b | if (a >= 0 && b < 0) => a > b ! 1542: bmi Lcmpdf$b$gt$a | if (a < 0 && b >= 0) => a < b ! 1543: 1: ! 1544: | If the signs are equal check for < 0 ! 1545: tstl d6 ! 1546: bpl 1f ! 1547: | If both are negative exchange them ! 1548: exg d0,d2 ! 1549: exg d1,d3 ! 1550: 1: ! 1551: | Now that they are positive we just compare them as longs (does this also ! 1552: | work for denormalized numbers?). ! 1553: cmpl d0,d2 ! 1554: bhi Lcmpdf$b$gt$a | |b| > |a| ! 1555: bne Lcmpdf$a$gt$b | |b| < |a| ! 1556: | If we got here d0 == d2, so we compare d1 and d3. ! 1557: cmpl d1,d3 ! 1558: bhi Lcmpdf$b$gt$a | |b| > |a| ! 1559: bne Lcmpdf$a$gt$b | |b| < |a| ! 1560: | If we got here a == b. ! 1561: movel IMM (EQUAL),d0 ! 1562: moveml sp@+,d2-d7 | put back the registers ! 1563: unlk a6 ! 1564: rts ! 1565: Lcmpdf$a$gt$b: ! 1566: movel IMM (GREATER),d0 ! 1567: moveml sp@+,d2-d7 | put back the registers ! 1568: unlk a6 ! 1569: rts ! 1570: Lcmpdf$b$gt$a: ! 1571: movel IMM (LESS),d0 ! 1572: moveml sp@+,d2-d7 | put back the registers ! 1573: unlk a6 ! 1574: rts ! 1575: ! 1576: Lcmpdf$a$0: ! 1577: bclr IMM (31),d6 ! 1578: bra Lcmpdf$0 ! 1579: Lcmpdf$b$0: ! 1580: bclr IMM (31),d7 ! 1581: bra Lcmpdf$1 ! 1582: ! 1583: Lcmpdf$a$nf: ! 1584: tstl d1 ! 1585: bne Ld$inop ! 1586: bra Lcmpdf$0 ! 1587: ! 1588: Lcmpdf$b$nf: ! 1589: tstl d3 ! 1590: bne Ld$inop ! 1591: bra Lcmpdf$1 ! 1592: ! 1593: |============================================================================= ! 1594: | rounding routines ! 1595: |============================================================================= ! 1596: ! 1597: | The rounding routines expect the number to be normalized in registers ! 1598: | d0-d1-d2-d3, with the exponent in register d4. They assume that the ! 1599: | exponent is larger or equal to 1. They return a properly normalized number ! 1600: | if possible, and a denormalized number otherwise. The exponent is returned ! 1601: | in d4. ! 1602: ! 1603: Lround$to$nearest: ! 1604: | We now normalize as suggested by D. Knuth ("Seminumerical Algorithms"): ! 1605: | Here we assume that the exponent is not too small (this should be checked ! 1606: | before entering the rounding routine), but the number could be denormalized. ! 1607: ! 1608: | Check for denormalized numbers: ! 1609: 1: btst IMM (DBL_MANT_DIG-32),d0 ! 1610: bne 2f | if set the number is normalized ! 1611: | Normalize shifting left until bit #DBL_MANT_DIG-32 is set or the exponent ! 1612: | is one (remember that a denormalized number corresponds to an ! 1613: | exponent of -D_BIAS+1). ! 1614: cmpw IMM (1),d4 | remember that the exponent is at least one ! 1615: beq 2f | an exponent of one means denormalized ! 1616: addl d3,d3 | else shift and adjust the exponent ! 1617: addxl d2,d2 | ! 1618: addxl d1,d1 | ! 1619: addxl d0,d0 | ! 1620: dbra d4,1b | ! 1621: 2: ! 1622: | Now round: we do it as follows: after the shifting we can write the ! 1623: | fraction part as f + delta, where 1 < f < 2^25, and 0 <= delta <= 2. ! 1624: | If delta < 1, do nothing. If delta > 1, add 1 to f. ! 1625: | If delta == 1, we make sure the rounded number will be even (odd?) ! 1626: | (after shifting). ! 1627: btst IMM (0),d1 | is delta < 1? ! 1628: beq 2f | if so, do not do anything ! 1629: orl d2,d3 | is delta == 1? ! 1630: bne 1f | if so round to even ! 1631: movel d1,d3 | ! 1632: andl IMM (2),d3 | bit 1 is the last significant bit ! 1633: movel IMM (0),d2 | ! 1634: addl d3,d1 | ! 1635: addxl d2,d0 | ! 1636: bra 2f | ! 1637: 1: movel IMM (1),d3 | else add 1 ! 1638: movel IMM (0),d2 | ! 1639: addl d3,d1 | ! 1640: addxl d2,d0 ! 1641: | Shift right once (because we used bit #DBL_MANT_DIG-32!). ! 1642: 2: lsrl IMM (1),d0 ! 1643: roxrl IMM (1),d1 ! 1644: ! 1645: | Now check again bit #DBL_MANT_DIG-32 (rounding could have produced a ! 1646: | 'fraction overflow' ...). ! 1647: btst IMM (DBL_MANT_DIG-32),d0 ! 1648: beq 1f ! 1649: lsrl IMM (1),d0 ! 1650: roxrl IMM (1),d1 ! 1651: addw IMM (1),d4 ! 1652: 1: ! 1653: | If bit #DBL_MANT_DIG-32-1 is clear we have a denormalized number, so we ! 1654: | have to put the exponent to zero and return a denormalized number. ! 1655: btst IMM (DBL_MANT_DIG-32-1),d0 ! 1656: beq 1f ! 1657: jmp a0@ ! 1658: 1: movel IMM (0),d4 ! 1659: jmp a0@ ! 1660: ! 1661: Lround$to$zero: ! 1662: Lround$to$plus: ! 1663: Lround$to$minus: ! 1664: jmp a0@ ! 1665: #endif /* L_double */ ! 1666: ! 1667: #ifdef L_float ! 1668: ! 1669: .globl SYM (_fpCCR) ! 1670: .globl $_exception_handler ! 1671: ! 1672: QUIET_NaN = 0xffffffff ! 1673: SIGNL_NaN = 0x7f800001 ! 1674: INFINITY = 0x7f800000 ! 1675: ! 1676: F_MAX_EXP = 0xff ! 1677: F_BIAS = 126 ! 1678: FLT_MAX_EXP = F_MAX_EXP - F_BIAS ! 1679: FLT_MIN_EXP = 1 - F_BIAS ! 1680: FLT_MANT_DIG = 24 ! 1681: ! 1682: INEXACT_RESULT = 0x0001 ! 1683: UNDERFLOW = 0x0002 ! 1684: OVERFLOW = 0x0004 ! 1685: DIVIDE_BY_ZERO = 0x0008 ! 1686: INVALID_OPERATION = 0x0010 ! 1687: ! 1688: SINGLE_FLOAT = 1 ! 1689: ! 1690: NOOP = 0 ! 1691: ADD = 1 ! 1692: MULTIPLY = 2 ! 1693: DIVIDE = 3 ! 1694: NEGATE = 4 ! 1695: COMPARE = 5 ! 1696: EXTENDSFDF = 6 ! 1697: TRUNCDFSF = 7 ! 1698: ! 1699: UNKNOWN = -1 ! 1700: ROUND_TO_NEAREST = 0 | round result to nearest representable value ! 1701: ROUND_TO_ZERO = 1 | round result towards zero ! 1702: ROUND_TO_PLUS = 2 | round result towards plus infinity ! 1703: ROUND_TO_MINUS = 3 | round result towards minus infinity ! 1704: ! 1705: | Entry points: ! 1706: ! 1707: .globl SYM (__addsf3) ! 1708: .globl SYM (__subsf3) ! 1709: .globl SYM (__mulsf3) ! 1710: .globl SYM (__divsf3) ! 1711: .globl SYM (__negsf2) ! 1712: .globl SYM (__cmpsf2) ! 1713: ! 1714: | These are common routines to return and signal exceptions. ! 1715: ! 1716: .text ! 1717: .even ! 1718: ! 1719: Lf$den: ! 1720: | Return and signal a denormalized number ! 1721: orl d7,d0 ! 1722: movew IMM (UNDERFLOW),d7 ! 1723: orw IMM (INEXACT_RESULT),d7 ! 1724: movew IMM (SINGLE_FLOAT),d6 ! 1725: jmp $_exception_handler ! 1726: ! 1727: Lf$infty: ! 1728: Lf$overflow: ! 1729: | Return a properly signed INFINITY and set the exception flags ! 1730: movel IMM (INFINITY),d0 ! 1731: orl d7,d0 ! 1732: movew IMM (OVERFLOW),d7 ! 1733: orw IMM (INEXACT_RESULT),d7 ! 1734: movew IMM (SINGLE_FLOAT),d6 ! 1735: jmp $_exception_handler ! 1736: ! 1737: Lf$underflow: ! 1738: | Return 0 and set the exception flags ! 1739: movel IMM (0),d0 ! 1740: movew IMM (UNDERFLOW),d7 ! 1741: orw IMM (INEXACT_RESULT),d7 ! 1742: movew IMM (SINGLE_FLOAT),d6 ! 1743: jmp $_exception_handler ! 1744: ! 1745: Lf$inop: ! 1746: | Return a quiet NaN and set the exception flags ! 1747: movel IMM (QUIET_NaN),d0 ! 1748: movew IMM (INVALID_OPERATION),d7 ! 1749: orw IMM (INEXACT_RESULT),d7 ! 1750: movew IMM (SINGLE_FLOAT),d6 ! 1751: jmp $_exception_handler ! 1752: ! 1753: Lf$div$0: ! 1754: | Return a properly signed INFINITY and set the exception flags ! 1755: movel IMM (INFINITY),d0 ! 1756: orl d7,d0 ! 1757: movew IMM (DIVIDE_BY_ZERO),d7 ! 1758: orw IMM (INEXACT_RESULT),d7 ! 1759: movew IMM (SINGLE_FLOAT),d6 ! 1760: jmp $_exception_handler ! 1761: ! 1762: |============================================================================= ! 1763: |============================================================================= ! 1764: | single precision routines ! 1765: |============================================================================= ! 1766: |============================================================================= ! 1767: ! 1768: | A single precision floating point number (float) has the format: ! 1769: | ! 1770: | struct _float { ! 1771: | unsigned int sign : 1; /* sign bit */ ! 1772: | unsigned int exponent : 8; /* exponent, shifted by 126 */ ! 1773: | unsigned int fraction : 23; /* fraction */ ! 1774: | } float; ! 1775: | ! 1776: | Thus sizeof(float) = 4 (32 bits). ! 1777: | ! 1778: | All the routines are callable from C programs, and return the result ! 1779: | in the single register d0. They also preserve all registers except ! 1780: | d0-d1 and a0-a1. ! 1781: ! 1782: |============================================================================= ! 1783: | __subsf3 ! 1784: |============================================================================= ! 1785: ! 1786: | float __subsf3(float, float); ! 1787: SYM (__subsf3): ! 1788: bchg IMM (31),sp@(8) | change sign of second operand ! 1789: | and fall through ! 1790: |============================================================================= ! 1791: | __addsf3 ! 1792: |============================================================================= ! 1793: ! 1794: | float __addsf3(float, float); ! 1795: SYM (__addsf3): ! 1796: link a6,IMM (0) | everything will be done in registers ! 1797: moveml d2-d7,sp@- | save all data registers but d0-d1 ! 1798: movel a6@(8),d0 | get first operand ! 1799: movel a6@(12),d1 | get second operand ! 1800: movel d0,d6 | get d0's sign bit ' ! 1801: addl d0,d0 | check and clear sign bit of a ! 1802: beq Laddsf$b | if zero return second operand ! 1803: movel d1,d7 | save b's sign bit ' ! 1804: addl d1,d1 | get rid of sign bit ! 1805: beq Laddsf$a | if zero return first operand ! 1806: ! 1807: movel d6,a0 | save signs in address registers ! 1808: movel d7,a1 | so we can use d6 and d7 ! 1809: ! 1810: | Get the exponents and check for denormalized and/or infinity. ! 1811: ! 1812: movel IMM (0x00ffffff),d4 | mask to get fraction ! 1813: movel IMM (0x01000000),d5 | mask to put hidden bit back ! 1814: ! 1815: movel d0,d6 | save a to get exponent ! 1816: andl d4,d0 | get fraction in d0 ! 1817: notl d4 | make d4 into a mask for the exponent ! 1818: andl d4,d6 | get exponent in d6 ! 1819: beq Laddsf$a$den | branch if a is denormalized ! 1820: cmpl d4,d6 | check for INFINITY or NaN ! 1821: beq Laddsf$nf ! 1822: swap d6 | put exponent into first word ! 1823: orl d5,d0 | and put hidden bit back ! 1824: Laddsf$1: ! 1825: | Now we have a's exponent in d6 (second byte) and the mantissa in d0. ' ! 1826: movel d1,d7 | get exponent in d7 ! 1827: andl d4,d7 | ! 1828: beq Laddsf$b$den | branch if b is denormalized ! 1829: cmpl d4,d7 | check for INFINITY or NaN ! 1830: beq Laddsf$nf ! 1831: swap d7 | put exponent into first word ! 1832: notl d4 | make d4 into a mask for the fraction ! 1833: andl d4,d1 | get fraction in d1 ! 1834: orl d5,d1 | and put hidden bit back ! 1835: Laddsf$2: ! 1836: | Now we have b's exponent in d7 (second byte) and the mantissa in d1. ' ! 1837: ! 1838: | Note that the hidden bit corresponds to bit #FLT_MANT_DIG-1, and we ! 1839: | shifted right once, so bit #FLT_MANT_DIG is set (so we have one extra ! 1840: | bit). ! 1841: ! 1842: movel d1,d2 | move b to d2, since we want to use ! 1843: | two registers to do the sum ! 1844: movel IMM (0),d1 | and clear the new ones ! 1845: movel d1,d3 | ! 1846: ! 1847: | Here we shift the numbers in registers d0 and d1 so the exponents are the ! 1848: | same, and put the largest exponent in d6. Note that we are using two ! 1849: | registers for each number (see the discussion by D. Knuth in "Seminumerical ! 1850: | Algorithms"). ! 1851: cmpw d6,d7 | compare exponents ! 1852: beq Laddsf$3 | if equal don't shift ' ! 1853: bhi 5f | branch if second exponent largest ! 1854: 1: ! 1855: subl d6,d7 | keep the largest exponent ! 1856: negl d7 ! 1857: lsrw IMM (8),d7 | put difference in lower byte ! 1858: | if difference is too large we don't shift (actually, we can just exit) ' ! 1859: cmpw IMM (FLT_MANT_DIG+2),d7 ! 1860: bge Laddsf$b$small ! 1861: cmpw IMM (16),d7 | if difference >= 16 swap ! 1862: bge 4f ! 1863: 2: ! 1864: subw IMM (1),d7 ! 1865: 3: lsrl IMM (1),d2 | shift right second operand ! 1866: roxrl IMM (1),d3 ! 1867: dbra d7,3b ! 1868: bra Laddsf$3 ! 1869: 4: ! 1870: movew d2,d3 ! 1871: swap d3 ! 1872: movew d3,d2 ! 1873: swap d2 ! 1874: subw IMM (16),d7 ! 1875: bne 2b | if still more bits, go back to normal case ! 1876: bra Laddsf$3 ! 1877: 5: ! 1878: exg d6,d7 | exchange the exponents ! 1879: subl d6,d7 | keep the largest exponent ! 1880: negl d7 | ! 1881: lsrw IMM (8),d7 | put difference in lower byte ! 1882: | if difference is too large we don't shift (and exit!) ' ! 1883: cmpw IMM (FLT_MANT_DIG+2),d7 ! 1884: bge Laddsf$a$small ! 1885: cmpw IMM (16),d7 | if difference >= 16 swap ! 1886: bge 8f ! 1887: 6: ! 1888: subw IMM (1),d7 ! 1889: 7: lsrl IMM (1),d0 | shift right first operand ! 1890: roxrl IMM (1),d1 ! 1891: dbra d7,7b ! 1892: bra Laddsf$3 ! 1893: 8: ! 1894: movew d0,d1 ! 1895: swap d1 ! 1896: movew d1,d0 ! 1897: swap d0 ! 1898: subw IMM (16),d7 ! 1899: bne 6b | if still more bits, go back to normal case ! 1900: | otherwise we fall through ! 1901: ! 1902: | Now we have a in d0-d1, b in d2-d3, and the largest exponent in d6 (the ! 1903: | signs are stored in a0 and a1). ! 1904: ! 1905: Laddsf$3: ! 1906: | Here we have to decide whether to add or subtract the numbers ! 1907: exg d6,a0 | get signs back ! 1908: exg d7,a1 | and save the exponents ! 1909: eorl d6,d7 | combine sign bits ! 1910: bmi Lsubsf$0 | if negative a and b have opposite ! 1911: | sign so we actually subtract the ! 1912: | numbers ! 1913: ! 1914: | Here we have both positive or both negative ! 1915: exg d6,a0 | now we have the exponent in d6 ! 1916: movel a0,d7 | and sign in d7 ! 1917: andl IMM (0x80000000),d7 ! 1918: | Here we do the addition. ! 1919: addl d3,d1 ! 1920: addxl d2,d0 ! 1921: | Note: now we have d2, d3, d4 and d5 to play with! ! 1922: ! 1923: | Put the exponent, in the first byte, in d2, to use the "standard" rounding ! 1924: | routines: ! 1925: movel d6,d2 ! 1926: lsrw IMM (8),d2 ! 1927: ! 1928: | Before rounding normalize so bit #FLT_MANT_DIG is set (we will consider ! 1929: | the case of denormalized numbers in the rounding routine itself). ! 1930: | As in the addition (not in the subtraction!) we could have set ! 1931: | one more bit we check this: ! 1932: btst IMM (FLT_MANT_DIG+1),d0 ! 1933: beq 1f ! 1934: lsrl IMM (1),d0 ! 1935: roxrl IMM (1),d1 ! 1936: addl IMM (1),d2 ! 1937: 1: ! 1938: lea Laddsf$4,a0 | to return from rounding routine ! 1939: lea SYM (_fpCCR),a1 | check the rounding mode ! 1940: movew a1@(6),d6 | rounding mode in d6 ! 1941: beq Lround$to$nearest ! 1942: cmpw IMM (ROUND_TO_PLUS),d6 ! 1943: bhi Lround$to$minus ! 1944: blt Lround$to$zero ! 1945: bra Lround$to$plus ! 1946: Laddsf$4: ! 1947: | Put back the exponent, but check for overflow. ! 1948: cmpw IMM (0xff),d2 ! 1949: bhi 1f ! 1950: bclr IMM (FLT_MANT_DIG-1),d0 ! 1951: lslw IMM (7),d2 ! 1952: swap d2 ! 1953: orl d2,d0 ! 1954: bra Laddsf$ret ! 1955: 1: ! 1956: movew IMM (ADD),d5 ! 1957: bra Lf$overflow ! 1958: ! 1959: Lsubsf$0: ! 1960: | We are here if a > 0 and b < 0 (sign bits cleared). ! 1961: | Here we do the subtraction. ! 1962: movel d6,d7 | put sign in d7 ! 1963: andl IMM (0x80000000),d7 ! 1964: ! 1965: subl d3,d1 | result in d0-d1 ! 1966: subxl d2,d0 | ! 1967: beq Laddsf$ret | if zero just exit ! 1968: bpl 1f | if positive skip the following ! 1969: bchg IMM (31),d7 | change sign bit in d7 ! 1970: negl d1 ! 1971: negxl d0 ! 1972: 1: ! 1973: exg d2,a0 | now we have the exponent in d2 ! 1974: lsrw IMM (8),d2 | put it in the first byte ! 1975: ! 1976: | Now d0-d1 is positive and the sign bit is in d7. ! 1977: ! 1978: | Note that we do not have to normalize, since in the subtraction bit ! 1979: | #FLT_MANT_DIG+1 is never set, and denormalized numbers are handled by ! 1980: | the rounding routines themselves. ! 1981: lea Lsubsf$1,a0 | to return from rounding routine ! 1982: lea SYM (_fpCCR),a1 | check the rounding mode ! 1983: movew a1@(6),d6 | rounding mode in d6 ! 1984: beq Lround$to$nearest ! 1985: cmpw IMM (ROUND_TO_PLUS),d6 ! 1986: bhi Lround$to$minus ! 1987: blt Lround$to$zero ! 1988: bra Lround$to$plus ! 1989: Lsubsf$1: ! 1990: | Put back the exponent (we can't have overflow!). ' ! 1991: bclr IMM (FLT_MANT_DIG-1),d0 ! 1992: lslw IMM (7),d2 ! 1993: swap d2 ! 1994: orl d2,d0 ! 1995: bra Laddsf$ret ! 1996: ! 1997: | If one of the numbers was too small (difference of exponents >= ! 1998: | FLT_MANT_DIG+2) we return the other (and now we don't have to ' ! 1999: | check for finiteness or zero). ! 2000: Laddsf$a$small: ! 2001: movel a6@(12),d0 ! 2002: lea SYM (_fpCCR),a0 ! 2003: movew IMM (0),a0@ ! 2004: moveml sp@+,d2-d7 | restore data registers ! 2005: unlk a6 | and return ! 2006: rts ! 2007: ! 2008: Laddsf$b$small: ! 2009: movel a6@(8),d0 ! 2010: lea SYM (_fpCCR),a0 ! 2011: movew IMM (0),a0@ ! 2012: moveml sp@+,d2-d7 | restore data registers ! 2013: unlk a6 | and return ! 2014: rts ! 2015: ! 2016: | If the numbers are denormalized remember to put exponent equal to 1. ! 2017: ! 2018: Laddsf$a$den: ! 2019: movel d5,d6 | d5 contains 0x01000000 ! 2020: swap d6 ! 2021: bra Laddsf$1 ! 2022: ! 2023: Laddsf$b$den: ! 2024: movel d5,d7 ! 2025: swap d7 ! 2026: notl d4 | make d4 into a mask for the fraction ! 2027: | (this was not executed after the jump) ! 2028: bra Laddsf$2 ! 2029: ! 2030: | The rest is mainly code for the different results which can be ! 2031: | returned (checking always for +/-INFINITY and NaN). ! 2032: ! 2033: Laddsf$b: ! 2034: | Return b (if a is zero). ! 2035: movel a6@(12),d0 ! 2036: bra 1f ! 2037: Laddsf$a: ! 2038: | Return a (if b is zero). ! 2039: movel a6@(8),d0 ! 2040: 1: ! 2041: movew IMM (ADD),d5 ! 2042: | We have to check for NaN and +/-infty. ! 2043: movel d0,d7 ! 2044: andl IMM (0x80000000),d7 | put sign in d7 ! 2045: bclr IMM (31),d0 | clear sign ! 2046: cmpl IMM (INFINITY),d0 | check for infty or NaN ! 2047: bge 2f ! 2048: movel d0,d0 | check for zero (we do this because we don't ' ! 2049: bne Laddsf$ret | want to return -0 by mistake ! 2050: bclr IMM (31),d7 | if zero be sure to clear sign ! 2051: bra Laddsf$ret | if everything OK just return ! 2052: 2: ! 2053: | The value to be returned is either +/-infty or NaN ! 2054: andl IMM (0x007fffff),d0 | check for NaN ! 2055: bne Lf$inop | if mantissa not zero is NaN ! 2056: bra Lf$infty ! 2057: ! 2058: Laddsf$ret: ! 2059: | Normal exit (a and b nonzero, result is not NaN nor +/-infty). ! 2060: | We have to clear the exception flags (just the exception type). ! 2061: lea SYM (_fpCCR),a0 ! 2062: movew IMM (0),a0@ ! 2063: orl d7,d0 | put sign bit ! 2064: moveml sp@+,d2-d7 | restore data registers ! 2065: unlk a6 | and return ! 2066: rts ! 2067: ! 2068: Laddsf$ret$den: ! 2069: | Return a denormalized number (for addition we don't signal underflow) ' ! 2070: lsrl IMM (1),d0 | remember to shift right back once ! 2071: bra Laddsf$ret | and return ! 2072: ! 2073: | Note: when adding two floats of the same sign if either one is ! 2074: | NaN we return NaN without regard to whether the other is finite or ! 2075: | not. When subtracting them (i.e., when adding two numbers of ! 2076: | opposite signs) things are more complicated: if both are INFINITY ! 2077: | we return NaN, if only one is INFINITY and the other is NaN we return ! 2078: | NaN, but if it is finite we return INFINITY with the corresponding sign. ! 2079: ! 2080: Laddsf$nf: ! 2081: movew IMM (ADD),d5 ! 2082: | This could be faster but it is not worth the effort, since it is not ! 2083: | executed very often. We sacrifice speed for clarity here. ! 2084: movel a6@(8),d0 | get the numbers back (remember that we ! 2085: movel a6@(12),d1 | did some processing already) ! 2086: movel IMM (INFINITY),d4 | useful constant (INFINITY) ! 2087: movel d0,d2 | save sign bits ! 2088: movel d1,d3 ! 2089: bclr IMM (31),d0 | clear sign bits ! 2090: bclr IMM (31),d1 ! 2091: | We know that one of them is either NaN of +/-INFINITY ! 2092: | Check for NaN (if either one is NaN return NaN) ! 2093: cmpl d4,d0 | check first a (d0) ! 2094: bhi Lf$inop ! 2095: cmpl d4,d1 | check now b (d1) ! 2096: bhi Lf$inop ! 2097: | Now comes the check for +/-INFINITY. We know that both are (maybe not ! 2098: | finite) numbers, but we have to check if both are infinite whether we ! 2099: | are adding or subtracting them. ! 2100: eorl d3,d2 | to check sign bits ! 2101: bmi 1f ! 2102: movel d0,d7 ! 2103: andl IMM (0x80000000),d7 | get (common) sign bit ! 2104: bra Lf$infty ! 2105: 1: ! 2106: | We know one (or both) are infinite, so we test for equality between the ! 2107: | two numbers (if they are equal they have to be infinite both, so we ! 2108: | return NaN). ! 2109: cmpl d1,d0 | are both infinite? ! 2110: beq Lf$inop | if so return NaN ! 2111: ! 2112: movel d0,d7 ! 2113: andl IMM (0x80000000),d7 | get a's sign bit ' ! 2114: cmpl d4,d0 | test now for infinity ! 2115: beq Lf$infty | if a is INFINITY return with this sign ! 2116: bchg IMM (31),d7 | else we know b is INFINITY and has ! 2117: bra Lf$infty | the opposite sign ! 2118: ! 2119: |============================================================================= ! 2120: | __mulsf3 ! 2121: |============================================================================= ! 2122: ! 2123: | float __mulsf3(float, float); ! 2124: SYM (__mulsf3): ! 2125: link a6,IMM (0) ! 2126: moveml d2-d7,sp@- ! 2127: movel a6@(8),d0 | get a into d0 ! 2128: movel a6@(12),d1 | and b into d1 ! 2129: movel d0,d7 | d7 will hold the sign of the product ! 2130: eorl d1,d7 | ! 2131: andl IMM (0x80000000),d7 ! 2132: movel IMM (INFINITY),d6 | useful constant (+INFINITY) ! 2133: movel d6,d5 | another (mask for fraction) ! 2134: notl d5 | ! 2135: movel IMM (0x00800000),d4 | this is to put hidden bit back ! 2136: bclr IMM (31),d0 | get rid of a's sign bit ' ! 2137: movel d0,d2 | ! 2138: beq Lmulsf$a$0 | branch if a is zero ! 2139: bclr IMM (31),d1 | get rid of b's sign bit ' ! 2140: movel d1,d3 | ! 2141: beq Lmulsf$b$0 | branch if b is zero ! 2142: cmpl d6,d0 | is a big? ! 2143: bhi Lmulsf$inop | if a is NaN return NaN ! 2144: beq Lmulsf$inf | if a is INFINITY we have to check b ! 2145: cmpl d6,d1 | now compare b with INFINITY ! 2146: bhi Lmulsf$inop | is b NaN? ! 2147: beq Lmulsf$overflow | is b INFINITY? ! 2148: | Here we have both numbers finite and nonzero (and with no sign bit). ! 2149: | Now we get the exponents into d2 and d3. ! 2150: andl d6,d2 | and isolate exponent in d2 ! 2151: beq Lmulsf$a$den | if exponent is zero we have a denormalized ! 2152: andl d5,d0 | and isolate fraction ! 2153: orl d4,d0 | and put hidden bit back ! 2154: swap d2 | I like exponents in the first byte ! 2155: lsrw IMM (7),d2 | ! 2156: Lmulsf$1: | number ! 2157: andl d6,d3 | ! 2158: beq Lmulsf$b$den | ! 2159: andl d5,d1 | ! 2160: orl d4,d1 | ! 2161: swap d3 | ! 2162: lsrw IMM (7),d3 | ! 2163: Lmulsf$2: | ! 2164: addw d3,d2 | add exponents ! 2165: subw IMM (F_BIAS+1),d2 | and subtract bias (plus one) ! 2166: ! 2167: | We are now ready to do the multiplication. The situation is as follows: ! 2168: | both a and b have bit FLT_MANT_DIG-1 set (even if they were ! 2169: | denormalized to start with!), which means that in the product ! 2170: | bit 2*(FLT_MANT_DIG-1) (that is, bit 2*FLT_MANT_DIG-2-32 of the ! 2171: | high long) is set. ! 2172: ! 2173: | To do the multiplication let us move the number a little bit around ... ! 2174: movel d1,d6 | second operand in d6 ! 2175: movel d0,d5 | first operand in d4-d5 ! 2176: movel IMM (0),d4 ! 2177: movel d4,d1 | the sums will go in d0-d1 ! 2178: movel d4,d0 ! 2179: ! 2180: | now bit FLT_MANT_DIG-1 becomes bit 31: ! 2181: lsll IMM (31-FLT_MANT_DIG+1),d6 ! 2182: ! 2183: | Start the loop (we loop #FLT_MANT_DIG times): ! 2184: movew IMM (FLT_MANT_DIG-1),d3 ! 2185: 1: addl d1,d1 | shift sum ! 2186: addxl d0,d0 ! 2187: lsll IMM (1),d6 | get bit bn ! 2188: bcc 2f | if not set skip sum ! 2189: addl d5,d1 | add a ! 2190: addxl d4,d0 ! 2191: 2: dbf d3,1b | loop back ! 2192: ! 2193: | Now we have the product in d0-d1, with bit (FLT_MANT_DIG - 1) + FLT_MANT_DIG ! 2194: | (mod 32) of d0 set. The first thing to do now is to normalize it so bit ! 2195: | FLT_MANT_DIG is set (to do the rounding). ! 2196: rorl IMM (6),d1 ! 2197: swap d1 ! 2198: movew d1,d3 ! 2199: andw IMM (0x03ff),d3 ! 2200: andw IMM (0xfd00),d1 ! 2201: lsll IMM (8),d0 ! 2202: addl d0,d0 ! 2203: addl d0,d0 ! 2204: orw d3,d0 ! 2205: ! 2206: movew IMM (MULTIPLY),d5 ! 2207: ! 2208: btst IMM (FLT_MANT_DIG+1),d0 ! 2209: beq Lround$exit ! 2210: lsrl IMM (1),d0 ! 2211: roxrl IMM (1),d1 ! 2212: addw IMM (1),d2 ! 2213: bra Lround$exit ! 2214: ! 2215: Lmulsf$inop: ! 2216: movew IMM (MULTIPLY),d5 ! 2217: bra Lf$inop ! 2218: ! 2219: Lmulsf$overflow: ! 2220: movew IMM (MULTIPLY),d5 ! 2221: bra Lf$overflow ! 2222: ! 2223: Lmulsf$inf: ! 2224: movew IMM (MULTIPLY),d5 ! 2225: | If either is NaN return NaN; else both are (maybe infinite) numbers, so ! 2226: | return INFINITY with the correct sign (which is in d7). ! 2227: cmpl d6,d1 | is b NaN? ! 2228: bhi Lf$inop | if so return NaN ! 2229: bra Lf$overflow | else return +/-INFINITY ! 2230: ! 2231: | If either number is zero return zero, unless the other is +/-INFINITY, ! 2232: | or NaN, in which case we return NaN. ! 2233: Lmulsf$b$0: ! 2234: | Here d1 (==b) is zero. ! 2235: movel d1,d0 | put b into d0 (just a zero) ! 2236: movel a6@(8),d1 | get a again to check for non-finiteness ! 2237: bra 1f ! 2238: Lmulsf$a$0: ! 2239: movel a6@(12),d1 | get b again to check for non-finiteness ! 2240: 1: bclr IMM (31),d1 | clear sign bit ! 2241: cmpl IMM (INFINITY),d1 | and check for a large exponent ! 2242: bge Lf$inop | if b is +/-INFINITY or NaN return NaN ! 2243: lea SYM (_fpCCR),a0 | else return zero ! 2244: movew IMM (0),a0@ | ! 2245: moveml sp@+,d2-d7 | ! 2246: unlk a6 | ! 2247: rts | ! 2248: ! 2249: | If a number is denormalized we put an exponent of 1 but do not put the ! 2250: | hidden bit back into the fraction; instead we shift left until bit 23 ! 2251: | (the hidden bit) is set, adjusting the exponent accordingly. We do this ! 2252: | to ensure that the product of the fractions is close to 1. ! 2253: Lmulsf$a$den: ! 2254: movel IMM (1),d2 ! 2255: andl d5,d0 ! 2256: 1: addl d0,d0 | shift a left (until bit 23 is set) ! 2257: subw IMM (1),d2 | and adjust exponent ! 2258: btst IMM (FLT_MANT_DIG-1),d0 ! 2259: bne Lmulsf$1 | ! 2260: bra 1b | else loop back ! 2261: ! 2262: Lmulsf$b$den: ! 2263: movel IMM (1),d3 ! 2264: andl d5,d1 ! 2265: 1: addl d1,d1 | shift b left until bit 23 is set ! 2266: subw IMM (1),d3 | and adjust exponent ! 2267: btst IMM (FLT_MANT_DIG-1),d1 ! 2268: bne Lmulsf$2 | ! 2269: bra 1b | else loop back ! 2270: ! 2271: |============================================================================= ! 2272: | __divsf3 ! 2273: |============================================================================= ! 2274: ! 2275: | float __divsf3(float, float); ! 2276: SYM (__divsf3): ! 2277: link a6,IMM (0) ! 2278: moveml d2-d7,sp@- ! 2279: movel a6@(8),d0 | get a into d0 ! 2280: movel a6@(12),d1 | and b into d1 ! 2281: movel d0,d7 | d7 will hold the sign of the result ! 2282: eorl d1,d7 | ! 2283: andl IMM (0x80000000),d7 | ! 2284: movel IMM (INFINITY),d6 | useful constant (+INFINITY) ! 2285: movel d6,d5 | another (mask for fraction) ! 2286: notl d5 | ! 2287: movel IMM (0x00800000),d4 | this is to put hidden bit back ! 2288: bclr IMM (31),d0 | get rid of a's sign bit ' ! 2289: movel d0,d2 | ! 2290: beq Ldivsf$a$0 | branch if a is zero ! 2291: bclr IMM (31),d1 | get rid of b's sign bit ' ! 2292: movel d1,d3 | ! 2293: beq Ldivsf$b$0 | branch if b is zero ! 2294: cmpl d6,d0 | is a big? ! 2295: bhi Ldivsf$inop | if a is NaN return NaN ! 2296: beq Ldivsf$inf | if a is INFINITY we have to check b ! 2297: cmpl d6,d1 | now compare b with INFINITY ! 2298: bhi Ldivsf$inop | if b is NaN return NaN ! 2299: beq Ldivsf$underflow ! 2300: | Here we have both numbers finite and nonzero (and with no sign bit). ! 2301: | Now we get the exponents into d2 and d3 and normalize the numbers to ! 2302: | ensure that the ratio of the fractions is close to 1. We do this by ! 2303: | making sure that bit #FLT_MANT_DIG-1 (hidden bit) is set. ! 2304: andl d6,d2 | and isolate exponent in d2 ! 2305: beq Ldivsf$a$den | if exponent is zero we have a denormalized ! 2306: andl d5,d0 | and isolate fraction ! 2307: orl d4,d0 | and put hidden bit back ! 2308: swap d2 | I like exponents in the first byte ! 2309: lsrw IMM (7),d2 | ! 2310: Ldivsf$1: | ! 2311: andl d6,d3 | ! 2312: beq Ldivsf$b$den | ! 2313: andl d5,d1 | ! 2314: orl d4,d1 | ! 2315: swap d3 | ! 2316: lsrw IMM (7),d3 | ! 2317: Ldivsf$2: | ! 2318: subw d3,d2 | subtract exponents ! 2319: addw IMM (F_BIAS),d2 | and add bias ! 2320: ! 2321: | We are now ready to do the division. We have prepared things in such a way ! 2322: | that the ratio of the fractions will be less than 2 but greater than 1/2. ! 2323: | At this point the registers in use are: ! 2324: | d0 holds a (first operand, bit FLT_MANT_DIG=0, bit FLT_MANT_DIG-1=1) ! 2325: | d1 holds b (second operand, bit FLT_MANT_DIG=1) ! 2326: | d2 holds the difference of the exponents, corrected by the bias ! 2327: | d7 holds the sign of the ratio ! 2328: | d4, d5, d6 hold some constants ! 2329: movel d7,a0 | d6-d7 will hold the ratio of the fractions ! 2330: movel IMM (0),d6 | ! 2331: movel d6,d7 ! 2332: ! 2333: movew IMM (FLT_MANT_DIG+1),d3 ! 2334: 1: cmpl d0,d1 | is a < b? ! 2335: bhi 2f | ! 2336: bset d3,d6 | set a bit in d6 ! 2337: subl d1,d0 | if a >= b a <-- a-b ! 2338: beq 3f | if a is zero, exit ! 2339: 2: addl d0,d0 | multiply a by 2 ! 2340: dbra d3,1b ! 2341: ! 2342: | Now we keep going to set the sticky bit ... ! 2343: movew IMM (FLT_MANT_DIG),d3 ! 2344: 1: cmpl d0,d1 ! 2345: ble 2f ! 2346: addl d0,d0 ! 2347: dbra d3,1b ! 2348: movel IMM (0),d1 ! 2349: bra 3f ! 2350: 2: movel IMM (0),d1 ! 2351: subw IMM (FLT_MANT_DIG),d3 ! 2352: addw IMM (31),d3 ! 2353: bset d3,d1 ! 2354: 3: ! 2355: movel d6,d0 | put the ratio in d0-d1 ! 2356: movel a0,d7 | get sign back ! 2357: ! 2358: | Because of the normalization we did before we are guaranteed that ! 2359: | d0 is smaller than 2^26 but larger than 2^24. Thus bit 26 is not set, ! 2360: | bit 25 could be set, and if it is not set then bit 24 is necessarily set. ! 2361: btst IMM (FLT_MANT_DIG+1),d0 ! 2362: beq 1f | if it is not set, then bit 24 is set ! 2363: lsrl IMM (1),d0 | ! 2364: addw IMM (1),d2 | ! 2365: 1: ! 2366: | Now round, check for over- and underflow, and exit. ! 2367: movew IMM (DIVIDE),d5 ! 2368: bra Lround$exit ! 2369: ! 2370: Ldivsf$inop: ! 2371: movew IMM (DIVIDE),d5 ! 2372: bra Lf$inop ! 2373: ! 2374: Ldivsf$overflow: ! 2375: movew IMM (DIVIDE),d5 ! 2376: bra Lf$overflow ! 2377: ! 2378: Ldivsf$underflow: ! 2379: movew IMM (DIVIDE),d5 ! 2380: bra Lf$underflow ! 2381: ! 2382: Ldivsf$a$0: ! 2383: movew IMM (DIVIDE),d5 ! 2384: | If a is zero check to see whether b is zero also. In that case return ! 2385: | NaN; then check if b is NaN, and return NaN also in that case. Else ! 2386: | return zero. ! 2387: andl IMM (0x7fffffff),d1 | clear sign bit and test b ! 2388: beq Lf$inop | if b is also zero return NaN ! 2389: cmpl IMM (INFINITY),d1 | check for NaN ! 2390: bhi Lf$inop | ! 2391: movel IMM (0),d0 | else return zero ! 2392: lea SYM (_fpCCR),a0 | ! 2393: movew IMM (0),a0@ | ! 2394: moveml sp@+,d2-d7 | ! 2395: unlk a6 | ! 2396: rts | ! 2397: ! 2398: Ldivsf$b$0: ! 2399: movew IMM (DIVIDE),d5 ! 2400: | If we got here a is not zero. Check if a is NaN; in that case return NaN, ! 2401: | else return +/-INFINITY. Remember that a is in d0 with the sign bit ! 2402: | cleared already. ! 2403: cmpl IMM (INFINITY),d0 | compare d0 with INFINITY ! 2404: bhi Lf$inop | if larger it is NaN ! 2405: bra Lf$div$0 | else signal DIVIDE_BY_ZERO ! 2406: ! 2407: Ldivsf$inf: ! 2408: movew IMM (DIVIDE),d5 ! 2409: | If a is INFINITY we have to check b ! 2410: cmpl IMM (INFINITY),d1 | compare b with INFINITY ! 2411: bge Lf$inop | if b is NaN or INFINITY return NaN ! 2412: bra Lf$overflow | else return overflow ! 2413: ! 2414: | If a number is denormalized we put an exponent of 1 but do not put the ! 2415: | bit back into the fraction. ! 2416: Ldivsf$a$den: ! 2417: movel IMM (1),d2 ! 2418: andl d5,d0 ! 2419: 1: addl d0,d0 | shift a left until bit FLT_MANT_DIG-1 is set ! 2420: subw IMM (1),d2 | and adjust exponent ! 2421: btst IMM (FLT_MANT_DIG-1),d0 ! 2422: bne Ldivsf$1 ! 2423: bra 1b ! 2424: ! 2425: Ldivsf$b$den: ! 2426: movel IMM (1),d3 ! 2427: andl d5,d1 ! 2428: 1: addl d1,d1 | shift b left until bit FLT_MANT_DIG is set ! 2429: subw IMM (1),d3 | and adjust exponent ! 2430: btst IMM (FLT_MANT_DIG-1),d1 ! 2431: bne Ldivsf$2 ! 2432: bra 1b ! 2433: ! 2434: Lround$exit: ! 2435: | This is a common exit point for __mulsf3 and __divsf3. ! 2436: ! 2437: | First check for underlow in the exponent: ! 2438: cmpw IMM (-FLT_MANT_DIG-1),d2 ! 2439: blt Lf$underflow ! 2440: | It could happen that the exponent is less than 1, in which case the ! 2441: | number is denormalized. In this case we shift right and adjust the ! 2442: | exponent until it becomes 1 or the fraction is zero (in the latter case ! 2443: | we signal underflow and return zero). ! 2444: movel IMM (0),d6 | d6 is used temporarily ! 2445: cmpw IMM (1),d2 | if the exponent is less than 1 we ! 2446: bge 2f | have to shift right (denormalize) ! 2447: 1: addw IMM (1),d2 | adjust the exponent ! 2448: lsrl IMM (1),d0 | shift right once ! 2449: roxrl IMM (1),d1 | ! 2450: roxrl IMM (1),d6 | d6 collect bits we would lose otherwise ! 2451: cmpw IMM (1),d2 | is the exponent 1 already? ! 2452: beq 2f | if not loop back ! 2453: bra 1b | ! 2454: bra Lf$underflow | safety check, shouldn't execute ' ! 2455: 2: orl d6,d1 | this is a trick so we don't lose ' ! 2456: | the extra bits which were flushed right ! 2457: | Now call the rounding routine (which takes care of denormalized numbers): ! 2458: lea Lround$0,a0 | to return from rounding routine ! 2459: lea SYM (_fpCCR),a1 | check the rounding mode ! 2460: movew a1@(6),d6 | rounding mode in d6 ! 2461: beq Lround$to$nearest ! 2462: cmpw IMM (ROUND_TO_PLUS),d6 ! 2463: bhi Lround$to$minus ! 2464: blt Lround$to$zero ! 2465: bra Lround$to$plus ! 2466: Lround$0: ! 2467: | Here we have a correctly rounded result (either normalized or denormalized). ! 2468: ! 2469: | Here we should have either a normalized number or a denormalized one, and ! 2470: | the exponent is necessarily larger or equal to 1 (so we don't have to ' ! 2471: | check again for underflow!). We have to check for overflow or for a ! 2472: | denormalized number (which also signals underflow). ! 2473: | Check for overflow (i.e., exponent >= 255). ! 2474: cmpw IMM (0x00ff),d2 ! 2475: bge Lf$overflow ! 2476: | Now check for a denormalized number (exponent==0). ! 2477: movew d2,d2 ! 2478: beq Lf$den ! 2479: 1: ! 2480: | Put back the exponents and sign and return. ! 2481: lslw IMM (7),d2 | exponent back to fourth byte ! 2482: bclr IMM (FLT_MANT_DIG-1),d0 ! 2483: swap d0 | and put back exponent ! 2484: orw d2,d0 | ! 2485: swap d0 | ! 2486: orl d7,d0 | and sign also ! 2487: ! 2488: lea SYM (_fpCCR),a0 ! 2489: movew IMM (0),a0@ ! 2490: moveml sp@+,d2-d7 ! 2491: unlk a6 ! 2492: rts ! 2493: ! 2494: |============================================================================= ! 2495: | __negsf2 ! 2496: |============================================================================= ! 2497: ! 2498: | This is trivial and could be shorter if we didn't bother checking for NaN ' ! 2499: | and +/-INFINITY. ! 2500: ! 2501: | float __negsf2(float); ! 2502: SYM (__negsf2): ! 2503: link a6,IMM (0) ! 2504: moveml d2-d7,sp@- ! 2505: movew IMM (NEGATE),d5 ! 2506: movel a6@(8),d0 | get number to negate in d0 ! 2507: bchg IMM (31),d0 | negate ! 2508: movel d0,d1 | make a positive copy ! 2509: bclr IMM (31),d1 | ! 2510: tstl d1 | check for zero ! 2511: beq 2f | if zero (either sign) return +zero ! 2512: cmpl IMM (INFINITY),d1 | compare to +INFINITY ! 2513: blt 1f | ! 2514: bhi Lf$inop | if larger (fraction not zero) is NaN ! 2515: movel d0,d7 | else get sign and return INFINITY ! 2516: andl IMM (0x80000000),d7 ! 2517: bra Lf$infty ! 2518: 1: lea SYM (_fpCCR),a0 ! 2519: movew IMM (0),a0@ ! 2520: moveml sp@+,d2-d7 ! 2521: unlk a6 ! 2522: rts ! 2523: 2: bclr IMM (31),d0 ! 2524: bra 1b ! 2525: ! 2526: |============================================================================= ! 2527: | __cmpsf2 ! 2528: |============================================================================= ! 2529: ! 2530: GREATER = 1 ! 2531: LESS = -1 ! 2532: EQUAL = 0 ! 2533: ! 2534: | int __cmpsf2(float, float); ! 2535: SYM (__cmpsf2): ! 2536: link a6,IMM (0) ! 2537: moveml d2-d7,sp@- | save registers ! 2538: movew IMM (COMPARE),d5 ! 2539: movel a6@(8),d0 | get first operand ! 2540: movel a6@(12),d1 | get second operand ! 2541: | Check if either is NaN, and in that case return garbage and signal ! 2542: | INVALID_OPERATION. Check also if either is zero, and clear the signs ! 2543: | if necessary. ! 2544: movel d0,d6 ! 2545: andl IMM (0x7fffffff),d0 ! 2546: beq Lcmpsf$a$0 ! 2547: cmpl IMM (0x7f800000),d0 ! 2548: bhi Lf$inop ! 2549: Lcmpsf$1: ! 2550: movel d1,d7 ! 2551: andl IMM (0x7fffffff),d1 ! 2552: beq Lcmpsf$b$0 ! 2553: cmpl IMM (0x7f800000),d1 ! 2554: bhi Lf$inop ! 2555: Lcmpsf$2: ! 2556: | Check the signs ! 2557: eorl d6,d7 ! 2558: bpl 1f ! 2559: | If the signs are not equal check if a >= 0 ! 2560: tstl d6 ! 2561: bpl Lcmpsf$a$gt$b | if (a >= 0 && b < 0) => a > b ! 2562: bmi Lcmpsf$b$gt$a | if (a < 0 && b >= 0) => a < b ! 2563: 1: ! 2564: | If the signs are equal check for < 0 ! 2565: tstl d6 ! 2566: bpl 1f ! 2567: | If both are negative exchange them ! 2568: exg d0,d1 ! 2569: 1: ! 2570: | Now that they are positive we just compare them as longs (does this also ! 2571: | work for denormalized numbers?). ! 2572: cmpl d0,d1 ! 2573: bhi Lcmpsf$b$gt$a | |b| > |a| ! 2574: bne Lcmpsf$a$gt$b | |b| < |a| ! 2575: | If we got here a == b. ! 2576: movel IMM (EQUAL),d0 ! 2577: moveml sp@+,d2-d7 | put back the registers ! 2578: unlk a6 ! 2579: rts ! 2580: Lcmpsf$a$gt$b: ! 2581: movel IMM (GREATER),d0 ! 2582: moveml sp@+,d2-d7 | put back the registers ! 2583: unlk a6 ! 2584: rts ! 2585: Lcmpsf$b$gt$a: ! 2586: movel IMM (LESS),d0 ! 2587: moveml sp@+,d2-d7 | put back the registers ! 2588: unlk a6 ! 2589: rts ! 2590: ! 2591: Lcmpsf$a$0: ! 2592: bclr IMM (31),d6 ! 2593: bra Lcmpsf$1 ! 2594: Lcmpsf$b$0: ! 2595: bclr IMM (31),d7 ! 2596: bra Lcmpsf$2 ! 2597: ! 2598: |============================================================================= ! 2599: | rounding routines ! 2600: |============================================================================= ! 2601: ! 2602: | The rounding routines expect the number to be normalized in registers ! 2603: | d0-d1, with the exponent in register d2. They assume that the ! 2604: | exponent is larger or equal to 1. They return a properly normalized number ! 2605: | if possible, and a denormalized number otherwise. The exponent is returned ! 2606: | in d2. ! 2607: ! 2608: Lround$to$nearest: ! 2609: | We now normalize as suggested by D. Knuth ("Seminumerical Algorithms"): ! 2610: | Here we assume that the exponent is not too small (this should be checked ! 2611: | before entering the rounding routine), but the number could be denormalized. ! 2612: ! 2613: | Check for denormalized numbers: ! 2614: 1: btst IMM (FLT_MANT_DIG),d0 ! 2615: bne 2f | if set the number is normalized ! 2616: | Normalize shifting left until bit #FLT_MANT_DIG is set or the exponent ! 2617: | is one (remember that a denormalized number corresponds to an ! 2618: | exponent of -F_BIAS+1). ! 2619: cmpw IMM (1),d2 | remember that the exponent is at least one ! 2620: beq 2f | an exponent of one means denormalized ! 2621: addl d1,d1 | else shift and adjust the exponent ! 2622: addxl d0,d0 | ! 2623: dbra d2,1b | ! 2624: 2: ! 2625: | Now round: we do it as follows: after the shifting we can write the ! 2626: | fraction part as f + delta, where 1 < f < 2^25, and 0 <= delta <= 2. ! 2627: | If delta < 1, do nothing. If delta > 1, add 1 to f. ! 2628: | If delta == 1, we make sure the rounded number will be even (odd?) ! 2629: | (after shifting). ! 2630: btst IMM (0),d0 | is delta < 1? ! 2631: beq 2f | if so, do not do anything ! 2632: tstl d1 | is delta == 1? ! 2633: bne 1f | if so round to even ! 2634: movel d0,d1 | ! 2635: andl IMM (2),d1 | bit 1 is the last significant bit ! 2636: addl d1,d0 | ! 2637: bra 2f | ! 2638: 1: movel IMM (1),d1 | else add 1 ! 2639: addl d1,d0 | ! 2640: | Shift right once (because we used bit #FLT_MANT_DIG!). ! 2641: 2: lsrl IMM (1),d0 ! 2642: | Now check again bit #FLT_MANT_DIG (rounding could have produced a ! 2643: | 'fraction overflow' ...). ! 2644: btst IMM (FLT_MANT_DIG),d0 ! 2645: beq 1f ! 2646: lsrl IMM (1),d0 ! 2647: addw IMM (1),d2 ! 2648: 1: ! 2649: | If bit #FLT_MANT_DIG-1 is clear we have a denormalized number, so we ! 2650: | have to put the exponent to zero and return a denormalized number. ! 2651: btst IMM (FLT_MANT_DIG-1),d0 ! 2652: beq 1f ! 2653: jmp a0@ ! 2654: 1: movel IMM (0),d2 ! 2655: jmp a0@ ! 2656: ! 2657: Lround$to$zero: ! 2658: Lround$to$plus: ! 2659: Lround$to$minus: ! 2660: jmp a0@ ! 2661: #endif /* L_float */ ! 2662: ! 2663: | gcc expects the routines __eqdf2, __nedf2, __gtdf2, __gedf2, ! 2664: | __ledf2, __ltdf2 to all return the same value as a direct call to ! 2665: | __cmpdf2 would. In this implementation, each of these routines ! 2666: | simply calls __cmpdf2. It would be more efficient to give the ! 2667: | __cmpdf2 routine several names, but separating them out will make it ! 2668: | easier to write efficient versions of these routines someday. ! 2669: ! 2670: #ifdef L_eqdf2 ! 2671: LL0: ! 2672: .text ! 2673: .proc ! 2674: |#PROC# 04 ! 2675: LF18 = 4 ! 2676: LS18 = 128 ! 2677: LFF18 = 0 ! 2678: LSS18 = 0 ! 2679: LV18 = 0 ! 2680: .text ! 2681: .globl SYM (__eqdf2) ! 2682: SYM (__eqdf2): ! 2683: |#PROLOGUE# 0 ! 2684: link a6,IMM (0) ! 2685: |#PROLOGUE# 1 ! 2686: movl a6@(20),sp@- ! 2687: movl a6@(16),sp@- ! 2688: movl a6@(12),sp@- ! 2689: movl a6@(8),sp@- ! 2690: jbsr SYM (__cmpdf2) ! 2691: |#PROLOGUE# 2 ! 2692: unlk a6 ! 2693: |#PROLOGUE# 3 ! 2694: rts ! 2695: #endif /* L_eqdf2 */ ! 2696: ! 2697: #ifdef L_nedf2 ! 2698: LL0: ! 2699: .text ! 2700: .proc ! 2701: |#PROC# 04 ! 2702: LF18 = 8 ! 2703: LS18 = 132 ! 2704: LFF18 = 0 ! 2705: LSS18 = 0 ! 2706: LV18 = 0 ! 2707: .text ! 2708: .globl SYM (__nedf2) ! 2709: SYM (__nedf2): ! 2710: |#PROLOGUE# 0 ! 2711: link a6,IMM (0) ! 2712: |#PROLOGUE# 1 ! 2713: movl a6@(20),sp@- ! 2714: movl a6@(16),sp@- ! 2715: movl a6@(12),sp@- ! 2716: movl a6@(8),sp@- ! 2717: jbsr SYM (__cmpdf2) ! 2718: |#PROLOGUE# 2 ! 2719: unlk a6 ! 2720: |#PROLOGUE# 3 ! 2721: rts ! 2722: #endif /* L_nedf2 */ ! 2723: ! 2724: #ifdef L_gtdf2 ! 2725: .text ! 2726: .proc ! 2727: |#PROC# 04 ! 2728: LF18 = 8 ! 2729: LS18 = 132 ! 2730: LFF18 = 0 ! 2731: LSS18 = 0 ! 2732: LV18 = 0 ! 2733: .text ! 2734: .globl SYM (__gtdf2) ! 2735: SYM (__gtdf2): ! 2736: |#PROLOGUE# 0 ! 2737: link a6,IMM (0) ! 2738: |#PROLOGUE# 1 ! 2739: movl a6@(20),sp@- ! 2740: movl a6@(16),sp@- ! 2741: movl a6@(12),sp@- ! 2742: movl a6@(8),sp@- ! 2743: jbsr SYM (__cmpdf2) ! 2744: |#PROLOGUE# 2 ! 2745: unlk a6 ! 2746: |#PROLOGUE# 3 ! 2747: rts ! 2748: #endif /* L_gtdf2 */ ! 2749: ! 2750: #ifdef L_gedf2 ! 2751: LL0: ! 2752: .text ! 2753: .proc ! 2754: |#PROC# 04 ! 2755: LF18 = 8 ! 2756: LS18 = 132 ! 2757: LFF18 = 0 ! 2758: LSS18 = 0 ! 2759: LV18 = 0 ! 2760: .text ! 2761: .globl SYM (__gedf2) ! 2762: SYM (__gedf2): ! 2763: |#PROLOGUE# 0 ! 2764: link a6,IMM (0) ! 2765: |#PROLOGUE# 1 ! 2766: movl a6@(20),sp@- ! 2767: movl a6@(16),sp@- ! 2768: movl a6@(12),sp@- ! 2769: movl a6@(8),sp@- ! 2770: jbsr SYM (__cmpdf2) ! 2771: |#PROLOGUE# 2 ! 2772: unlk a6 ! 2773: |#PROLOGUE# 3 ! 2774: rts ! 2775: #endif /* L_gedf2 */ ! 2776: ! 2777: #ifdef L_ltdf2 ! 2778: LL0: ! 2779: .text ! 2780: .proc ! 2781: |#PROC# 04 ! 2782: LF18 = 8 ! 2783: LS18 = 132 ! 2784: LFF18 = 0 ! 2785: LSS18 = 0 ! 2786: LV18 = 0 ! 2787: .text ! 2788: .globl SYM (__ltdf2) ! 2789: SYM (__ltdf2): ! 2790: |#PROLOGUE# 0 ! 2791: link a6,IMM (0) ! 2792: |#PROLOGUE# 1 ! 2793: movl a6@(20),sp@- ! 2794: movl a6@(16),sp@- ! 2795: movl a6@(12),sp@- ! 2796: movl a6@(8),sp@- ! 2797: jbsr SYM (__cmpdf2) ! 2798: |#PROLOGUE# 2 ! 2799: unlk a6 ! 2800: |#PROLOGUE# 3 ! 2801: rts ! 2802: #endif /* L_ltdf2 */ ! 2803: ! 2804: #ifdef L_ledf2 ! 2805: .text ! 2806: .proc ! 2807: |#PROC# 04 ! 2808: LF18 = 8 ! 2809: LS18 = 132 ! 2810: LFF18 = 0 ! 2811: LSS18 = 0 ! 2812: LV18 = 0 ! 2813: .text ! 2814: .globl SYM (__ledf2) ! 2815: SYM (__ledf2): ! 2816: |#PROLOGUE# 0 ! 2817: link a6,IMM (0) ! 2818: |#PROLOGUE# 1 ! 2819: movl a6@(20),sp@- ! 2820: movl a6@(16),sp@- ! 2821: movl a6@(12),sp@- ! 2822: movl a6@(8),sp@- ! 2823: jbsr SYM (__cmpdf2) ! 2824: |#PROLOGUE# 2 ! 2825: unlk a6 ! 2826: |#PROLOGUE# 3 ! 2827: rts ! 2828: #endif /* L_ledf2 */ ! 2829: ! 2830: | The comments above about __eqdf2, et. al., also apply to __eqsf2, ! 2831: | et. al., except that the latter call __cmpsf2 rather than __cmpdf2. ! 2832: ! 2833: #ifdef L_eqsf2 ! 2834: .text ! 2835: .proc ! 2836: |#PROC# 04 ! 2837: LF18 = 4 ! 2838: LS18 = 128 ! 2839: LFF18 = 0 ! 2840: LSS18 = 0 ! 2841: LV18 = 0 ! 2842: .text ! 2843: .globl SYM (__eqsf2) ! 2844: SYM (__eqsf2): ! 2845: |#PROLOGUE# 0 ! 2846: link a6,IMM (0) ! 2847: |#PROLOGUE# 1 ! 2848: movl a6@(12),sp@- ! 2849: movl a6@(8),sp@- ! 2850: jbsr SYM (__cmpsf2) ! 2851: |#PROLOGUE# 2 ! 2852: unlk a6 ! 2853: |#PROLOGUE# 3 ! 2854: rts ! 2855: #endif /* L_eqsf2 */ ! 2856: ! 2857: #ifdef L_nesf2 ! 2858: .text ! 2859: .proc ! 2860: |#PROC# 04 ! 2861: LF18 = 8 ! 2862: LS18 = 132 ! 2863: LFF18 = 0 ! 2864: LSS18 = 0 ! 2865: LV18 = 0 ! 2866: .text ! 2867: .globl SYM (__nesf2) ! 2868: SYM (__nesf2): ! 2869: |#PROLOGUE# 0 ! 2870: link a6,IMM (0) ! 2871: |#PROLOGUE# 1 ! 2872: movl a6@(12),sp@- ! 2873: movl a6@(8),sp@- ! 2874: jbsr SYM (__cmpsf2) ! 2875: |#PROLOGUE# 2 ! 2876: unlk a6 ! 2877: |#PROLOGUE# 3 ! 2878: rts ! 2879: #endif /* L_nesf2 */ ! 2880: ! 2881: #ifdef L_gtsf2 ! 2882: .text ! 2883: .proc ! 2884: |#PROC# 04 ! 2885: LF18 = 8 ! 2886: LS18 = 132 ! 2887: LFF18 = 0 ! 2888: LSS18 = 0 ! 2889: LV18 = 0 ! 2890: .text ! 2891: .globl SYM (__gtsf2) ! 2892: SYM (__gtsf2): ! 2893: |#PROLOGUE# 0 ! 2894: link a6,IMM (0) ! 2895: |#PROLOGUE# 1 ! 2896: movl a6@(12),sp@- ! 2897: movl a6@(8),sp@- ! 2898: jbsr SYM (__cmpsf2) ! 2899: |#PROLOGUE# 2 ! 2900: unlk a6 ! 2901: |#PROLOGUE# 3 ! 2902: rts ! 2903: #endif /* L_gtsf2 */ ! 2904: ! 2905: #ifdef L_gesf2 ! 2906: .text ! 2907: .proc ! 2908: |#PROC# 04 ! 2909: LF18 = 8 ! 2910: LS18 = 132 ! 2911: LFF18 = 0 ! 2912: LSS18 = 0 ! 2913: LV18 = 0 ! 2914: .text ! 2915: .globl SYM (__gesf2) ! 2916: SYM (__gesf2): ! 2917: |#PROLOGUE# 0 ! 2918: link a6,IMM (0) ! 2919: |#PROLOGUE# 1 ! 2920: movl a6@(12),sp@- ! 2921: movl a6@(8),sp@- ! 2922: jbsr SYM (__cmpsf2) ! 2923: |#PROLOGUE# 2 ! 2924: unlk a6 ! 2925: |#PROLOGUE# 3 ! 2926: rts ! 2927: #endif /* L_gesf2 */ ! 2928: ! 2929: #ifdef L_ltsf2 ! 2930: .text ! 2931: .proc ! 2932: |#PROC# 04 ! 2933: LF18 = 8 ! 2934: LS18 = 132 ! 2935: LFF18 = 0 ! 2936: LSS18 = 0 ! 2937: LV18 = 0 ! 2938: .text ! 2939: .globl SYM (__ltsf2) ! 2940: SYM (__ltsf2): ! 2941: |#PROLOGUE# 0 ! 2942: link a6,IMM (0) ! 2943: |#PROLOGUE# 1 ! 2944: movl a6@(12),sp@- ! 2945: movl a6@(8),sp@- ! 2946: jbsr SYM (__cmpsf2) ! 2947: |#PROLOGUE# 2 ! 2948: unlk a6 ! 2949: |#PROLOGUE# 3 ! 2950: rts ! 2951: #endif /* L_ltsf2 */ ! 2952: ! 2953: #ifdef L_lesf2 ! 2954: .text ! 2955: .proc ! 2956: |#PROC# 04 ! 2957: LF18 = 8 ! 2958: LS18 = 132 ! 2959: LFF18 = 0 ! 2960: LSS18 = 0 ! 2961: LV18 = 0 ! 2962: .text ! 2963: .globl SYM (__lesf2) ! 2964: SYM (__lesf2): ! 2965: |#PROLOGUE# 0 ! 2966: link a6,IMM (0) ! 2967: |#PROLOGUE# 1 ! 2968: movl a6@(12),sp@- ! 2969: movl a6@(8),sp@- ! 2970: jbsr SYM (__cmpsf2) ! 2971: |#PROLOGUE# 2 ! 2972: unlk a6 ! 2973: |#PROLOGUE# 3 ! 2974: rts ! 2975: #endif /* L_lesf2 */ ! 2976:
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.