lb1sf68.S revision 1.1.1.6
1 /* libgcc routines for 68000 w/o floating-point hardware. 2 Copyright (C) 1994-2017 Free Software Foundation, Inc. 3 4 This file is part of GCC. 5 6 GCC is free software; you can redistribute it and/or modify it 7 under the terms of the GNU General Public License as published by the 8 Free Software Foundation; either version 3, or (at your option) any 9 later version. 10 11 This file is distributed in the hope that it will be useful, but 12 WITHOUT ANY WARRANTY; without even the implied warranty of 13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 14 General Public License for more details. 15 16 Under Section 7 of GPL version 3, you are granted additional 17 permissions described in the GCC Runtime Library Exception, version 18 3.1, as published by the Free Software Foundation. 19 20 You should have received a copy of the GNU General Public License and 21 a copy of the GCC Runtime Library Exception along with this program; 22 see the files COPYING3 and COPYING.RUNTIME respectively. If not, see 23 <http://www.gnu.org/licenses/>. */ 24 25 /* Use this one for any 680x0; assumes no floating point hardware. 26 The trailing " '" appearing on some lines is for ANSI preprocessors. Yuk. 27 Some of this code comes from MINIX, via the folks at ericsson. 28 D. V. Henkel-Wallace (gumby (at) cygnus.com) Fete Bastille, 1992 29 */ 30 31 /* These are predefined by new versions of GNU cpp. */ 32 33 #ifndef __USER_LABEL_PREFIX__ 34 #define __USER_LABEL_PREFIX__ _ 35 #endif 36 37 #ifndef __REGISTER_PREFIX__ 38 #define __REGISTER_PREFIX__ 39 #endif 40 41 #ifndef __IMMEDIATE_PREFIX__ 42 #define __IMMEDIATE_PREFIX__ # 43 #endif 44 45 /* ANSI concatenation macros. */ 46 47 #define CONCAT1(a, b) CONCAT2(a, b) 48 #define CONCAT2(a, b) a ## b 49 50 /* Use the right prefix for global labels. */ 51 52 #define SYM(x) CONCAT1 (__USER_LABEL_PREFIX__, x) 53 54 /* Note that X is a function. */ 55 56 #ifdef __ELF__ 57 #define FUNC(x) .type SYM(x),function 58 #else 59 /* The .proc pseudo-op is accepted, but ignored, by GAS. We could just 60 define this to the empty string for non-ELF systems, but defining it 61 to .proc means that the information is available to the assembler if 62 the need arises. */ 63 #define FUNC(x) .proc 64 #endif 65 66 /* Use the right prefix for registers. */ 67 68 #define REG(x) CONCAT1 (__REGISTER_PREFIX__, x) 69 70 /* Use the right prefix for immediate values. */ 71 72 #define IMM(x) CONCAT1 (__IMMEDIATE_PREFIX__, x) 73 74 #define d0 REG (d0) 75 #define d1 REG (d1) 76 #define d2 REG (d2) 77 #define d3 REG (d3) 78 #define d4 REG (d4) 79 #define d5 REG (d5) 80 #define d6 REG (d6) 81 #define d7 REG (d7) 82 #define a0 REG (a0) 83 #define a1 REG (a1) 84 #define a2 REG (a2) 85 #define a3 REG (a3) 86 #define a4 REG (a4) 87 #define a5 REG (a5) 88 #define a6 REG (a6) 89 #define fp REG (fp) 90 #define sp REG (sp) 91 #define pc REG (pc) 92 93 /* Provide a few macros to allow for PIC code support. 94 * With PIC, data is stored A5 relative so we've got to take a bit of special 95 * care to ensure that all loads of global data is via A5. PIC also requires 96 * jumps and subroutine calls to be PC relative rather than absolute. We cheat 97 * a little on this and in the PIC case, we use short offset branches and 98 * hope that the final object code is within range (which it should be). 99 */ 100 #ifndef __PIC__ 101 102 /* Non PIC (absolute/relocatable) versions */ 103 104 .macro PICCALL addr 105 jbsr \addr 106 .endm 107 108 .macro PICJUMP addr 109 jmp \addr 110 .endm 111 112 .macro PICLEA sym, reg 113 lea \sym, \reg 114 .endm 115 116 .macro PICPEA sym, areg 117 pea \sym 118 .endm 119 120 #else /* __PIC__ */ 121 122 # if defined (__uClinux__) 123 124 /* Versions for uClinux */ 125 126 # if defined(__ID_SHARED_LIBRARY__) 127 128 /* -mid-shared-library versions */ 129 130 .macro PICLEA sym, reg 131 movel a5@(_current_shared_library_a5_offset_), \reg 132 movel \sym@GOT(\reg), \reg 133 .endm 134 135 .macro PICPEA sym, areg 136 movel a5@(_current_shared_library_a5_offset_), \areg 137 movel \sym@GOT(\areg), sp@- 138 .endm 139 140 .macro PICCALL addr 141 PICLEA \addr,a0 142 jsr a0@ 143 .endm 144 145 .macro PICJUMP addr 146 PICLEA \addr,a0 147 jmp a0@ 148 .endm 149 150 # else /* !__ID_SHARED_LIBRARY__ */ 151 152 /* Versions for -msep-data */ 153 154 .macro PICLEA sym, reg 155 movel \sym@GOT(a5), \reg 156 .endm 157 158 .macro PICPEA sym, areg 159 movel \sym@GOT(a5), sp@- 160 .endm 161 162 .macro PICCALL addr 163 #if defined (__mcoldfire__) && !defined (__mcfisab__) && !defined (__mcfisac__) 164 lea \addr-.-8,a0 165 jsr pc@(a0) 166 #else 167 jbsr \addr 168 #endif 169 .endm 170 171 .macro PICJUMP addr 172 /* ISA C has no bra.l instruction, and since this assembly file 173 gets assembled into multiple object files, we avoid the 174 bra instruction entirely. */ 175 #if defined (__mcoldfire__) && !defined (__mcfisab__) 176 lea \addr-.-8,a0 177 jmp pc@(a0) 178 #else 179 bra \addr 180 #endif 181 .endm 182 183 # endif 184 185 # else /* !__uClinux__ */ 186 187 /* Versions for Linux */ 188 189 .macro PICLEA sym, reg 190 movel #_GLOBAL_OFFSET_TABLE_@GOTPC, \reg 191 lea (-6, pc, \reg), \reg 192 movel \sym@GOT(\reg), \reg 193 .endm 194 195 .macro PICPEA sym, areg 196 movel #_GLOBAL_OFFSET_TABLE_@GOTPC, \areg 197 lea (-6, pc, \areg), \areg 198 movel \sym@GOT(\areg), sp@- 199 .endm 200 201 .macro PICCALL addr 202 #if defined (__mcoldfire__) && !defined (__mcfisab__) && !defined (__mcfisac__) 203 lea \addr-.-8,a0 204 jsr pc@(a0) 205 #else 206 jbsr \addr 207 #endif 208 .endm 209 210 .macro PICJUMP addr 211 /* ISA C has no bra.l instruction, and since this assembly file 212 gets assembled into multiple object files, we avoid the 213 bra instruction entirely. */ 214 #if defined (__mcoldfire__) && !defined (__mcfisab__) 215 lea \addr-.-8,a0 216 jmp pc@(a0) 217 #else 218 bra \addr 219 #endif 220 .endm 221 222 # endif 223 #endif /* __PIC__ */ 224 225 226 #ifdef L_floatex 227 228 | This is an attempt at a decent floating point (single, double and 229 | extended double) code for the GNU C compiler. It should be easy to 230 | adapt to other compilers (but beware of the local labels!). 231 232 | Starting date: 21 October, 1990 233 234 | It is convenient to introduce the notation (s,e,f) for a floating point 235 | number, where s=sign, e=exponent, f=fraction. We will call a floating 236 | point number fpn to abbreviate, independently of the precision. 237 | Let MAX_EXP be in each case the maximum exponent (255 for floats, 1023 238 | for doubles and 16383 for long doubles). We then have the following 239 | different cases: 240 | 1. Normalized fpns have 0 < e < MAX_EXP. They correspond to 241 | (-1)^s x 1.f x 2^(e-bias-1). 242 | 2. Denormalized fpns have e=0. They correspond to numbers of the form 243 | (-1)^s x 0.f x 2^(-bias). 244 | 3. +/-INFINITY have e=MAX_EXP, f=0. 245 | 4. Quiet NaN (Not a Number) have all bits set. 246 | 5. Signaling NaN (Not a Number) have s=0, e=MAX_EXP, f=1. 247 248 |============================================================================= 249 | exceptions 250 |============================================================================= 251 252 | This is the floating point condition code register (_fpCCR): 253 | 254 | struct { 255 | short _exception_bits; 256 | short _trap_enable_bits; 257 | short _sticky_bits; 258 | short _rounding_mode; 259 | short _format; 260 | short _last_operation; 261 | union { 262 | float sf; 263 | double df; 264 | } _operand1; 265 | union { 266 | float sf; 267 | double df; 268 | } _operand2; 269 | } _fpCCR; 270 271 .data 272 .even 273 274 .globl SYM (_fpCCR) 275 276 SYM (_fpCCR): 277 __exception_bits: 278 .word 0 279 __trap_enable_bits: 280 .word 0 281 __sticky_bits: 282 .word 0 283 __rounding_mode: 284 .word ROUND_TO_NEAREST 285 __format: 286 .word NIL 287 __last_operation: 288 .word NOOP 289 __operand1: 290 .long 0 291 .long 0 292 __operand2: 293 .long 0 294 .long 0 295 296 | Offsets: 297 EBITS = __exception_bits - SYM (_fpCCR) 298 TRAPE = __trap_enable_bits - SYM (_fpCCR) 299 STICK = __sticky_bits - SYM (_fpCCR) 300 ROUND = __rounding_mode - SYM (_fpCCR) 301 FORMT = __format - SYM (_fpCCR) 302 LASTO = __last_operation - SYM (_fpCCR) 303 OPER1 = __operand1 - SYM (_fpCCR) 304 OPER2 = __operand2 - SYM (_fpCCR) 305 306 | The following exception types are supported: 307 INEXACT_RESULT = 0x0001 308 UNDERFLOW = 0x0002 309 OVERFLOW = 0x0004 310 DIVIDE_BY_ZERO = 0x0008 311 INVALID_OPERATION = 0x0010 312 313 | The allowed rounding modes are: 314 UNKNOWN = -1 315 ROUND_TO_NEAREST = 0 | round result to nearest representable value 316 ROUND_TO_ZERO = 1 | round result towards zero 317 ROUND_TO_PLUS = 2 | round result towards plus infinity 318 ROUND_TO_MINUS = 3 | round result towards minus infinity 319 320 | The allowed values of format are: 321 NIL = 0 322 SINGLE_FLOAT = 1 323 DOUBLE_FLOAT = 2 324 LONG_FLOAT = 3 325 326 | The allowed values for the last operation are: 327 NOOP = 0 328 ADD = 1 329 MULTIPLY = 2 330 DIVIDE = 3 331 NEGATE = 4 332 COMPARE = 5 333 EXTENDSFDF = 6 334 TRUNCDFSF = 7 335 336 |============================================================================= 337 | __clear_sticky_bits 338 |============================================================================= 339 340 | The sticky bits are normally not cleared (thus the name), whereas the 341 | exception type and exception value reflect the last computation. 342 | This routine is provided to clear them (you can also write to _fpCCR, 343 | since it is globally visible). 344 345 .globl SYM (__clear_sticky_bit) 346 347 .text 348 .even 349 350 | void __clear_sticky_bits(void); 351 SYM (__clear_sticky_bit): 352 PICLEA SYM (_fpCCR),a0 353 #ifndef __mcoldfire__ 354 movew IMM (0),a0@(STICK) 355 #else 356 clr.w a0@(STICK) 357 #endif 358 rts 359 360 |============================================================================= 361 | $_exception_handler 362 |============================================================================= 363 364 .globl $_exception_handler 365 366 .text 367 .even 368 369 | This is the common exit point if an exception occurs. 370 | NOTE: it is NOT callable from C! 371 | It expects the exception type in d7, the format (SINGLE_FLOAT, 372 | DOUBLE_FLOAT or LONG_FLOAT) in d6, and the last operation code in d5. 373 | It sets the corresponding exception and sticky bits, and the format. 374 | Depending on the format if fills the corresponding slots for the 375 | operands which produced the exception (all this information is provided 376 | so if you write your own exception handlers you have enough information 377 | to deal with the problem). 378 | Then checks to see if the corresponding exception is trap-enabled, 379 | in which case it pushes the address of _fpCCR and traps through 380 | trap FPTRAP (15 for the moment). 381 382 FPTRAP = 15 383 384 $_exception_handler: 385 PICLEA SYM (_fpCCR),a0 386 movew d7,a0@(EBITS) | set __exception_bits 387 #ifndef __mcoldfire__ 388 orw d7,a0@(STICK) | and __sticky_bits 389 #else 390 movew a0@(STICK),d4 391 orl d7,d4 392 movew d4,a0@(STICK) 393 #endif 394 movew d6,a0@(FORMT) | and __format 395 movew d5,a0@(LASTO) | and __last_operation 396 397 | Now put the operands in place: 398 #ifndef __mcoldfire__ 399 cmpw IMM (SINGLE_FLOAT),d6 400 #else 401 cmpl IMM (SINGLE_FLOAT),d6 402 #endif 403 beq 1f 404 movel a6@(8),a0@(OPER1) 405 movel a6@(12),a0@(OPER1+4) 406 movel a6@(16),a0@(OPER2) 407 movel a6@(20),a0@(OPER2+4) 408 bra 2f 409 1: movel a6@(8),a0@(OPER1) 410 movel a6@(12),a0@(OPER2) 411 2: 412 | And check whether the exception is trap-enabled: 413 #ifndef __mcoldfire__ 414 andw a0@(TRAPE),d7 | is exception trap-enabled? 415 #else 416 clrl d6 417 movew a0@(TRAPE),d6 418 andl d6,d7 419 #endif 420 beq 1f | no, exit 421 PICPEA SYM (_fpCCR),a1 | yes, push address of _fpCCR 422 trap IMM (FPTRAP) | and trap 423 #ifndef __mcoldfire__ 424 1: moveml sp@+,d2-d7 | restore data registers 425 #else 426 1: moveml sp@,d2-d7 427 | XXX if frame pointer is ever removed, stack pointer must 428 | be adjusted here. 429 #endif 430 unlk a6 | and return 431 rts 432 #endif /* L_floatex */ 433 434 #ifdef L_mulsi3 435 .text 436 FUNC(__mulsi3) 437 .globl SYM (__mulsi3) 438 SYM (__mulsi3): 439 movew sp@(4), d0 /* x0 -> d0 */ 440 muluw sp@(10), d0 /* x0*y1 */ 441 movew sp@(6), d1 /* x1 -> d1 */ 442 muluw sp@(8), d1 /* x1*y0 */ 443 #ifndef __mcoldfire__ 444 addw d1, d0 445 #else 446 addl d1, d0 447 #endif 448 swap d0 449 clrw d0 450 movew sp@(6), d1 /* x1 -> d1 */ 451 muluw sp@(10), d1 /* x1*y1 */ 452 addl d1, d0 453 454 rts 455 #endif /* L_mulsi3 */ 456 457 #ifdef L_udivsi3 458 .text 459 FUNC(__udivsi3) 460 .globl SYM (__udivsi3) 461 SYM (__udivsi3): 462 #ifndef __mcoldfire__ 463 movel d2, sp@- 464 movel sp@(12), d1 /* d1 = divisor */ 465 movel sp@(8), d0 /* d0 = dividend */ 466 467 cmpl IMM (0x10000), d1 /* divisor >= 2 ^ 16 ? */ 468 jcc L3 /* then try next algorithm */ 469 movel d0, d2 470 clrw d2 471 swap d2 472 divu d1, d2 /* high quotient in lower word */ 473 movew d2, d0 /* save high quotient */ 474 swap d0 475 movew sp@(10), d2 /* get low dividend + high rest */ 476 divu d1, d2 /* low quotient */ 477 movew d2, d0 478 jra L6 479 480 L3: movel d1, d2 /* use d2 as divisor backup */ 481 L4: lsrl IMM (1), d1 /* shift divisor */ 482 lsrl IMM (1), d0 /* shift dividend */ 483 cmpl IMM (0x10000), d1 /* still divisor >= 2 ^ 16 ? */ 484 jcc L4 485 divu d1, d0 /* now we have 16-bit divisor */ 486 andl IMM (0xffff), d0 /* mask out divisor, ignore remainder */ 487 488 /* Multiply the 16-bit tentative quotient with the 32-bit divisor. Because of 489 the operand ranges, this might give a 33-bit product. If this product is 490 greater than the dividend, the tentative quotient was too large. */ 491 movel d2, d1 492 mulu d0, d1 /* low part, 32 bits */ 493 swap d2 494 mulu d0, d2 /* high part, at most 17 bits */ 495 swap d2 /* align high part with low part */ 496 tstw d2 /* high part 17 bits? */ 497 jne L5 /* if 17 bits, quotient was too large */ 498 addl d2, d1 /* add parts */ 499 jcs L5 /* if sum is 33 bits, quotient was too large */ 500 cmpl sp@(8), d1 /* compare the sum with the dividend */ 501 jls L6 /* if sum > dividend, quotient was too large */ 502 L5: subql IMM (1), d0 /* adjust quotient */ 503 504 L6: movel sp@+, d2 505 rts 506 507 #else /* __mcoldfire__ */ 508 509 /* ColdFire implementation of non-restoring division algorithm from 510 Hennessy & Patterson, Appendix A. */ 511 link a6,IMM (-12) 512 moveml d2-d4,sp@ 513 movel a6@(8),d0 514 movel a6@(12),d1 515 clrl d2 | clear p 516 moveq IMM (31),d4 517 L1: addl d0,d0 | shift reg pair (p,a) one bit left 518 addxl d2,d2 519 movl d2,d3 | subtract b from p, store in tmp. 520 subl d1,d3 521 jcs L2 | if no carry, 522 bset IMM (0),d0 | set the low order bit of a to 1, 523 movl d3,d2 | and store tmp in p. 524 L2: subql IMM (1),d4 525 jcc L1 526 moveml sp@,d2-d4 | restore data registers 527 unlk a6 | and return 528 rts 529 #endif /* __mcoldfire__ */ 530 531 #endif /* L_udivsi3 */ 532 533 #ifdef L_divsi3 534 .text 535 FUNC(__divsi3) 536 .globl SYM (__divsi3) 537 SYM (__divsi3): 538 movel d2, sp@- 539 540 moveq IMM (1), d2 /* sign of result stored in d2 (=1 or =-1) */ 541 movel sp@(12), d1 /* d1 = divisor */ 542 jpl L1 543 negl d1 544 #ifndef __mcoldfire__ 545 negb d2 /* change sign because divisor <0 */ 546 #else 547 negl d2 /* change sign because divisor <0 */ 548 #endif 549 L1: movel sp@(8), d0 /* d0 = dividend */ 550 jpl L2 551 negl d0 552 #ifndef __mcoldfire__ 553 negb d2 554 #else 555 negl d2 556 #endif 557 558 L2: movel d1, sp@- 559 movel d0, sp@- 560 PICCALL SYM (__udivsi3) /* divide abs(dividend) by abs(divisor) */ 561 addql IMM (8), sp 562 563 tstb d2 564 jpl L3 565 negl d0 566 567 L3: movel sp@+, d2 568 rts 569 #endif /* L_divsi3 */ 570 571 #ifdef L_umodsi3 572 .text 573 FUNC(__umodsi3) 574 .globl SYM (__umodsi3) 575 SYM (__umodsi3): 576 movel sp@(8), d1 /* d1 = divisor */ 577 movel sp@(4), d0 /* d0 = dividend */ 578 movel d1, sp@- 579 movel d0, sp@- 580 PICCALL SYM (__udivsi3) 581 addql IMM (8), sp 582 movel sp@(8), d1 /* d1 = divisor */ 583 #ifndef __mcoldfire__ 584 movel d1, sp@- 585 movel d0, sp@- 586 PICCALL SYM (__mulsi3) /* d0 = (a/b)*b */ 587 addql IMM (8), sp 588 #else 589 mulsl d1,d0 590 #endif 591 movel sp@(4), d1 /* d1 = dividend */ 592 subl d0, d1 /* d1 = a - (a/b)*b */ 593 movel d1, d0 594 rts 595 #endif /* L_umodsi3 */ 596 597 #ifdef L_modsi3 598 .text 599 FUNC(__modsi3) 600 .globl SYM (__modsi3) 601 SYM (__modsi3): 602 movel sp@(8), d1 /* d1 = divisor */ 603 movel sp@(4), d0 /* d0 = dividend */ 604 movel d1, sp@- 605 movel d0, sp@- 606 PICCALL SYM (__divsi3) 607 addql IMM (8), sp 608 movel sp@(8), d1 /* d1 = divisor */ 609 #ifndef __mcoldfire__ 610 movel d1, sp@- 611 movel d0, sp@- 612 PICCALL SYM (__mulsi3) /* d0 = (a/b)*b */ 613 addql IMM (8), sp 614 #else 615 mulsl d1,d0 616 #endif 617 movel sp@(4), d1 /* d1 = dividend */ 618 subl d0, d1 /* d1 = a - (a/b)*b */ 619 movel d1, d0 620 rts 621 #endif /* L_modsi3 */ 622 623 624 #ifdef L_double 625 626 .globl SYM (_fpCCR) 627 .globl $_exception_handler 628 629 QUIET_NaN = 0xffffffff 630 631 D_MAX_EXP = 0x07ff 632 D_BIAS = 1022 633 DBL_MAX_EXP = D_MAX_EXP - D_BIAS 634 DBL_MIN_EXP = 1 - D_BIAS 635 DBL_MANT_DIG = 53 636 637 INEXACT_RESULT = 0x0001 638 UNDERFLOW = 0x0002 639 OVERFLOW = 0x0004 640 DIVIDE_BY_ZERO = 0x0008 641 INVALID_OPERATION = 0x0010 642 643 DOUBLE_FLOAT = 2 644 645 NOOP = 0 646 ADD = 1 647 MULTIPLY = 2 648 DIVIDE = 3 649 NEGATE = 4 650 COMPARE = 5 651 EXTENDSFDF = 6 652 TRUNCDFSF = 7 653 654 UNKNOWN = -1 655 ROUND_TO_NEAREST = 0 | round result to nearest representable value 656 ROUND_TO_ZERO = 1 | round result towards zero 657 ROUND_TO_PLUS = 2 | round result towards plus infinity 658 ROUND_TO_MINUS = 3 | round result towards minus infinity 659 660 | Entry points: 661 662 .globl SYM (__adddf3) 663 .globl SYM (__subdf3) 664 .globl SYM (__muldf3) 665 .globl SYM (__divdf3) 666 .globl SYM (__negdf2) 667 .globl SYM (__cmpdf2) 668 .globl SYM (__cmpdf2_internal) 669 .hidden SYM (__cmpdf2_internal) 670 671 .text 672 .even 673 674 | These are common routines to return and signal exceptions. 675 676 Ld$den: 677 | Return and signal a denormalized number 678 orl d7,d0 679 movew IMM (INEXACT_RESULT+UNDERFLOW),d7 680 moveq IMM (DOUBLE_FLOAT),d6 681 PICJUMP $_exception_handler 682 683 Ld$infty: 684 Ld$overflow: 685 | Return a properly signed INFINITY and set the exception flags 686 movel IMM (0x7ff00000),d0 687 movel IMM (0),d1 688 orl d7,d0 689 movew IMM (INEXACT_RESULT+OVERFLOW),d7 690 moveq IMM (DOUBLE_FLOAT),d6 691 PICJUMP $_exception_handler 692 693 Ld$underflow: 694 | Return 0 and set the exception flags 695 movel IMM (0),d0 696 movel d0,d1 697 movew IMM (INEXACT_RESULT+UNDERFLOW),d7 698 moveq IMM (DOUBLE_FLOAT),d6 699 PICJUMP $_exception_handler 700 701 Ld$inop: 702 | Return a quiet NaN and set the exception flags 703 movel IMM (QUIET_NaN),d0 704 movel d0,d1 705 movew IMM (INEXACT_RESULT+INVALID_OPERATION),d7 706 moveq IMM (DOUBLE_FLOAT),d6 707 PICJUMP $_exception_handler 708 709 Ld$div$0: 710 | Return a properly signed INFINITY and set the exception flags 711 movel IMM (0x7ff00000),d0 712 movel IMM (0),d1 713 orl d7,d0 714 movew IMM (INEXACT_RESULT+DIVIDE_BY_ZERO),d7 715 moveq IMM (DOUBLE_FLOAT),d6 716 PICJUMP $_exception_handler 717 718 |============================================================================= 719 |============================================================================= 720 | double precision routines 721 |============================================================================= 722 |============================================================================= 723 724 | A double precision floating point number (double) has the format: 725 | 726 | struct _double { 727 | unsigned int sign : 1; /* sign bit */ 728 | unsigned int exponent : 11; /* exponent, shifted by 126 */ 729 | unsigned int fraction : 52; /* fraction */ 730 | } double; 731 | 732 | Thus sizeof(double) = 8 (64 bits). 733 | 734 | All the routines are callable from C programs, and return the result 735 | in the register pair d0-d1. They also preserve all registers except 736 | d0-d1 and a0-a1. 737 738 |============================================================================= 739 | __subdf3 740 |============================================================================= 741 742 | double __subdf3(double, double); 743 FUNC(__subdf3) 744 SYM (__subdf3): 745 bchg IMM (31),sp@(12) | change sign of second operand 746 | and fall through, so we always add 747 |============================================================================= 748 | __adddf3 749 |============================================================================= 750 751 | double __adddf3(double, double); 752 FUNC(__adddf3) 753 SYM (__adddf3): 754 #ifndef __mcoldfire__ 755 link a6,IMM (0) | everything will be done in registers 756 moveml d2-d7,sp@- | save all data registers and a2 (but d0-d1) 757 #else 758 link a6,IMM (-24) 759 moveml d2-d7,sp@ 760 #endif 761 movel a6@(8),d0 | get first operand 762 movel a6@(12),d1 | 763 movel a6@(16),d2 | get second operand 764 movel a6@(20),d3 | 765 766 movel d0,d7 | get d0's sign bit in d7 ' 767 addl d1,d1 | check and clear sign bit of a, and gain one 768 addxl d0,d0 | bit of extra precision 769 beq Ladddf$b | if zero return second operand 770 771 movel d2,d6 | save sign in d6 772 addl d3,d3 | get rid of sign bit and gain one bit of 773 addxl d2,d2 | extra precision 774 beq Ladddf$a | if zero return first operand 775 776 andl IMM (0x80000000),d7 | isolate a's sign bit ' 777 swap d6 | and also b's sign bit ' 778 #ifndef __mcoldfire__ 779 andw IMM (0x8000),d6 | 780 orw d6,d7 | and combine them into d7, so that a's sign ' 781 | bit is in the high word and b's is in the ' 782 | low word, so d6 is free to be used 783 #else 784 andl IMM (0x8000),d6 785 orl d6,d7 786 #endif 787 movel d7,a0 | now save d7 into a0, so d7 is free to 788 | be used also 789 790 | Get the exponents and check for denormalized and/or infinity. 791 792 movel IMM (0x001fffff),d6 | mask for the fraction 793 movel IMM (0x00200000),d7 | mask to put hidden bit back 794 795 movel d0,d4 | 796 andl d6,d0 | get fraction in d0 797 notl d6 | make d6 into mask for the exponent 798 andl d6,d4 | get exponent in d4 799 beq Ladddf$a$den | branch if a is denormalized 800 cmpl d6,d4 | check for INFINITY or NaN 801 beq Ladddf$nf | 802 orl d7,d0 | and put hidden bit back 803 Ladddf$1: 804 swap d4 | shift right exponent so that it starts 805 #ifndef __mcoldfire__ 806 lsrw IMM (5),d4 | in bit 0 and not bit 20 807 #else 808 lsrl IMM (5),d4 | in bit 0 and not bit 20 809 #endif 810 | Now we have a's exponent in d4 and fraction in d0-d1 ' 811 movel d2,d5 | save b to get exponent 812 andl d6,d5 | get exponent in d5 813 beq Ladddf$b$den | branch if b is denormalized 814 cmpl d6,d5 | check for INFINITY or NaN 815 beq Ladddf$nf 816 notl d6 | make d6 into mask for the fraction again 817 andl d6,d2 | and get fraction in d2 818 orl d7,d2 | and put hidden bit back 819 Ladddf$2: 820 swap d5 | shift right exponent so that it starts 821 #ifndef __mcoldfire__ 822 lsrw IMM (5),d5 | in bit 0 and not bit 20 823 #else 824 lsrl IMM (5),d5 | in bit 0 and not bit 20 825 #endif 826 827 | Now we have b's exponent in d5 and fraction in d2-d3. ' 828 829 | The situation now is as follows: the signs are combined in a0, the 830 | numbers are in d0-d1 (a) and d2-d3 (b), and the exponents in d4 (a) 831 | and d5 (b). To do the rounding correctly we need to keep all the 832 | bits until the end, so we need to use d0-d1-d2-d3 for the first number 833 | and d4-d5-d6-d7 for the second. To do this we store (temporarily) the 834 | exponents in a2-a3. 835 836 #ifndef __mcoldfire__ 837 moveml a2-a3,sp@- | save the address registers 838 #else 839 movel a2,sp@- 840 movel a3,sp@- 841 movel a4,sp@- 842 #endif 843 844 movel d4,a2 | save the exponents 845 movel d5,a3 | 846 847 movel IMM (0),d7 | and move the numbers around 848 movel d7,d6 | 849 movel d3,d5 | 850 movel d2,d4 | 851 movel d7,d3 | 852 movel d7,d2 | 853 854 | Here we shift the numbers until the exponents are the same, and put 855 | the largest exponent in a2. 856 #ifndef __mcoldfire__ 857 exg d4,a2 | get exponents back 858 exg d5,a3 | 859 cmpw d4,d5 | compare the exponents 860 #else 861 movel d4,a4 | get exponents back 862 movel a2,d4 863 movel a4,a2 864 movel d5,a4 865 movel a3,d5 866 movel a4,a3 867 cmpl d4,d5 | compare the exponents 868 #endif 869 beq Ladddf$3 | if equal don't shift ' 870 bhi 9f | branch if second exponent is higher 871 872 | Here we have a's exponent larger than b's, so we have to shift b. We do 873 | this by using as counter d2: 874 1: movew d4,d2 | move largest exponent to d2 875 #ifndef __mcoldfire__ 876 subw d5,d2 | and subtract second exponent 877 exg d4,a2 | get back the longs we saved 878 exg d5,a3 | 879 #else 880 subl d5,d2 | and subtract second exponent 881 movel d4,a4 | get back the longs we saved 882 movel a2,d4 883 movel a4,a2 884 movel d5,a4 885 movel a3,d5 886 movel a4,a3 887 #endif 888 | if difference is too large we don't shift (actually, we can just exit) ' 889 #ifndef __mcoldfire__ 890 cmpw IMM (DBL_MANT_DIG+2),d2 891 #else 892 cmpl IMM (DBL_MANT_DIG+2),d2 893 #endif 894 bge Ladddf$b$small 895 #ifndef __mcoldfire__ 896 cmpw IMM (32),d2 | if difference >= 32, shift by longs 897 #else 898 cmpl IMM (32),d2 | if difference >= 32, shift by longs 899 #endif 900 bge 5f 901 2: 902 #ifndef __mcoldfire__ 903 cmpw IMM (16),d2 | if difference >= 16, shift by words 904 #else 905 cmpl IMM (16),d2 | if difference >= 16, shift by words 906 #endif 907 bge 6f 908 bra 3f | enter dbra loop 909 910 4: 911 #ifndef __mcoldfire__ 912 lsrl IMM (1),d4 913 roxrl IMM (1),d5 914 roxrl IMM (1),d6 915 roxrl IMM (1),d7 916 #else 917 lsrl IMM (1),d7 918 btst IMM (0),d6 919 beq 10f 920 bset IMM (31),d7 921 10: lsrl IMM (1),d6 922 btst IMM (0),d5 923 beq 11f 924 bset IMM (31),d6 925 11: lsrl IMM (1),d5 926 btst IMM (0),d4 927 beq 12f 928 bset IMM (31),d5 929 12: lsrl IMM (1),d4 930 #endif 931 3: 932 #ifndef __mcoldfire__ 933 dbra d2,4b 934 #else 935 subql IMM (1),d2 936 bpl 4b 937 #endif 938 movel IMM (0),d2 939 movel d2,d3 940 bra Ladddf$4 941 5: 942 movel d6,d7 943 movel d5,d6 944 movel d4,d5 945 movel IMM (0),d4 946 #ifndef __mcoldfire__ 947 subw IMM (32),d2 948 #else 949 subl IMM (32),d2 950 #endif 951 bra 2b 952 6: 953 movew d6,d7 954 swap d7 955 movew d5,d6 956 swap d6 957 movew d4,d5 958 swap d5 959 movew IMM (0),d4 960 swap d4 961 #ifndef __mcoldfire__ 962 subw IMM (16),d2 963 #else 964 subl IMM (16),d2 965 #endif 966 bra 3b 967 968 9: 969 #ifndef __mcoldfire__ 970 exg d4,d5 971 movew d4,d6 972 subw d5,d6 | keep d5 (largest exponent) in d4 973 exg d4,a2 974 exg d5,a3 975 #else 976 movel d5,d6 977 movel d4,d5 978 movel d6,d4 979 subl d5,d6 980 movel d4,a4 981 movel a2,d4 982 movel a4,a2 983 movel d5,a4 984 movel a3,d5 985 movel a4,a3 986 #endif 987 | if difference is too large we don't shift (actually, we can just exit) ' 988 #ifndef __mcoldfire__ 989 cmpw IMM (DBL_MANT_DIG+2),d6 990 #else 991 cmpl IMM (DBL_MANT_DIG+2),d6 992 #endif 993 bge Ladddf$a$small 994 #ifndef __mcoldfire__ 995 cmpw IMM (32),d6 | if difference >= 32, shift by longs 996 #else 997 cmpl IMM (32),d6 | if difference >= 32, shift by longs 998 #endif 999 bge 5f 1000 2: 1001 #ifndef __mcoldfire__ 1002 cmpw IMM (16),d6 | if difference >= 16, shift by words 1003 #else 1004 cmpl IMM (16),d6 | if difference >= 16, shift by words 1005 #endif 1006 bge 6f 1007 bra 3f | enter dbra loop 1008 1009 4: 1010 #ifndef __mcoldfire__ 1011 lsrl IMM (1),d0 1012 roxrl IMM (1),d1 1013 roxrl IMM (1),d2 1014 roxrl IMM (1),d3 1015 #else 1016 lsrl IMM (1),d3 1017 btst IMM (0),d2 1018 beq 10f 1019 bset IMM (31),d3 1020 10: lsrl IMM (1),d2 1021 btst IMM (0),d1 1022 beq 11f 1023 bset IMM (31),d2 1024 11: lsrl IMM (1),d1 1025 btst IMM (0),d0 1026 beq 12f 1027 bset IMM (31),d1 1028 12: lsrl IMM (1),d0 1029 #endif 1030 3: 1031 #ifndef __mcoldfire__ 1032 dbra d6,4b 1033 #else 1034 subql IMM (1),d6 1035 bpl 4b 1036 #endif 1037 movel IMM (0),d7 1038 movel d7,d6 1039 bra Ladddf$4 1040 5: 1041 movel d2,d3 1042 movel d1,d2 1043 movel d0,d1 1044 movel IMM (0),d0 1045 #ifndef __mcoldfire__ 1046 subw IMM (32),d6 1047 #else 1048 subl IMM (32),d6 1049 #endif 1050 bra 2b 1051 6: 1052 movew d2,d3 1053 swap d3 1054 movew d1,d2 1055 swap d2 1056 movew d0,d1 1057 swap d1 1058 movew IMM (0),d0 1059 swap d0 1060 #ifndef __mcoldfire__ 1061 subw IMM (16),d6 1062 #else 1063 subl IMM (16),d6 1064 #endif 1065 bra 3b 1066 Ladddf$3: 1067 #ifndef __mcoldfire__ 1068 exg d4,a2 1069 exg d5,a3 1070 #else 1071 movel d4,a4 1072 movel a2,d4 1073 movel a4,a2 1074 movel d5,a4 1075 movel a3,d5 1076 movel a4,a3 1077 #endif 1078 Ladddf$4: 1079 | Now we have the numbers in d0--d3 and d4--d7, the exponent in a2, and 1080 | the signs in a4. 1081 1082 | Here we have to decide whether to add or subtract the numbers: 1083 #ifndef __mcoldfire__ 1084 exg d7,a0 | get the signs 1085 exg d6,a3 | a3 is free to be used 1086 #else 1087 movel d7,a4 1088 movel a0,d7 1089 movel a4,a0 1090 movel d6,a4 1091 movel a3,d6 1092 movel a4,a3 1093 #endif 1094 movel d7,d6 | 1095 movew IMM (0),d7 | get a's sign in d7 ' 1096 swap d6 | 1097 movew IMM (0),d6 | and b's sign in d6 ' 1098 eorl d7,d6 | compare the signs 1099 bmi Lsubdf$0 | if the signs are different we have 1100 | to subtract 1101 #ifndef __mcoldfire__ 1102 exg d7,a0 | else we add the numbers 1103 exg d6,a3 | 1104 #else 1105 movel d7,a4 1106 movel a0,d7 1107 movel a4,a0 1108 movel d6,a4 1109 movel a3,d6 1110 movel a4,a3 1111 #endif 1112 addl d7,d3 | 1113 addxl d6,d2 | 1114 addxl d5,d1 | 1115 addxl d4,d0 | 1116 1117 movel a2,d4 | return exponent to d4 1118 movel a0,d7 | 1119 andl IMM (0x80000000),d7 | d7 now has the sign 1120 1121 #ifndef __mcoldfire__ 1122 moveml sp@+,a2-a3 1123 #else 1124 movel sp@+,a4 1125 movel sp@+,a3 1126 movel sp@+,a2 1127 #endif 1128 1129 | Before rounding normalize so bit #DBL_MANT_DIG is set (we will consider 1130 | the case of denormalized numbers in the rounding routine itself). 1131 | As in the addition (not in the subtraction!) we could have set 1132 | one more bit we check this: 1133 btst IMM (DBL_MANT_DIG+1),d0 1134 beq 1f 1135 #ifndef __mcoldfire__ 1136 lsrl IMM (1),d0 1137 roxrl IMM (1),d1 1138 roxrl IMM (1),d2 1139 roxrl IMM (1),d3 1140 addw IMM (1),d4 1141 #else 1142 lsrl IMM (1),d3 1143 btst IMM (0),d2 1144 beq 10f 1145 bset IMM (31),d3 1146 10: lsrl IMM (1),d2 1147 btst IMM (0),d1 1148 beq 11f 1149 bset IMM (31),d2 1150 11: lsrl IMM (1),d1 1151 btst IMM (0),d0 1152 beq 12f 1153 bset IMM (31),d1 1154 12: lsrl IMM (1),d0 1155 addl IMM (1),d4 1156 #endif 1157 1: 1158 lea pc@(Ladddf$5),a0 | to return from rounding routine 1159 PICLEA SYM (_fpCCR),a1 | check the rounding mode 1160 #ifdef __mcoldfire__ 1161 clrl d6 1162 #endif 1163 movew a1@(6),d6 | rounding mode in d6 1164 beq Lround$to$nearest 1165 #ifndef __mcoldfire__ 1166 cmpw IMM (ROUND_TO_PLUS),d6 1167 #else 1168 cmpl IMM (ROUND_TO_PLUS),d6 1169 #endif 1170 bhi Lround$to$minus 1171 blt Lround$to$zero 1172 bra Lround$to$plus 1173 Ladddf$5: 1174 | Put back the exponent and check for overflow 1175 #ifndef __mcoldfire__ 1176 cmpw IMM (0x7ff),d4 | is the exponent big? 1177 #else 1178 cmpl IMM (0x7ff),d4 | is the exponent big? 1179 #endif 1180 bge 1f 1181 bclr IMM (DBL_MANT_DIG-1),d0 1182 #ifndef __mcoldfire__ 1183 lslw IMM (4),d4 | put exponent back into position 1184 #else 1185 lsll IMM (4),d4 | put exponent back into position 1186 #endif 1187 swap d0 | 1188 #ifndef __mcoldfire__ 1189 orw d4,d0 | 1190 #else 1191 orl d4,d0 | 1192 #endif 1193 swap d0 | 1194 bra Ladddf$ret 1195 1: 1196 moveq IMM (ADD),d5 1197 bra Ld$overflow 1198 1199 Lsubdf$0: 1200 | Here we do the subtraction. 1201 #ifndef __mcoldfire__ 1202 exg d7,a0 | put sign back in a0 1203 exg d6,a3 | 1204 #else 1205 movel d7,a4 1206 movel a0,d7 1207 movel a4,a0 1208 movel d6,a4 1209 movel a3,d6 1210 movel a4,a3 1211 #endif 1212 subl d7,d3 | 1213 subxl d6,d2 | 1214 subxl d5,d1 | 1215 subxl d4,d0 | 1216 beq Ladddf$ret$1 | if zero just exit 1217 bpl 1f | if positive skip the following 1218 movel a0,d7 | 1219 bchg IMM (31),d7 | change sign bit in d7 1220 movel d7,a0 | 1221 negl d3 | 1222 negxl d2 | 1223 negxl d1 | and negate result 1224 negxl d0 | 1225 1: 1226 movel a2,d4 | return exponent to d4 1227 movel a0,d7 1228 andl IMM (0x80000000),d7 | isolate sign bit 1229 #ifndef __mcoldfire__ 1230 moveml sp@+,a2-a3 | 1231 #else 1232 movel sp@+,a4 1233 movel sp@+,a3 1234 movel sp@+,a2 1235 #endif 1236 1237 | Before rounding normalize so bit #DBL_MANT_DIG is set (we will consider 1238 | the case of denormalized numbers in the rounding routine itself). 1239 | As in the addition (not in the subtraction!) we could have set 1240 | one more bit we check this: 1241 btst IMM (DBL_MANT_DIG+1),d0 1242 beq 1f 1243 #ifndef __mcoldfire__ 1244 lsrl IMM (1),d0 1245 roxrl IMM (1),d1 1246 roxrl IMM (1),d2 1247 roxrl IMM (1),d3 1248 addw IMM (1),d4 1249 #else 1250 lsrl IMM (1),d3 1251 btst IMM (0),d2 1252 beq 10f 1253 bset IMM (31),d3 1254 10: lsrl IMM (1),d2 1255 btst IMM (0),d1 1256 beq 11f 1257 bset IMM (31),d2 1258 11: lsrl IMM (1),d1 1259 btst IMM (0),d0 1260 beq 12f 1261 bset IMM (31),d1 1262 12: lsrl IMM (1),d0 1263 addl IMM (1),d4 1264 #endif 1265 1: 1266 lea pc@(Lsubdf$1),a0 | to return from rounding routine 1267 PICLEA SYM (_fpCCR),a1 | check the rounding mode 1268 #ifdef __mcoldfire__ 1269 clrl d6 1270 #endif 1271 movew a1@(6),d6 | rounding mode in d6 1272 beq Lround$to$nearest 1273 #ifndef __mcoldfire__ 1274 cmpw IMM (ROUND_TO_PLUS),d6 1275 #else 1276 cmpl IMM (ROUND_TO_PLUS),d6 1277 #endif 1278 bhi Lround$to$minus 1279 blt Lround$to$zero 1280 bra Lround$to$plus 1281 Lsubdf$1: 1282 | Put back the exponent and sign (we don't have overflow). ' 1283 bclr IMM (DBL_MANT_DIG-1),d0 1284 #ifndef __mcoldfire__ 1285 lslw IMM (4),d4 | put exponent back into position 1286 #else 1287 lsll IMM (4),d4 | put exponent back into position 1288 #endif 1289 swap d0 | 1290 #ifndef __mcoldfire__ 1291 orw d4,d0 | 1292 #else 1293 orl d4,d0 | 1294 #endif 1295 swap d0 | 1296 bra Ladddf$ret 1297 1298 | If one of the numbers was too small (difference of exponents >= 1299 | DBL_MANT_DIG+1) we return the other (and now we don't have to ' 1300 | check for finiteness or zero). 1301 Ladddf$a$small: 1302 #ifndef __mcoldfire__ 1303 moveml sp@+,a2-a3 1304 #else 1305 movel sp@+,a4 1306 movel sp@+,a3 1307 movel sp@+,a2 1308 #endif 1309 movel a6@(16),d0 1310 movel a6@(20),d1 1311 PICLEA SYM (_fpCCR),a0 1312 movew IMM (0),a0@ 1313 #ifndef __mcoldfire__ 1314 moveml sp@+,d2-d7 | restore data registers 1315 #else 1316 moveml sp@,d2-d7 1317 | XXX if frame pointer is ever removed, stack pointer must 1318 | be adjusted here. 1319 #endif 1320 unlk a6 | and return 1321 rts 1322 1323 Ladddf$b$small: 1324 #ifndef __mcoldfire__ 1325 moveml sp@+,a2-a3 1326 #else 1327 movel sp@+,a4 1328 movel sp@+,a3 1329 movel sp@+,a2 1330 #endif 1331 movel a6@(8),d0 1332 movel a6@(12),d1 1333 PICLEA SYM (_fpCCR),a0 1334 movew IMM (0),a0@ 1335 #ifndef __mcoldfire__ 1336 moveml sp@+,d2-d7 | restore data registers 1337 #else 1338 moveml sp@,d2-d7 1339 | XXX if frame pointer is ever removed, stack pointer must 1340 | be adjusted here. 1341 #endif 1342 unlk a6 | and return 1343 rts 1344 1345 Ladddf$a$den: 1346 movel d7,d4 | d7 contains 0x00200000 1347 bra Ladddf$1 1348 1349 Ladddf$b$den: 1350 movel d7,d5 | d7 contains 0x00200000 1351 notl d6 1352 bra Ladddf$2 1353 1354 Ladddf$b: 1355 | Return b (if a is zero) 1356 movel d2,d0 1357 movel d3,d1 1358 bne 1f | Check if b is -0 1359 cmpl IMM (0x80000000),d0 1360 bne 1f 1361 andl IMM (0x80000000),d7 | Use the sign of a 1362 clrl d0 1363 bra Ladddf$ret 1364 Ladddf$a: 1365 movel a6@(8),d0 1366 movel a6@(12),d1 1367 1: 1368 moveq IMM (ADD),d5 1369 | Check for NaN and +/-INFINITY. 1370 movel d0,d7 | 1371 andl IMM (0x80000000),d7 | 1372 bclr IMM (31),d0 | 1373 cmpl IMM (0x7ff00000),d0 | 1374 bge 2f | 1375 movel d0,d0 | check for zero, since we don't ' 1376 bne Ladddf$ret | want to return -0 by mistake 1377 bclr IMM (31),d7 | 1378 bra Ladddf$ret | 1379 2: 1380 andl IMM (0x000fffff),d0 | check for NaN (nonzero fraction) 1381 orl d1,d0 | 1382 bne Ld$inop | 1383 bra Ld$infty | 1384 1385 Ladddf$ret$1: 1386 #ifndef __mcoldfire__ 1387 moveml sp@+,a2-a3 | restore regs and exit 1388 #else 1389 movel sp@+,a4 1390 movel sp@+,a3 1391 movel sp@+,a2 1392 #endif 1393 1394 Ladddf$ret: 1395 | Normal exit. 1396 PICLEA SYM (_fpCCR),a0 1397 movew IMM (0),a0@ 1398 orl d7,d0 | put sign bit back 1399 #ifndef __mcoldfire__ 1400 moveml sp@+,d2-d7 1401 #else 1402 moveml sp@,d2-d7 1403 | XXX if frame pointer is ever removed, stack pointer must 1404 | be adjusted here. 1405 #endif 1406 unlk a6 1407 rts 1408 1409 Ladddf$ret$den: 1410 | Return a denormalized number. 1411 #ifndef __mcoldfire__ 1412 lsrl IMM (1),d0 | shift right once more 1413 roxrl IMM (1),d1 | 1414 #else 1415 lsrl IMM (1),d1 1416 btst IMM (0),d0 1417 beq 10f 1418 bset IMM (31),d1 1419 10: lsrl IMM (1),d0 1420 #endif 1421 bra Ladddf$ret 1422 1423 Ladddf$nf: 1424 moveq IMM (ADD),d5 1425 | This could be faster but it is not worth the effort, since it is not 1426 | executed very often. We sacrifice speed for clarity here. 1427 movel a6@(8),d0 | get the numbers back (remember that we 1428 movel a6@(12),d1 | did some processing already) 1429 movel a6@(16),d2 | 1430 movel a6@(20),d3 | 1431 movel IMM (0x7ff00000),d4 | useful constant (INFINITY) 1432 movel d0,d7 | save sign bits 1433 movel d2,d6 | 1434 bclr IMM (31),d0 | clear sign bits 1435 bclr IMM (31),d2 | 1436 | We know that one of them is either NaN of +/-INFINITY 1437 | Check for NaN (if either one is NaN return NaN) 1438 cmpl d4,d0 | check first a (d0) 1439 bhi Ld$inop | if d0 > 0x7ff00000 or equal and 1440 bne 2f 1441 tstl d1 | d1 > 0, a is NaN 1442 bne Ld$inop | 1443 2: cmpl d4,d2 | check now b (d1) 1444 bhi Ld$inop | 1445 bne 3f 1446 tstl d3 | 1447 bne Ld$inop | 1448 3: 1449 | Now comes the check for +/-INFINITY. We know that both are (maybe not 1450 | finite) numbers, but we have to check if both are infinite whether we 1451 | are adding or subtracting them. 1452 eorl d7,d6 | to check sign bits 1453 bmi 1f 1454 andl IMM (0x80000000),d7 | get (common) sign bit 1455 bra Ld$infty 1456 1: 1457 | We know one (or both) are infinite, so we test for equality between the 1458 | two numbers (if they are equal they have to be infinite both, so we 1459 | return NaN). 1460 cmpl d2,d0 | are both infinite? 1461 bne 1f | if d0 <> d2 they are not equal 1462 cmpl d3,d1 | if d0 == d2 test d3 and d1 1463 beq Ld$inop | if equal return NaN 1464 1: 1465 andl IMM (0x80000000),d7 | get a's sign bit ' 1466 cmpl d4,d0 | test now for infinity 1467 beq Ld$infty | if a is INFINITY return with this sign 1468 bchg IMM (31),d7 | else we know b is INFINITY and has 1469 bra Ld$infty | the opposite sign 1470 1471 |============================================================================= 1472 | __muldf3 1473 |============================================================================= 1474 1475 | double __muldf3(double, double); 1476 FUNC(__muldf3) 1477 SYM (__muldf3): 1478 #ifndef __mcoldfire__ 1479 link a6,IMM (0) 1480 moveml d2-d7,sp@- 1481 #else 1482 link a6,IMM (-24) 1483 moveml d2-d7,sp@ 1484 #endif 1485 movel a6@(8),d0 | get a into d0-d1 1486 movel a6@(12),d1 | 1487 movel a6@(16),d2 | and b into d2-d3 1488 movel a6@(20),d3 | 1489 movel d0,d7 | d7 will hold the sign of the product 1490 eorl d2,d7 | 1491 andl IMM (0x80000000),d7 | 1492 movel d7,a0 | save sign bit into a0 1493 movel IMM (0x7ff00000),d7 | useful constant (+INFINITY) 1494 movel d7,d6 | another (mask for fraction) 1495 notl d6 | 1496 bclr IMM (31),d0 | get rid of a's sign bit ' 1497 movel d0,d4 | 1498 orl d1,d4 | 1499 beq Lmuldf$a$0 | branch if a is zero 1500 movel d0,d4 | 1501 bclr IMM (31),d2 | get rid of b's sign bit ' 1502 movel d2,d5 | 1503 orl d3,d5 | 1504 beq Lmuldf$b$0 | branch if b is zero 1505 movel d2,d5 | 1506 cmpl d7,d0 | is a big? 1507 bhi Lmuldf$inop | if a is NaN return NaN 1508 beq Lmuldf$a$nf | we still have to check d1 and b ... 1509 cmpl d7,d2 | now compare b with INFINITY 1510 bhi Lmuldf$inop | is b NaN? 1511 beq Lmuldf$b$nf | we still have to check d3 ... 1512 | Here we have both numbers finite and nonzero (and with no sign bit). 1513 | Now we get the exponents into d4 and d5. 1514 andl d7,d4 | isolate exponent in d4 1515 beq Lmuldf$a$den | if exponent zero, have denormalized 1516 andl d6,d0 | isolate fraction 1517 orl IMM (0x00100000),d0 | and put hidden bit back 1518 swap d4 | I like exponents in the first byte 1519 #ifndef __mcoldfire__ 1520 lsrw IMM (4),d4 | 1521 #else 1522 lsrl IMM (4),d4 | 1523 #endif 1524 Lmuldf$1: 1525 andl d7,d5 | 1526 beq Lmuldf$b$den | 1527 andl d6,d2 | 1528 orl IMM (0x00100000),d2 | and put hidden bit back 1529 swap d5 | 1530 #ifndef __mcoldfire__ 1531 lsrw IMM (4),d5 | 1532 #else 1533 lsrl IMM (4),d5 | 1534 #endif 1535 Lmuldf$2: | 1536 #ifndef __mcoldfire__ 1537 addw d5,d4 | add exponents 1538 subw IMM (D_BIAS+1),d4 | and subtract bias (plus one) 1539 #else 1540 addl d5,d4 | add exponents 1541 subl IMM (D_BIAS+1),d4 | and subtract bias (plus one) 1542 #endif 1543 1544 | We are now ready to do the multiplication. The situation is as follows: 1545 | both a and b have bit 52 ( bit 20 of d0 and d2) set (even if they were 1546 | denormalized to start with!), which means that in the product bit 104 1547 | (which will correspond to bit 8 of the fourth long) is set. 1548 1549 | Here we have to do the product. 1550 | To do it we have to juggle the registers back and forth, as there are not 1551 | enough to keep everything in them. So we use the address registers to keep 1552 | some intermediate data. 1553 1554 #ifndef __mcoldfire__ 1555 moveml a2-a3,sp@- | save a2 and a3 for temporary use 1556 #else 1557 movel a2,sp@- 1558 movel a3,sp@- 1559 movel a4,sp@- 1560 #endif 1561 movel IMM (0),a2 | a2 is a null register 1562 movel d4,a3 | and a3 will preserve the exponent 1563 1564 | First, shift d2-d3 so bit 20 becomes bit 31: 1565 #ifndef __mcoldfire__ 1566 rorl IMM (5),d2 | rotate d2 5 places right 1567 swap d2 | and swap it 1568 rorl IMM (5),d3 | do the same thing with d3 1569 swap d3 | 1570 movew d3,d6 | get the rightmost 11 bits of d3 1571 andw IMM (0x07ff),d6 | 1572 orw d6,d2 | and put them into d2 1573 andw IMM (0xf800),d3 | clear those bits in d3 1574 #else 1575 moveq IMM (11),d7 | left shift d2 11 bits 1576 lsll d7,d2 1577 movel d3,d6 | get a copy of d3 1578 lsll d7,d3 | left shift d3 11 bits 1579 andl IMM (0xffe00000),d6 | get the top 11 bits of d3 1580 moveq IMM (21),d7 | right shift them 21 bits 1581 lsrl d7,d6 1582 orl d6,d2 | stick them at the end of d2 1583 #endif 1584 1585 movel d2,d6 | move b into d6-d7 1586 movel d3,d7 | move a into d4-d5 1587 movel d0,d4 | and clear d0-d1-d2-d3 (to put result) 1588 movel d1,d5 | 1589 movel IMM (0),d3 | 1590 movel d3,d2 | 1591 movel d3,d1 | 1592 movel d3,d0 | 1593 1594 | We use a1 as counter: 1595 movel IMM (DBL_MANT_DIG-1),a1 1596 #ifndef __mcoldfire__ 1597 exg d7,a1 1598 #else 1599 movel d7,a4 1600 movel a1,d7 1601 movel a4,a1 1602 #endif 1603 1604 1: 1605 #ifndef __mcoldfire__ 1606 exg d7,a1 | put counter back in a1 1607 #else 1608 movel d7,a4 1609 movel a1,d7 1610 movel a4,a1 1611 #endif 1612 addl d3,d3 | shift sum once left 1613 addxl d2,d2 | 1614 addxl d1,d1 | 1615 addxl d0,d0 | 1616 addl d7,d7 | 1617 addxl d6,d6 | 1618 bcc 2f | if bit clear skip the following 1619 #ifndef __mcoldfire__ 1620 exg d7,a2 | 1621 #else 1622 movel d7,a4 1623 movel a2,d7 1624 movel a4,a2 1625 #endif 1626 addl d5,d3 | else add a to the sum 1627 addxl d4,d2 | 1628 addxl d7,d1 | 1629 addxl d7,d0 | 1630 #ifndef __mcoldfire__ 1631 exg d7,a2 | 1632 #else 1633 movel d7,a4 1634 movel a2,d7 1635 movel a4,a2 1636 #endif 1637 2: 1638 #ifndef __mcoldfire__ 1639 exg d7,a1 | put counter in d7 1640 dbf d7,1b | decrement and branch 1641 #else 1642 movel d7,a4 1643 movel a1,d7 1644 movel a4,a1 1645 subql IMM (1),d7 1646 bpl 1b 1647 #endif 1648 1649 movel a3,d4 | restore exponent 1650 #ifndef __mcoldfire__ 1651 moveml sp@+,a2-a3 1652 #else 1653 movel sp@+,a4 1654 movel sp@+,a3 1655 movel sp@+,a2 1656 #endif 1657 1658 | Now we have the product in d0-d1-d2-d3, with bit 8 of d0 set. The 1659 | first thing to do now is to normalize it so bit 8 becomes bit 1660 | DBL_MANT_DIG-32 (to do the rounding); later we will shift right. 1661 swap d0 1662 swap d1 1663 movew d1,d0 1664 swap d2 1665 movew d2,d1 1666 swap d3 1667 movew d3,d2 1668 movew IMM (0),d3 1669 #ifndef __mcoldfire__ 1670 lsrl IMM (1),d0 1671 roxrl IMM (1),d1 1672 roxrl IMM (1),d2 1673 roxrl IMM (1),d3 1674 lsrl IMM (1),d0 1675 roxrl IMM (1),d1 1676 roxrl IMM (1),d2 1677 roxrl IMM (1),d3 1678 lsrl IMM (1),d0 1679 roxrl IMM (1),d1 1680 roxrl IMM (1),d2 1681 roxrl IMM (1),d3 1682 #else 1683 moveq IMM (29),d6 1684 lsrl IMM (3),d3 1685 movel d2,d7 1686 lsll d6,d7 1687 orl d7,d3 1688 lsrl IMM (3),d2 1689 movel d1,d7 1690 lsll d6,d7 1691 orl d7,d2 1692 lsrl IMM (3),d1 1693 movel d0,d7 1694 lsll d6,d7 1695 orl d7,d1 1696 lsrl IMM (3),d0 1697 #endif 1698 1699 | Now round, check for over- and underflow, and exit. 1700 movel a0,d7 | get sign bit back into d7 1701 moveq IMM (MULTIPLY),d5 1702 1703 btst IMM (DBL_MANT_DIG+1-32),d0 1704 beq Lround$exit 1705 #ifndef __mcoldfire__ 1706 lsrl IMM (1),d0 1707 roxrl IMM (1),d1 1708 addw IMM (1),d4 1709 #else 1710 lsrl IMM (1),d1 1711 btst IMM (0),d0 1712 beq 10f 1713 bset IMM (31),d1 1714 10: lsrl IMM (1),d0 1715 addl IMM (1),d4 1716 #endif 1717 bra Lround$exit 1718 1719 Lmuldf$inop: 1720 moveq IMM (MULTIPLY),d5 1721 bra Ld$inop 1722 1723 Lmuldf$b$nf: 1724 moveq IMM (MULTIPLY),d5 1725 movel a0,d7 | get sign bit back into d7 1726 tstl d3 | we know d2 == 0x7ff00000, so check d3 1727 bne Ld$inop | if d3 <> 0 b is NaN 1728 bra Ld$overflow | else we have overflow (since a is finite) 1729 1730 Lmuldf$a$nf: 1731 moveq IMM (MULTIPLY),d5 1732 movel a0,d7 | get sign bit back into d7 1733 tstl d1 | we know d0 == 0x7ff00000, so check d1 1734 bne Ld$inop | if d1 <> 0 a is NaN 1735 bra Ld$overflow | else signal overflow 1736 1737 | If either number is zero return zero, unless the other is +/-INFINITY or 1738 | NaN, in which case we return NaN. 1739 Lmuldf$b$0: 1740 moveq IMM (MULTIPLY),d5 1741 #ifndef __mcoldfire__ 1742 exg d2,d0 | put b (==0) into d0-d1 1743 exg d3,d1 | and a (with sign bit cleared) into d2-d3 1744 movel a0,d0 | set result sign 1745 #else 1746 movel d0,d2 | put a into d2-d3 1747 movel d1,d3 1748 movel a0,d0 | put result zero into d0-d1 1749 movq IMM(0),d1 1750 #endif 1751 bra 1f 1752 Lmuldf$a$0: 1753 movel a0,d0 | set result sign 1754 movel a6@(16),d2 | put b into d2-d3 again 1755 movel a6@(20),d3 | 1756 bclr IMM (31),d2 | clear sign bit 1757 1: cmpl IMM (0x7ff00000),d2 | check for non-finiteness 1758 bge Ld$inop | in case NaN or +/-INFINITY return NaN 1759 PICLEA SYM (_fpCCR),a0 1760 movew IMM (0),a0@ 1761 #ifndef __mcoldfire__ 1762 moveml sp@+,d2-d7 1763 #else 1764 moveml sp@,d2-d7 1765 | XXX if frame pointer is ever removed, stack pointer must 1766 | be adjusted here. 1767 #endif 1768 unlk a6 1769 rts 1770 1771 | If a number is denormalized we put an exponent of 1 but do not put the 1772 | hidden bit back into the fraction; instead we shift left until bit 21 1773 | (the hidden bit) is set, adjusting the exponent accordingly. We do this 1774 | to ensure that the product of the fractions is close to 1. 1775 Lmuldf$a$den: 1776 movel IMM (1),d4 1777 andl d6,d0 1778 1: addl d1,d1 | shift a left until bit 20 is set 1779 addxl d0,d0 | 1780 #ifndef __mcoldfire__ 1781 subw IMM (1),d4 | and adjust exponent 1782 #else 1783 subl IMM (1),d4 | and adjust exponent 1784 #endif 1785 btst IMM (20),d0 | 1786 bne Lmuldf$1 | 1787 bra 1b 1788 1789 Lmuldf$b$den: 1790 movel IMM (1),d5 1791 andl d6,d2 1792 1: addl d3,d3 | shift b left until bit 20 is set 1793 addxl d2,d2 | 1794 #ifndef __mcoldfire__ 1795 subw IMM (1),d5 | and adjust exponent 1796 #else 1797 subql IMM (1),d5 | and adjust exponent 1798 #endif 1799 btst IMM (20),d2 | 1800 bne Lmuldf$2 | 1801 bra 1b 1802 1803 1804 |============================================================================= 1805 | __divdf3 1806 |============================================================================= 1807 1808 | double __divdf3(double, double); 1809 FUNC(__divdf3) 1810 SYM (__divdf3): 1811 #ifndef __mcoldfire__ 1812 link a6,IMM (0) 1813 moveml d2-d7,sp@- 1814 #else 1815 link a6,IMM (-24) 1816 moveml d2-d7,sp@ 1817 #endif 1818 movel a6@(8),d0 | get a into d0-d1 1819 movel a6@(12),d1 | 1820 movel a6@(16),d2 | and b into d2-d3 1821 movel a6@(20),d3 | 1822 movel d0,d7 | d7 will hold the sign of the result 1823 eorl d2,d7 | 1824 andl IMM (0x80000000),d7 1825 movel d7,a0 | save sign into a0 1826 movel IMM (0x7ff00000),d7 | useful constant (+INFINITY) 1827 movel d7,d6 | another (mask for fraction) 1828 notl d6 | 1829 bclr IMM (31),d0 | get rid of a's sign bit ' 1830 movel d0,d4 | 1831 orl d1,d4 | 1832 beq Ldivdf$a$0 | branch if a is zero 1833 movel d0,d4 | 1834 bclr IMM (31),d2 | get rid of b's sign bit ' 1835 movel d2,d5 | 1836 orl d3,d5 | 1837 beq Ldivdf$b$0 | branch if b is zero 1838 movel d2,d5 1839 cmpl d7,d0 | is a big? 1840 bhi Ldivdf$inop | if a is NaN return NaN 1841 beq Ldivdf$a$nf | if d0 == 0x7ff00000 we check d1 1842 cmpl d7,d2 | now compare b with INFINITY 1843 bhi Ldivdf$inop | if b is NaN return NaN 1844 beq Ldivdf$b$nf | if d2 == 0x7ff00000 we check d3 1845 | Here we have both numbers finite and nonzero (and with no sign bit). 1846 | Now we get the exponents into d4 and d5 and normalize the numbers to 1847 | ensure that the ratio of the fractions is around 1. We do this by 1848 | making sure that both numbers have bit #DBL_MANT_DIG-32-1 (hidden bit) 1849 | set, even if they were denormalized to start with. 1850 | Thus, the result will satisfy: 2 > result > 1/2. 1851 andl d7,d4 | and isolate exponent in d4 1852 beq Ldivdf$a$den | if exponent is zero we have a denormalized 1853 andl d6,d0 | and isolate fraction 1854 orl IMM (0x00100000),d0 | and put hidden bit back 1855 swap d4 | I like exponents in the first byte 1856 #ifndef __mcoldfire__ 1857 lsrw IMM (4),d4 | 1858 #else 1859 lsrl IMM (4),d4 | 1860 #endif 1861 Ldivdf$1: | 1862 andl d7,d5 | 1863 beq Ldivdf$b$den | 1864 andl d6,d2 | 1865 orl IMM (0x00100000),d2 1866 swap d5 | 1867 #ifndef __mcoldfire__ 1868 lsrw IMM (4),d5 | 1869 #else 1870 lsrl IMM (4),d5 | 1871 #endif 1872 Ldivdf$2: | 1873 #ifndef __mcoldfire__ 1874 subw d5,d4 | subtract exponents 1875 addw IMM (D_BIAS),d4 | and add bias 1876 #else 1877 subl d5,d4 | subtract exponents 1878 addl IMM (D_BIAS),d4 | and add bias 1879 #endif 1880 1881 | We are now ready to do the division. We have prepared things in such a way 1882 | that the ratio of the fractions will be less than 2 but greater than 1/2. 1883 | At this point the registers in use are: 1884 | d0-d1 hold a (first operand, bit DBL_MANT_DIG-32=0, bit 1885 | DBL_MANT_DIG-1-32=1) 1886 | d2-d3 hold b (second operand, bit DBL_MANT_DIG-32=1) 1887 | d4 holds the difference of the exponents, corrected by the bias 1888 | a0 holds the sign of the ratio 1889 1890 | To do the rounding correctly we need to keep information about the 1891 | nonsignificant bits. One way to do this would be to do the division 1892 | using four registers; another is to use two registers (as originally 1893 | I did), but use a sticky bit to preserve information about the 1894 | fractional part. Note that we can keep that info in a1, which is not 1895 | used. 1896 movel IMM (0),d6 | d6-d7 will hold the result 1897 movel d6,d7 | 1898 movel IMM (0),a1 | and a1 will hold the sticky bit 1899 1900 movel IMM (DBL_MANT_DIG-32+1),d5 1901 1902 1: cmpl d0,d2 | is a < b? 1903 bhi 3f | if b > a skip the following 1904 beq 4f | if d0==d2 check d1 and d3 1905 2: subl d3,d1 | 1906 subxl d2,d0 | a <-- a - b 1907 bset d5,d6 | set the corresponding bit in d6 1908 3: addl d1,d1 | shift a by 1 1909 addxl d0,d0 | 1910 #ifndef __mcoldfire__ 1911 dbra d5,1b | and branch back 1912 #else 1913 subql IMM (1), d5 1914 bpl 1b 1915 #endif 1916 bra 5f 1917 4: cmpl d1,d3 | here d0==d2, so check d1 and d3 1918 bhi 3b | if d1 > d2 skip the subtraction 1919 bra 2b | else go do it 1920 5: 1921 | Here we have to start setting the bits in the second long. 1922 movel IMM (31),d5 | again d5 is counter 1923 1924 1: cmpl d0,d2 | is a < b? 1925 bhi 3f | if b > a skip the following 1926 beq 4f | if d0==d2 check d1 and d3 1927 2: subl d3,d1 | 1928 subxl d2,d0 | a <-- a - b 1929 bset d5,d7 | set the corresponding bit in d7 1930 3: addl d1,d1 | shift a by 1 1931 addxl d0,d0 | 1932 #ifndef __mcoldfire__ 1933 dbra d5,1b | and branch back 1934 #else 1935 subql IMM (1), d5 1936 bpl 1b 1937 #endif 1938 bra 5f 1939 4: cmpl d1,d3 | here d0==d2, so check d1 and d3 1940 bhi 3b | if d1 > d2 skip the subtraction 1941 bra 2b | else go do it 1942 5: 1943 | Now go ahead checking until we hit a one, which we store in d2. 1944 movel IMM (DBL_MANT_DIG),d5 1945 1: cmpl d2,d0 | is a < b? 1946 bhi 4f | if b < a, exit 1947 beq 3f | if d0==d2 check d1 and d3 1948 2: addl d1,d1 | shift a by 1 1949 addxl d0,d0 | 1950 #ifndef __mcoldfire__ 1951 dbra d5,1b | and branch back 1952 #else 1953 subql IMM (1), d5 1954 bpl 1b 1955 #endif 1956 movel IMM (0),d2 | here no sticky bit was found 1957 movel d2,d3 1958 bra 5f 1959 3: cmpl d1,d3 | here d0==d2, so check d1 and d3 1960 bhi 2b | if d1 > d2 go back 1961 4: 1962 | Here put the sticky bit in d2-d3 (in the position which actually corresponds 1963 | to it; if you don't do this the algorithm loses in some cases). ' 1964 movel IMM (0),d2 1965 movel d2,d3 1966 #ifndef __mcoldfire__ 1967 subw IMM (DBL_MANT_DIG),d5 1968 addw IMM (63),d5 1969 cmpw IMM (31),d5 1970 #else 1971 subl IMM (DBL_MANT_DIG),d5 1972 addl IMM (63),d5 1973 cmpl IMM (31),d5 1974 #endif 1975 bhi 2f 1976 1: bset d5,d3 1977 bra 5f 1978 #ifndef __mcoldfire__ 1979 subw IMM (32),d5 1980 #else 1981 subl IMM (32),d5 1982 #endif 1983 2: bset d5,d2 1984 5: 1985 | Finally we are finished! Move the longs in the address registers to 1986 | their final destination: 1987 movel d6,d0 1988 movel d7,d1 1989 movel IMM (0),d3 1990 1991 | Here we have finished the division, with the result in d0-d1-d2-d3, with 1992 | 2^21 <= d6 < 2^23. Thus bit 23 is not set, but bit 22 could be set. 1993 | If it is not, then definitely bit 21 is set. Normalize so bit 22 is 1994 | not set: 1995 btst IMM (DBL_MANT_DIG-32+1),d0 1996 beq 1f 1997 #ifndef __mcoldfire__ 1998 lsrl IMM (1),d0 1999 roxrl IMM (1),d1 2000 roxrl IMM (1),d2 2001 roxrl IMM (1),d3 2002 addw IMM (1),d4 2003 #else 2004 lsrl IMM (1),d3 2005 btst IMM (0),d2 2006 beq 10f 2007 bset IMM (31),d3 2008 10: lsrl IMM (1),d2 2009 btst IMM (0),d1 2010 beq 11f 2011 bset IMM (31),d2 2012 11: lsrl IMM (1),d1 2013 btst IMM (0),d0 2014 beq 12f 2015 bset IMM (31),d1 2016 12: lsrl IMM (1),d0 2017 addl IMM (1),d4 2018 #endif 2019 1: 2020 | Now round, check for over- and underflow, and exit. 2021 movel a0,d7 | restore sign bit to d7 2022 moveq IMM (DIVIDE),d5 2023 bra Lround$exit 2024 2025 Ldivdf$inop: 2026 moveq IMM (DIVIDE),d5 2027 bra Ld$inop 2028 2029 Ldivdf$a$0: 2030 | If a is zero check to see whether b is zero also. In that case return 2031 | NaN; then check if b is NaN, and return NaN also in that case. Else 2032 | return a properly signed zero. 2033 moveq IMM (DIVIDE),d5 2034 bclr IMM (31),d2 | 2035 movel d2,d4 | 2036 orl d3,d4 | 2037 beq Ld$inop | if b is also zero return NaN 2038 cmpl IMM (0x7ff00000),d2 | check for NaN 2039 bhi Ld$inop | 2040 blt 1f | 2041 tstl d3 | 2042 bne Ld$inop | 2043 1: movel a0,d0 | else return signed zero 2044 moveq IMM(0),d1 | 2045 PICLEA SYM (_fpCCR),a0 | clear exception flags 2046 movew IMM (0),a0@ | 2047 #ifndef __mcoldfire__ 2048 moveml sp@+,d2-d7 | 2049 #else 2050 moveml sp@,d2-d7 | 2051 | XXX if frame pointer is ever removed, stack pointer must 2052 | be adjusted here. 2053 #endif 2054 unlk a6 | 2055 rts | 2056 2057 Ldivdf$b$0: 2058 moveq IMM (DIVIDE),d5 2059 | If we got here a is not zero. Check if a is NaN; in that case return NaN, 2060 | else return +/-INFINITY. Remember that a is in d0 with the sign bit 2061 | cleared already. 2062 movel a0,d7 | put a's sign bit back in d7 ' 2063 cmpl IMM (0x7ff00000),d0 | compare d0 with INFINITY 2064 bhi Ld$inop | if larger it is NaN 2065 tstl d1 | 2066 bne Ld$inop | 2067 bra Ld$div$0 | else signal DIVIDE_BY_ZERO 2068 2069 Ldivdf$b$nf: 2070 moveq IMM (DIVIDE),d5 2071 | If d2 == 0x7ff00000 we have to check d3. 2072 tstl d3 | 2073 bne Ld$inop | if d3 <> 0, b is NaN 2074 bra Ld$underflow | else b is +/-INFINITY, so signal underflow 2075 2076 Ldivdf$a$nf: 2077 moveq IMM (DIVIDE),d5 2078 | If d0 == 0x7ff00000 we have to check d1. 2079 tstl d1 | 2080 bne Ld$inop | if d1 <> 0, a is NaN 2081 | If a is INFINITY we have to check b 2082 cmpl d7,d2 | compare b with INFINITY 2083 bge Ld$inop | if b is NaN or INFINITY return NaN 2084 tstl d3 | 2085 bne Ld$inop | 2086 bra Ld$overflow | else return overflow 2087 2088 | If a number is denormalized we put an exponent of 1 but do not put the 2089 | bit back into the fraction. 2090 Ldivdf$a$den: 2091 movel IMM (1),d4 2092 andl d6,d0 2093 1: addl d1,d1 | shift a left until bit 20 is set 2094 addxl d0,d0 2095 #ifndef __mcoldfire__ 2096 subw IMM (1),d4 | and adjust exponent 2097 #else 2098 subl IMM (1),d4 | and adjust exponent 2099 #endif 2100 btst IMM (DBL_MANT_DIG-32-1),d0 2101 bne Ldivdf$1 2102 bra 1b 2103 2104 Ldivdf$b$den: 2105 movel IMM (1),d5 2106 andl d6,d2 2107 1: addl d3,d3 | shift b left until bit 20 is set 2108 addxl d2,d2 2109 #ifndef __mcoldfire__ 2110 subw IMM (1),d5 | and adjust exponent 2111 #else 2112 subql IMM (1),d5 | and adjust exponent 2113 #endif 2114 btst IMM (DBL_MANT_DIG-32-1),d2 2115 bne Ldivdf$2 2116 bra 1b 2117 2118 Lround$exit: 2119 | This is a common exit point for __muldf3 and __divdf3. When they enter 2120 | this point the sign of the result is in d7, the result in d0-d1, normalized 2121 | so that 2^21 <= d0 < 2^22, and the exponent is in the lower byte of d4. 2122 2123 | First check for underlow in the exponent: 2124 #ifndef __mcoldfire__ 2125 cmpw IMM (-DBL_MANT_DIG-1),d4 2126 #else 2127 cmpl IMM (-DBL_MANT_DIG-1),d4 2128 #endif 2129 blt Ld$underflow 2130 | It could happen that the exponent is less than 1, in which case the 2131 | number is denormalized. In this case we shift right and adjust the 2132 | exponent until it becomes 1 or the fraction is zero (in the latter case 2133 | we signal underflow and return zero). 2134 movel d7,a0 | 2135 movel IMM (0),d6 | use d6-d7 to collect bits flushed right 2136 movel d6,d7 | use d6-d7 to collect bits flushed right 2137 #ifndef __mcoldfire__ 2138 cmpw IMM (1),d4 | if the exponent is less than 1 we 2139 #else 2140 cmpl IMM (1),d4 | if the exponent is less than 1 we 2141 #endif 2142 bge 2f | have to shift right (denormalize) 2143 1: 2144 #ifndef __mcoldfire__ 2145 addw IMM (1),d4 | adjust the exponent 2146 lsrl IMM (1),d0 | shift right once 2147 roxrl IMM (1),d1 | 2148 roxrl IMM (1),d2 | 2149 roxrl IMM (1),d3 | 2150 roxrl IMM (1),d6 | 2151 roxrl IMM (1),d7 | 2152 cmpw IMM (1),d4 | is the exponent 1 already? 2153 #else 2154 addl IMM (1),d4 | adjust the exponent 2155 lsrl IMM (1),d7 2156 btst IMM (0),d6 2157 beq 13f 2158 bset IMM (31),d7 2159 13: lsrl IMM (1),d6 2160 btst IMM (0),d3 2161 beq 14f 2162 bset IMM (31),d6 2163 14: lsrl IMM (1),d3 2164 btst IMM (0),d2 2165 beq 10f 2166 bset IMM (31),d3 2167 10: lsrl IMM (1),d2 2168 btst IMM (0),d1 2169 beq 11f 2170 bset IMM (31),d2 2171 11: lsrl IMM (1),d1 2172 btst IMM (0),d0 2173 beq 12f 2174 bset IMM (31),d1 2175 12: lsrl IMM (1),d0 2176 cmpl IMM (1),d4 | is the exponent 1 already? 2177 #endif 2178 beq 2f | if not loop back 2179 bra 1b | 2180 bra Ld$underflow | safety check, shouldn't execute ' 2181 2: orl d6,d2 | this is a trick so we don't lose ' 2182 orl d7,d3 | the bits which were flushed right 2183 movel a0,d7 | get back sign bit into d7 2184 | Now call the rounding routine (which takes care of denormalized numbers): 2185 lea pc@(Lround$0),a0 | to return from rounding routine 2186 PICLEA SYM (_fpCCR),a1 | check the rounding mode 2187 #ifdef __mcoldfire__ 2188 clrl d6 2189 #endif 2190 movew a1@(6),d6 | rounding mode in d6 2191 beq Lround$to$nearest 2192 #ifndef __mcoldfire__ 2193 cmpw IMM (ROUND_TO_PLUS),d6 2194 #else 2195 cmpl IMM (ROUND_TO_PLUS),d6 2196 #endif 2197 bhi Lround$to$minus 2198 blt Lround$to$zero 2199 bra Lround$to$plus 2200 Lround$0: 2201 | Here we have a correctly rounded result (either normalized or denormalized). 2202 2203 | Here we should have either a normalized number or a denormalized one, and 2204 | the exponent is necessarily larger or equal to 1 (so we don't have to ' 2205 | check again for underflow!). We have to check for overflow or for a 2206 | denormalized number (which also signals underflow). 2207 | Check for overflow (i.e., exponent >= 0x7ff). 2208 #ifndef __mcoldfire__ 2209 cmpw IMM (0x07ff),d4 2210 #else 2211 cmpl IMM (0x07ff),d4 2212 #endif 2213 bge Ld$overflow 2214 | Now check for a denormalized number (exponent==0): 2215 movew d4,d4 2216 beq Ld$den 2217 1: 2218 | Put back the exponents and sign and return. 2219 #ifndef __mcoldfire__ 2220 lslw IMM (4),d4 | exponent back to fourth byte 2221 #else 2222 lsll IMM (4),d4 | exponent back to fourth byte 2223 #endif 2224 bclr IMM (DBL_MANT_DIG-32-1),d0 2225 swap d0 | and put back exponent 2226 #ifndef __mcoldfire__ 2227 orw d4,d0 | 2228 #else 2229 orl d4,d0 | 2230 #endif 2231 swap d0 | 2232 orl d7,d0 | and sign also 2233 2234 PICLEA SYM (_fpCCR),a0 2235 movew IMM (0),a0@ 2236 #ifndef __mcoldfire__ 2237 moveml sp@+,d2-d7 2238 #else 2239 moveml sp@,d2-d7 2240 | XXX if frame pointer is ever removed, stack pointer must 2241 | be adjusted here. 2242 #endif 2243 unlk a6 2244 rts 2245 2246 |============================================================================= 2247 | __negdf2 2248 |============================================================================= 2249 2250 | double __negdf2(double, double); 2251 FUNC(__negdf2) 2252 SYM (__negdf2): 2253 #ifndef __mcoldfire__ 2254 link a6,IMM (0) 2255 moveml d2-d7,sp@- 2256 #else 2257 link a6,IMM (-24) 2258 moveml d2-d7,sp@ 2259 #endif 2260 moveq IMM (NEGATE),d5 2261 movel a6@(8),d0 | get number to negate in d0-d1 2262 movel a6@(12),d1 | 2263 bchg IMM (31),d0 | negate 2264 movel d0,d2 | make a positive copy (for the tests) 2265 bclr IMM (31),d2 | 2266 movel d2,d4 | check for zero 2267 orl d1,d4 | 2268 beq 2f | if zero (either sign) return +zero 2269 cmpl IMM (0x7ff00000),d2 | compare to +INFINITY 2270 blt 1f | if finite, return 2271 bhi Ld$inop | if larger (fraction not zero) is NaN 2272 tstl d1 | if d2 == 0x7ff00000 check d1 2273 bne Ld$inop | 2274 movel d0,d7 | else get sign and return INFINITY 2275 andl IMM (0x80000000),d7 2276 bra Ld$infty 2277 1: PICLEA SYM (_fpCCR),a0 2278 movew IMM (0),a0@ 2279 #ifndef __mcoldfire__ 2280 moveml sp@+,d2-d7 2281 #else 2282 moveml sp@,d2-d7 2283 | XXX if frame pointer is ever removed, stack pointer must 2284 | be adjusted here. 2285 #endif 2286 unlk a6 2287 rts 2288 2: bclr IMM (31),d0 2289 bra 1b 2290 2291 |============================================================================= 2292 | __cmpdf2 2293 |============================================================================= 2294 2295 GREATER = 1 2296 LESS = -1 2297 EQUAL = 0 2298 2299 | int __cmpdf2_internal(double, double, int); 2300 SYM (__cmpdf2_internal): 2301 #ifndef __mcoldfire__ 2302 link a6,IMM (0) 2303 moveml d2-d7,sp@- | save registers 2304 #else 2305 link a6,IMM (-24) 2306 moveml d2-d7,sp@ 2307 #endif 2308 moveq IMM (COMPARE),d5 2309 movel a6@(8),d0 | get first operand 2310 movel a6@(12),d1 | 2311 movel a6@(16),d2 | get second operand 2312 movel a6@(20),d3 | 2313 | First check if a and/or b are (+/-) zero and in that case clear 2314 | the sign bit. 2315 movel d0,d6 | copy signs into d6 (a) and d7(b) 2316 bclr IMM (31),d0 | and clear signs in d0 and d2 2317 movel d2,d7 | 2318 bclr IMM (31),d2 | 2319 cmpl IMM (0x7ff00000),d0 | check for a == NaN 2320 bhi Lcmpd$inop | if d0 > 0x7ff00000, a is NaN 2321 beq Lcmpdf$a$nf | if equal can be INFINITY, so check d1 2322 movel d0,d4 | copy into d4 to test for zero 2323 orl d1,d4 | 2324 beq Lcmpdf$a$0 | 2325 Lcmpdf$0: 2326 cmpl IMM (0x7ff00000),d2 | check for b == NaN 2327 bhi Lcmpd$inop | if d2 > 0x7ff00000, b is NaN 2328 beq Lcmpdf$b$nf | if equal can be INFINITY, so check d3 2329 movel d2,d4 | 2330 orl d3,d4 | 2331 beq Lcmpdf$b$0 | 2332 Lcmpdf$1: 2333 | Check the signs 2334 eorl d6,d7 2335 bpl 1f 2336 | If the signs are not equal check if a >= 0 2337 tstl d6 2338 bpl Lcmpdf$a$gt$b | if (a >= 0 && b < 0) => a > b 2339 bmi Lcmpdf$b$gt$a | if (a < 0 && b >= 0) => a < b 2340 1: 2341 | If the signs are equal check for < 0 2342 tstl d6 2343 bpl 1f 2344 | If both are negative exchange them 2345 #ifndef __mcoldfire__ 2346 exg d0,d2 2347 exg d1,d3 2348 #else 2349 movel d0,d7 2350 movel d2,d0 2351 movel d7,d2 2352 movel d1,d7 2353 movel d3,d1 2354 movel d7,d3 2355 #endif 2356 1: 2357 | Now that they are positive we just compare them as longs (does this also 2358 | work for denormalized numbers?). 2359 cmpl d0,d2 2360 bhi Lcmpdf$b$gt$a | |b| > |a| 2361 bne Lcmpdf$a$gt$b | |b| < |a| 2362 | If we got here d0 == d2, so we compare d1 and d3. 2363 cmpl d1,d3 2364 bhi Lcmpdf$b$gt$a | |b| > |a| 2365 bne Lcmpdf$a$gt$b | |b| < |a| 2366 | If we got here a == b. 2367 movel IMM (EQUAL),d0 2368 #ifndef __mcoldfire__ 2369 moveml sp@+,d2-d7 | put back the registers 2370 #else 2371 moveml sp@,d2-d7 2372 | XXX if frame pointer is ever removed, stack pointer must 2373 | be adjusted here. 2374 #endif 2375 unlk a6 2376 rts 2377 Lcmpdf$a$gt$b: 2378 movel IMM (GREATER),d0 2379 #ifndef __mcoldfire__ 2380 moveml sp@+,d2-d7 | put back the registers 2381 #else 2382 moveml sp@,d2-d7 2383 | XXX if frame pointer is ever removed, stack pointer must 2384 | be adjusted here. 2385 #endif 2386 unlk a6 2387 rts 2388 Lcmpdf$b$gt$a: 2389 movel IMM (LESS),d0 2390 #ifndef __mcoldfire__ 2391 moveml sp@+,d2-d7 | put back the registers 2392 #else 2393 moveml sp@,d2-d7 2394 | XXX if frame pointer is ever removed, stack pointer must 2395 | be adjusted here. 2396 #endif 2397 unlk a6 2398 rts 2399 2400 Lcmpdf$a$0: 2401 bclr IMM (31),d6 2402 bra Lcmpdf$0 2403 Lcmpdf$b$0: 2404 bclr IMM (31),d7 2405 bra Lcmpdf$1 2406 2407 Lcmpdf$a$nf: 2408 tstl d1 2409 bne Ld$inop 2410 bra Lcmpdf$0 2411 2412 Lcmpdf$b$nf: 2413 tstl d3 2414 bne Ld$inop 2415 bra Lcmpdf$1 2416 2417 Lcmpd$inop: 2418 movl a6@(24),d0 2419 moveq IMM (INEXACT_RESULT+INVALID_OPERATION),d7 2420 moveq IMM (DOUBLE_FLOAT),d6 2421 PICJUMP $_exception_handler 2422 2423 | int __cmpdf2(double, double); 2424 FUNC(__cmpdf2) 2425 SYM (__cmpdf2): 2426 link a6,IMM (0) 2427 pea 1 2428 movl a6@(20),sp@- 2429 movl a6@(16),sp@- 2430 movl a6@(12),sp@- 2431 movl a6@(8),sp@- 2432 PICCALL SYM (__cmpdf2_internal) 2433 unlk a6 2434 rts 2435 2436 |============================================================================= 2437 | rounding routines 2438 |============================================================================= 2439 2440 | The rounding routines expect the number to be normalized in registers 2441 | d0-d1-d2-d3, with the exponent in register d4. They assume that the 2442 | exponent is larger or equal to 1. They return a properly normalized number 2443 | if possible, and a denormalized number otherwise. The exponent is returned 2444 | in d4. 2445 2446 Lround$to$nearest: 2447 | We now normalize as suggested by D. Knuth ("Seminumerical Algorithms"): 2448 | Here we assume that the exponent is not too small (this should be checked 2449 | before entering the rounding routine), but the number could be denormalized. 2450 2451 | Check for denormalized numbers: 2452 1: btst IMM (DBL_MANT_DIG-32),d0 2453 bne 2f | if set the number is normalized 2454 | Normalize shifting left until bit #DBL_MANT_DIG-32 is set or the exponent 2455 | is one (remember that a denormalized number corresponds to an 2456 | exponent of -D_BIAS+1). 2457 #ifndef __mcoldfire__ 2458 cmpw IMM (1),d4 | remember that the exponent is at least one 2459 #else 2460 cmpl IMM (1),d4 | remember that the exponent is at least one 2461 #endif 2462 beq 2f | an exponent of one means denormalized 2463 addl d3,d3 | else shift and adjust the exponent 2464 addxl d2,d2 | 2465 addxl d1,d1 | 2466 addxl d0,d0 | 2467 #ifndef __mcoldfire__ 2468 dbra d4,1b | 2469 #else 2470 subql IMM (1), d4 2471 bpl 1b 2472 #endif 2473 2: 2474 | Now round: we do it as follows: after the shifting we can write the 2475 | fraction part as f + delta, where 1 < f < 2^25, and 0 <= delta <= 2. 2476 | If delta < 1, do nothing. If delta > 1, add 1 to f. 2477 | If delta == 1, we make sure the rounded number will be even (odd?) 2478 | (after shifting). 2479 btst IMM (0),d1 | is delta < 1? 2480 beq 2f | if so, do not do anything 2481 orl d2,d3 | is delta == 1? 2482 bne 1f | if so round to even 2483 movel d1,d3 | 2484 andl IMM (2),d3 | bit 1 is the last significant bit 2485 movel IMM (0),d2 | 2486 addl d3,d1 | 2487 addxl d2,d0 | 2488 bra 2f | 2489 1: movel IMM (1),d3 | else add 1 2490 movel IMM (0),d2 | 2491 addl d3,d1 | 2492 addxl d2,d0 2493 | Shift right once (because we used bit #DBL_MANT_DIG-32!). 2494 2: 2495 #ifndef __mcoldfire__ 2496 lsrl IMM (1),d0 2497 roxrl IMM (1),d1 2498 #else 2499 lsrl IMM (1),d1 2500 btst IMM (0),d0 2501 beq 10f 2502 bset IMM (31),d1 2503 10: lsrl IMM (1),d0 2504 #endif 2505 2506 | Now check again bit #DBL_MANT_DIG-32 (rounding could have produced a 2507 | 'fraction overflow' ...). 2508 btst IMM (DBL_MANT_DIG-32),d0 2509 beq 1f 2510 #ifndef __mcoldfire__ 2511 lsrl IMM (1),d0 2512 roxrl IMM (1),d1 2513 addw IMM (1),d4 2514 #else 2515 lsrl IMM (1),d1 2516 btst IMM (0),d0 2517 beq 10f 2518 bset IMM (31),d1 2519 10: lsrl IMM (1),d0 2520 addl IMM (1),d4 2521 #endif 2522 1: 2523 | If bit #DBL_MANT_DIG-32-1 is clear we have a denormalized number, so we 2524 | have to put the exponent to zero and return a denormalized number. 2525 btst IMM (DBL_MANT_DIG-32-1),d0 2526 beq 1f 2527 jmp a0@ 2528 1: movel IMM (0),d4 2529 jmp a0@ 2530 2531 Lround$to$zero: 2532 Lround$to$plus: 2533 Lround$to$minus: 2534 jmp a0@ 2535 #endif /* L_double */ 2536 2537 #ifdef L_float 2538 2539 .globl SYM (_fpCCR) 2540 .globl $_exception_handler 2541 2542 QUIET_NaN = 0xffffffff 2543 SIGNL_NaN = 0x7f800001 2544 INFINITY = 0x7f800000 2545 2546 F_MAX_EXP = 0xff 2547 F_BIAS = 126 2548 FLT_MAX_EXP = F_MAX_EXP - F_BIAS 2549 FLT_MIN_EXP = 1 - F_BIAS 2550 FLT_MANT_DIG = 24 2551 2552 INEXACT_RESULT = 0x0001 2553 UNDERFLOW = 0x0002 2554 OVERFLOW = 0x0004 2555 DIVIDE_BY_ZERO = 0x0008 2556 INVALID_OPERATION = 0x0010 2557 2558 SINGLE_FLOAT = 1 2559 2560 NOOP = 0 2561 ADD = 1 2562 MULTIPLY = 2 2563 DIVIDE = 3 2564 NEGATE = 4 2565 COMPARE = 5 2566 EXTENDSFDF = 6 2567 TRUNCDFSF = 7 2568 2569 UNKNOWN = -1 2570 ROUND_TO_NEAREST = 0 | round result to nearest representable value 2571 ROUND_TO_ZERO = 1 | round result towards zero 2572 ROUND_TO_PLUS = 2 | round result towards plus infinity 2573 ROUND_TO_MINUS = 3 | round result towards minus infinity 2574 2575 | Entry points: 2576 2577 .globl SYM (__addsf3) 2578 .globl SYM (__subsf3) 2579 .globl SYM (__mulsf3) 2580 .globl SYM (__divsf3) 2581 .globl SYM (__negsf2) 2582 .globl SYM (__cmpsf2) 2583 .globl SYM (__cmpsf2_internal) 2584 .hidden SYM (__cmpsf2_internal) 2585 2586 | These are common routines to return and signal exceptions. 2587 2588 .text 2589 .even 2590 2591 Lf$den: 2592 | Return and signal a denormalized number 2593 orl d7,d0 2594 moveq IMM (INEXACT_RESULT+UNDERFLOW),d7 2595 moveq IMM (SINGLE_FLOAT),d6 2596 PICJUMP $_exception_handler 2597 2598 Lf$infty: 2599 Lf$overflow: 2600 | Return a properly signed INFINITY and set the exception flags 2601 movel IMM (INFINITY),d0 2602 orl d7,d0 2603 moveq IMM (INEXACT_RESULT+OVERFLOW),d7 2604 moveq IMM (SINGLE_FLOAT),d6 2605 PICJUMP $_exception_handler 2606 2607 Lf$underflow: 2608 | Return 0 and set the exception flags 2609 moveq IMM (0),d0 2610 moveq IMM (INEXACT_RESULT+UNDERFLOW),d7 2611 moveq IMM (SINGLE_FLOAT),d6 2612 PICJUMP $_exception_handler 2613 2614 Lf$inop: 2615 | Return a quiet NaN and set the exception flags 2616 movel IMM (QUIET_NaN),d0 2617 moveq IMM (INEXACT_RESULT+INVALID_OPERATION),d7 2618 moveq IMM (SINGLE_FLOAT),d6 2619 PICJUMP $_exception_handler 2620 2621 Lf$div$0: 2622 | Return a properly signed INFINITY and set the exception flags 2623 movel IMM (INFINITY),d0 2624 orl d7,d0 2625 moveq IMM (INEXACT_RESULT+DIVIDE_BY_ZERO),d7 2626 moveq IMM (SINGLE_FLOAT),d6 2627 PICJUMP $_exception_handler 2628 2629 |============================================================================= 2630 |============================================================================= 2631 | single precision routines 2632 |============================================================================= 2633 |============================================================================= 2634 2635 | A single precision floating point number (float) has the format: 2636 | 2637 | struct _float { 2638 | unsigned int sign : 1; /* sign bit */ 2639 | unsigned int exponent : 8; /* exponent, shifted by 126 */ 2640 | unsigned int fraction : 23; /* fraction */ 2641 | } float; 2642 | 2643 | Thus sizeof(float) = 4 (32 bits). 2644 | 2645 | All the routines are callable from C programs, and return the result 2646 | in the single register d0. They also preserve all registers except 2647 | d0-d1 and a0-a1. 2648 2649 |============================================================================= 2650 | __subsf3 2651 |============================================================================= 2652 2653 | float __subsf3(float, float); 2654 FUNC(__subsf3) 2655 SYM (__subsf3): 2656 bchg IMM (31),sp@(8) | change sign of second operand 2657 | and fall through 2658 |============================================================================= 2659 | __addsf3 2660 |============================================================================= 2661 2662 | float __addsf3(float, float); 2663 FUNC(__addsf3) 2664 SYM (__addsf3): 2665 #ifndef __mcoldfire__ 2666 link a6,IMM (0) | everything will be done in registers 2667 moveml d2-d7,sp@- | save all data registers but d0-d1 2668 #else 2669 link a6,IMM (-24) 2670 moveml d2-d7,sp@ 2671 #endif 2672 movel a6@(8),d0 | get first operand 2673 movel a6@(12),d1 | get second operand 2674 movel d0,a0 | get d0's sign bit ' 2675 addl d0,d0 | check and clear sign bit of a 2676 beq Laddsf$b | if zero return second operand 2677 movel d1,a1 | save b's sign bit ' 2678 addl d1,d1 | get rid of sign bit 2679 beq Laddsf$a | if zero return first operand 2680 2681 | Get the exponents and check for denormalized and/or infinity. 2682 2683 movel IMM (0x00ffffff),d4 | mask to get fraction 2684 movel IMM (0x01000000),d5 | mask to put hidden bit back 2685 2686 movel d0,d6 | save a to get exponent 2687 andl d4,d0 | get fraction in d0 2688 notl d4 | make d4 into a mask for the exponent 2689 andl d4,d6 | get exponent in d6 2690 beq Laddsf$a$den | branch if a is denormalized 2691 cmpl d4,d6 | check for INFINITY or NaN 2692 beq Laddsf$nf 2693 swap d6 | put exponent into first word 2694 orl d5,d0 | and put hidden bit back 2695 Laddsf$1: 2696 | Now we have a's exponent in d6 (second byte) and the mantissa in d0. ' 2697 movel d1,d7 | get exponent in d7 2698 andl d4,d7 | 2699 beq Laddsf$b$den | branch if b is denormalized 2700 cmpl d4,d7 | check for INFINITY or NaN 2701 beq Laddsf$nf 2702 swap d7 | put exponent into first word 2703 notl d4 | make d4 into a mask for the fraction 2704 andl d4,d1 | get fraction in d1 2705 orl d5,d1 | and put hidden bit back 2706 Laddsf$2: 2707 | Now we have b's exponent in d7 (second byte) and the mantissa in d1. ' 2708 2709 | Note that the hidden bit corresponds to bit #FLT_MANT_DIG-1, and we 2710 | shifted right once, so bit #FLT_MANT_DIG is set (so we have one extra 2711 | bit). 2712 2713 movel d1,d2 | move b to d2, since we want to use 2714 | two registers to do the sum 2715 movel IMM (0),d1 | and clear the new ones 2716 movel d1,d3 | 2717 2718 | Here we shift the numbers in registers d0 and d1 so the exponents are the 2719 | same, and put the largest exponent in d6. Note that we are using two 2720 | registers for each number (see the discussion by D. Knuth in "Seminumerical 2721 | Algorithms"). 2722 #ifndef __mcoldfire__ 2723 cmpw d6,d7 | compare exponents 2724 #else 2725 cmpl d6,d7 | compare exponents 2726 #endif 2727 beq Laddsf$3 | if equal don't shift ' 2728 bhi 5f | branch if second exponent largest 2729 1: 2730 subl d6,d7 | keep the largest exponent 2731 negl d7 2732 #ifndef __mcoldfire__ 2733 lsrw IMM (8),d7 | put difference in lower byte 2734 #else 2735 lsrl IMM (8),d7 | put difference in lower byte 2736 #endif 2737 | if difference is too large we don't shift (actually, we can just exit) ' 2738 #ifndef __mcoldfire__ 2739 cmpw IMM (FLT_MANT_DIG+2),d7 2740 #else 2741 cmpl IMM (FLT_MANT_DIG+2),d7 2742 #endif 2743 bge Laddsf$b$small 2744 #ifndef __mcoldfire__ 2745 cmpw IMM (16),d7 | if difference >= 16 swap 2746 #else 2747 cmpl IMM (16),d7 | if difference >= 16 swap 2748 #endif 2749 bge 4f 2750 2: 2751 #ifndef __mcoldfire__ 2752 subw IMM (1),d7 2753 #else 2754 subql IMM (1), d7 2755 #endif 2756 3: 2757 #ifndef __mcoldfire__ 2758 lsrl IMM (1),d2 | shift right second operand 2759 roxrl IMM (1),d3 2760 dbra d7,3b 2761 #else 2762 lsrl IMM (1),d3 2763 btst IMM (0),d2 2764 beq 10f 2765 bset IMM (31),d3 2766 10: lsrl IMM (1),d2 2767 subql IMM (1), d7 2768 bpl 3b 2769 #endif 2770 bra Laddsf$3 2771 4: 2772 movew d2,d3 2773 swap d3 2774 movew d3,d2 2775 swap d2 2776 #ifndef __mcoldfire__ 2777 subw IMM (16),d7 2778 #else 2779 subl IMM (16),d7 2780 #endif 2781 bne 2b | if still more bits, go back to normal case 2782 bra Laddsf$3 2783 5: 2784 #ifndef __mcoldfire__ 2785 exg d6,d7 | exchange the exponents 2786 #else 2787 eorl d6,d7 2788 eorl d7,d6 2789 eorl d6,d7 2790 #endif 2791 subl d6,d7 | keep the largest exponent 2792 negl d7 | 2793 #ifndef __mcoldfire__ 2794 lsrw IMM (8),d7 | put difference in lower byte 2795 #else 2796 lsrl IMM (8),d7 | put difference in lower byte 2797 #endif 2798 | if difference is too large we don't shift (and exit!) ' 2799 #ifndef __mcoldfire__ 2800 cmpw IMM (FLT_MANT_DIG+2),d7 2801 #else 2802 cmpl IMM (FLT_MANT_DIG+2),d7 2803 #endif 2804 bge Laddsf$a$small 2805 #ifndef __mcoldfire__ 2806 cmpw IMM (16),d7 | if difference >= 16 swap 2807 #else 2808 cmpl IMM (16),d7 | if difference >= 16 swap 2809 #endif 2810 bge 8f 2811 6: 2812 #ifndef __mcoldfire__ 2813 subw IMM (1),d7 2814 #else 2815 subl IMM (1),d7 2816 #endif 2817 7: 2818 #ifndef __mcoldfire__ 2819 lsrl IMM (1),d0 | shift right first operand 2820 roxrl IMM (1),d1 2821 dbra d7,7b 2822 #else 2823 lsrl IMM (1),d1 2824 btst IMM (0),d0 2825 beq 10f 2826 bset IMM (31),d1 2827 10: lsrl IMM (1),d0 2828 subql IMM (1),d7 2829 bpl 7b 2830 #endif 2831 bra Laddsf$3 2832 8: 2833 movew d0,d1 2834 swap d1 2835 movew d1,d0 2836 swap d0 2837 #ifndef __mcoldfire__ 2838 subw IMM (16),d7 2839 #else 2840 subl IMM (16),d7 2841 #endif 2842 bne 6b | if still more bits, go back to normal case 2843 | otherwise we fall through 2844 2845 | Now we have a in d0-d1, b in d2-d3, and the largest exponent in d6 (the 2846 | signs are stored in a0 and a1). 2847 2848 Laddsf$3: 2849 | Here we have to decide whether to add or subtract the numbers 2850 #ifndef __mcoldfire__ 2851 exg d6,a0 | get signs back 2852 exg d7,a1 | and save the exponents 2853 #else 2854 movel d6,d4 2855 movel a0,d6 2856 movel d4,a0 2857 movel d7,d4 2858 movel a1,d7 2859 movel d4,a1 2860 #endif 2861 eorl d6,d7 | combine sign bits 2862 bmi Lsubsf$0 | if negative a and b have opposite 2863 | sign so we actually subtract the 2864 | numbers 2865 2866 | Here we have both positive or both negative 2867 #ifndef __mcoldfire__ 2868 exg d6,a0 | now we have the exponent in d6 2869 #else 2870 movel d6,d4 2871 movel a0,d6 2872 movel d4,a0 2873 #endif 2874 movel a0,d7 | and sign in d7 2875 andl IMM (0x80000000),d7 2876 | Here we do the addition. 2877 addl d3,d1 2878 addxl d2,d0 2879 | Note: now we have d2, d3, d4 and d5 to play with! 2880 2881 | Put the exponent, in the first byte, in d2, to use the "standard" rounding 2882 | routines: 2883 movel d6,d2 2884 #ifndef __mcoldfire__ 2885 lsrw IMM (8),d2 2886 #else 2887 lsrl IMM (8),d2 2888 #endif 2889 2890 | Before rounding normalize so bit #FLT_MANT_DIG is set (we will consider 2891 | the case of denormalized numbers in the rounding routine itself). 2892 | As in the addition (not in the subtraction!) we could have set 2893 | one more bit we check this: 2894 btst IMM (FLT_MANT_DIG+1),d0 2895 beq 1f 2896 #ifndef __mcoldfire__ 2897 lsrl IMM (1),d0 2898 roxrl IMM (1),d1 2899 #else 2900 lsrl IMM (1),d1 2901 btst IMM (0),d0 2902 beq 10f 2903 bset IMM (31),d1 2904 10: lsrl IMM (1),d0 2905 #endif 2906 addl IMM (1),d2 2907 1: 2908 lea pc@(Laddsf$4),a0 | to return from rounding routine 2909 PICLEA SYM (_fpCCR),a1 | check the rounding mode 2910 #ifdef __mcoldfire__ 2911 clrl d6 2912 #endif 2913 movew a1@(6),d6 | rounding mode in d6 2914 beq Lround$to$nearest 2915 #ifndef __mcoldfire__ 2916 cmpw IMM (ROUND_TO_PLUS),d6 2917 #else 2918 cmpl IMM (ROUND_TO_PLUS),d6 2919 #endif 2920 bhi Lround$to$minus 2921 blt Lround$to$zero 2922 bra Lround$to$plus 2923 Laddsf$4: 2924 | Put back the exponent, but check for overflow. 2925 #ifndef __mcoldfire__ 2926 cmpw IMM (0xff),d2 2927 #else 2928 cmpl IMM (0xff),d2 2929 #endif 2930 bhi 1f 2931 bclr IMM (FLT_MANT_DIG-1),d0 2932 #ifndef __mcoldfire__ 2933 lslw IMM (7),d2 2934 #else 2935 lsll IMM (7),d2 2936 #endif 2937 swap d2 2938 orl d2,d0 2939 bra Laddsf$ret 2940 1: 2941 moveq IMM (ADD),d5 2942 bra Lf$overflow 2943 2944 Lsubsf$0: 2945 | We are here if a > 0 and b < 0 (sign bits cleared). 2946 | Here we do the subtraction. 2947 movel d6,d7 | put sign in d7 2948 andl IMM (0x80000000),d7 2949 2950 subl d3,d1 | result in d0-d1 2951 subxl d2,d0 | 2952 beq Laddsf$ret | if zero just exit 2953 bpl 1f | if positive skip the following 2954 bchg IMM (31),d7 | change sign bit in d7 2955 negl d1 2956 negxl d0 2957 1: 2958 #ifndef __mcoldfire__ 2959 exg d2,a0 | now we have the exponent in d2 2960 lsrw IMM (8),d2 | put it in the first byte 2961 #else 2962 movel d2,d4 2963 movel a0,d2 2964 movel d4,a0 2965 lsrl IMM (8),d2 | put it in the first byte 2966 #endif 2967 2968 | Now d0-d1 is positive and the sign bit is in d7. 2969 2970 | Note that we do not have to normalize, since in the subtraction bit 2971 | #FLT_MANT_DIG+1 is never set, and denormalized numbers are handled by 2972 | the rounding routines themselves. 2973 lea pc@(Lsubsf$1),a0 | to return from rounding routine 2974 PICLEA SYM (_fpCCR),a1 | check the rounding mode 2975 #ifdef __mcoldfire__ 2976 clrl d6 2977 #endif 2978 movew a1@(6),d6 | rounding mode in d6 2979 beq Lround$to$nearest 2980 #ifndef __mcoldfire__ 2981 cmpw IMM (ROUND_TO_PLUS),d6 2982 #else 2983 cmpl IMM (ROUND_TO_PLUS),d6 2984 #endif 2985 bhi Lround$to$minus 2986 blt Lround$to$zero 2987 bra Lround$to$plus 2988 Lsubsf$1: 2989 | Put back the exponent (we can't have overflow!). ' 2990 bclr IMM (FLT_MANT_DIG-1),d0 2991 #ifndef __mcoldfire__ 2992 lslw IMM (7),d2 2993 #else 2994 lsll IMM (7),d2 2995 #endif 2996 swap d2 2997 orl d2,d0 2998 bra Laddsf$ret 2999 3000 | If one of the numbers was too small (difference of exponents >= 3001 | FLT_MANT_DIG+2) we return the other (and now we don't have to ' 3002 | check for finiteness or zero). 3003 Laddsf$a$small: 3004 movel a6@(12),d0 3005 PICLEA SYM (_fpCCR),a0 3006 movew IMM (0),a0@ 3007 #ifndef __mcoldfire__ 3008 moveml sp@+,d2-d7 | restore data registers 3009 #else 3010 moveml sp@,d2-d7 3011 | XXX if frame pointer is ever removed, stack pointer must 3012 | be adjusted here. 3013 #endif 3014 unlk a6 | and return 3015 rts 3016 3017 Laddsf$b$small: 3018 movel a6@(8),d0 3019 PICLEA SYM (_fpCCR),a0 3020 movew IMM (0),a0@ 3021 #ifndef __mcoldfire__ 3022 moveml sp@+,d2-d7 | restore data registers 3023 #else 3024 moveml sp@,d2-d7 3025 | XXX if frame pointer is ever removed, stack pointer must 3026 | be adjusted here. 3027 #endif 3028 unlk a6 | and return 3029 rts 3030 3031 | If the numbers are denormalized remember to put exponent equal to 1. 3032 3033 Laddsf$a$den: 3034 movel d5,d6 | d5 contains 0x01000000 3035 swap d6 3036 bra Laddsf$1 3037 3038 Laddsf$b$den: 3039 movel d5,d7 3040 swap d7 3041 notl d4 | make d4 into a mask for the fraction 3042 | (this was not executed after the jump) 3043 bra Laddsf$2 3044 3045 | The rest is mainly code for the different results which can be 3046 | returned (checking always for +/-INFINITY and NaN). 3047 3048 Laddsf$b: 3049 | Return b (if a is zero). 3050 movel a6@(12),d0 3051 cmpl IMM (0x80000000),d0 | Check if b is -0 3052 bne 1f 3053 movel a0,d7 3054 andl IMM (0x80000000),d7 | Use the sign of a 3055 clrl d0 3056 bra Laddsf$ret 3057 Laddsf$a: 3058 | Return a (if b is zero). 3059 movel a6@(8),d0 3060 1: 3061 moveq IMM (ADD),d5 3062 | We have to check for NaN and +/-infty. 3063 movel d0,d7 3064 andl IMM (0x80000000),d7 | put sign in d7 3065 bclr IMM (31),d0 | clear sign 3066 cmpl IMM (INFINITY),d0 | check for infty or NaN 3067 bge 2f 3068 movel d0,d0 | check for zero (we do this because we don't ' 3069 bne Laddsf$ret | want to return -0 by mistake 3070 bclr IMM (31),d7 | if zero be sure to clear sign 3071 bra Laddsf$ret | if everything OK just return 3072 2: 3073 | The value to be returned is either +/-infty or NaN 3074 andl IMM (0x007fffff),d0 | check for NaN 3075 bne Lf$inop | if mantissa not zero is NaN 3076 bra Lf$infty 3077 3078 Laddsf$ret: 3079 | Normal exit (a and b nonzero, result is not NaN nor +/-infty). 3080 | We have to clear the exception flags (just the exception type). 3081 PICLEA SYM (_fpCCR),a0 3082 movew IMM (0),a0@ 3083 orl d7,d0 | put sign bit 3084 #ifndef __mcoldfire__ 3085 moveml sp@+,d2-d7 | restore data registers 3086 #else 3087 moveml sp@,d2-d7 3088 | XXX if frame pointer is ever removed, stack pointer must 3089 | be adjusted here. 3090 #endif 3091 unlk a6 | and return 3092 rts 3093 3094 Laddsf$ret$den: 3095 | Return a denormalized number (for addition we don't signal underflow) ' 3096 lsrl IMM (1),d0 | remember to shift right back once 3097 bra Laddsf$ret | and return 3098 3099 | Note: when adding two floats of the same sign if either one is 3100 | NaN we return NaN without regard to whether the other is finite or 3101 | not. When subtracting them (i.e., when adding two numbers of 3102 | opposite signs) things are more complicated: if both are INFINITY 3103 | we return NaN, if only one is INFINITY and the other is NaN we return 3104 | NaN, but if it is finite we return INFINITY with the corresponding sign. 3105 3106 Laddsf$nf: 3107 moveq IMM (ADD),d5 3108 | This could be faster but it is not worth the effort, since it is not 3109 | executed very often. We sacrifice speed for clarity here. 3110 movel a6@(8),d0 | get the numbers back (remember that we 3111 movel a6@(12),d1 | did some processing already) 3112 movel IMM (INFINITY),d4 | useful constant (INFINITY) 3113 movel d0,d2 | save sign bits 3114 movel d1,d3 3115 bclr IMM (31),d0 | clear sign bits 3116 bclr IMM (31),d1 3117 | We know that one of them is either NaN of +/-INFINITY 3118 | Check for NaN (if either one is NaN return NaN) 3119 cmpl d4,d0 | check first a (d0) 3120 bhi Lf$inop 3121 cmpl d4,d1 | check now b (d1) 3122 bhi Lf$inop 3123 | Now comes the check for +/-INFINITY. We know that both are (maybe not 3124 | finite) numbers, but we have to check if both are infinite whether we 3125 | are adding or subtracting them. 3126 eorl d3,d2 | to check sign bits 3127 bmi 1f 3128 movel d0,d7 3129 andl IMM (0x80000000),d7 | get (common) sign bit 3130 bra Lf$infty 3131 1: 3132 | We know one (or both) are infinite, so we test for equality between the 3133 | two numbers (if they are equal they have to be infinite both, so we 3134 | return NaN). 3135 cmpl d1,d0 | are both infinite? 3136 beq Lf$inop | if so return NaN 3137 3138 movel d0,d7 3139 andl IMM (0x80000000),d7 | get a's sign bit ' 3140 cmpl d4,d0 | test now for infinity 3141 beq Lf$infty | if a is INFINITY return with this sign 3142 bchg IMM (31),d7 | else we know b is INFINITY and has 3143 bra Lf$infty | the opposite sign 3144 3145 |============================================================================= 3146 | __mulsf3 3147 |============================================================================= 3148 3149 | float __mulsf3(float, float); 3150 FUNC(__mulsf3) 3151 SYM (__mulsf3): 3152 #ifndef __mcoldfire__ 3153 link a6,IMM (0) 3154 moveml d2-d7,sp@- 3155 #else 3156 link a6,IMM (-24) 3157 moveml d2-d7,sp@ 3158 #endif 3159 movel a6@(8),d0 | get a into d0 3160 movel a6@(12),d1 | and b into d1 3161 movel d0,d7 | d7 will hold the sign of the product 3162 eorl d1,d7 | 3163 andl IMM (0x80000000),d7 3164 movel IMM (INFINITY),d6 | useful constant (+INFINITY) 3165 movel d6,d5 | another (mask for fraction) 3166 notl d5 | 3167 movel IMM (0x00800000),d4 | this is to put hidden bit back 3168 bclr IMM (31),d0 | get rid of a's sign bit ' 3169 movel d0,d2 | 3170 beq Lmulsf$a$0 | branch if a is zero 3171 bclr IMM (31),d1 | get rid of b's sign bit ' 3172 movel d1,d3 | 3173 beq Lmulsf$b$0 | branch if b is zero 3174 cmpl d6,d0 | is a big? 3175 bhi Lmulsf$inop | if a is NaN return NaN 3176 beq Lmulsf$inf | if a is INFINITY we have to check b 3177 cmpl d6,d1 | now compare b with INFINITY 3178 bhi Lmulsf$inop | is b NaN? 3179 beq Lmulsf$overflow | is b INFINITY? 3180 | Here we have both numbers finite and nonzero (and with no sign bit). 3181 | Now we get the exponents into d2 and d3. 3182 andl d6,d2 | and isolate exponent in d2 3183 beq Lmulsf$a$den | if exponent is zero we have a denormalized 3184 andl d5,d0 | and isolate fraction 3185 orl d4,d0 | and put hidden bit back 3186 swap d2 | I like exponents in the first byte 3187 #ifndef __mcoldfire__ 3188 lsrw IMM (7),d2 | 3189 #else 3190 lsrl IMM (7),d2 | 3191 #endif 3192 Lmulsf$1: | number 3193 andl d6,d3 | 3194 beq Lmulsf$b$den | 3195 andl d5,d1 | 3196 orl d4,d1 | 3197 swap d3 | 3198 #ifndef __mcoldfire__ 3199 lsrw IMM (7),d3 | 3200 #else 3201 lsrl IMM (7),d3 | 3202 #endif 3203 Lmulsf$2: | 3204 #ifndef __mcoldfire__ 3205 addw d3,d2 | add exponents 3206 subw IMM (F_BIAS+1),d2 | and subtract bias (plus one) 3207 #else 3208 addl d3,d2 | add exponents 3209 subl IMM (F_BIAS+1),d2 | and subtract bias (plus one) 3210 #endif 3211 3212 | We are now ready to do the multiplication. The situation is as follows: 3213 | both a and b have bit FLT_MANT_DIG-1 set (even if they were 3214 | denormalized to start with!), which means that in the product 3215 | bit 2*(FLT_MANT_DIG-1) (that is, bit 2*FLT_MANT_DIG-2-32 of the 3216 | high long) is set. 3217 3218 | To do the multiplication let us move the number a little bit around ... 3219 movel d1,d6 | second operand in d6 3220 movel d0,d5 | first operand in d4-d5 3221 movel IMM (0),d4 3222 movel d4,d1 | the sums will go in d0-d1 3223 movel d4,d0 3224 3225 | now bit FLT_MANT_DIG-1 becomes bit 31: 3226 lsll IMM (31-FLT_MANT_DIG+1),d6 3227 3228 | Start the loop (we loop #FLT_MANT_DIG times): 3229 moveq IMM (FLT_MANT_DIG-1),d3 3230 1: addl d1,d1 | shift sum 3231 addxl d0,d0 3232 lsll IMM (1),d6 | get bit bn 3233 bcc 2f | if not set skip sum 3234 addl d5,d1 | add a 3235 addxl d4,d0 3236 2: 3237 #ifndef __mcoldfire__ 3238 dbf d3,1b | loop back 3239 #else 3240 subql IMM (1),d3 3241 bpl 1b 3242 #endif 3243 3244 | Now we have the product in d0-d1, with bit (FLT_MANT_DIG - 1) + FLT_MANT_DIG 3245 | (mod 32) of d0 set. The first thing to do now is to normalize it so bit 3246 | FLT_MANT_DIG is set (to do the rounding). 3247 #ifndef __mcoldfire__ 3248 rorl IMM (6),d1 3249 swap d1 3250 movew d1,d3 3251 andw IMM (0x03ff),d3 3252 andw IMM (0xfd00),d1 3253 #else 3254 movel d1,d3 3255 lsll IMM (8),d1 3256 addl d1,d1 3257 addl d1,d1 3258 moveq IMM (22),d5 3259 lsrl d5,d3 3260 orl d3,d1 3261 andl IMM (0xfffffd00),d1 3262 #endif 3263 lsll IMM (8),d0 3264 addl d0,d0 3265 addl d0,d0 3266 #ifndef __mcoldfire__ 3267 orw d3,d0 3268 #else 3269 orl d3,d0 3270 #endif 3271 3272 moveq IMM (MULTIPLY),d5 3273 3274 btst IMM (FLT_MANT_DIG+1),d0 3275 beq Lround$exit 3276 #ifndef __mcoldfire__ 3277 lsrl IMM (1),d0 3278 roxrl IMM (1),d1 3279 addw IMM (1),d2 3280 #else 3281 lsrl IMM (1),d1 3282 btst IMM (0),d0 3283 beq 10f 3284 bset IMM (31),d1 3285 10: lsrl IMM (1),d0 3286 addql IMM (1),d2 3287 #endif 3288 bra Lround$exit 3289 3290 Lmulsf$inop: 3291 moveq IMM (MULTIPLY),d5 3292 bra Lf$inop 3293 3294 Lmulsf$overflow: 3295 moveq IMM (MULTIPLY),d5 3296 bra Lf$overflow 3297 3298 Lmulsf$inf: 3299 moveq IMM (MULTIPLY),d5 3300 | If either is NaN return NaN; else both are (maybe infinite) numbers, so 3301 | return INFINITY with the correct sign (which is in d7). 3302 cmpl d6,d1 | is b NaN? 3303 bhi Lf$inop | if so return NaN 3304 bra Lf$overflow | else return +/-INFINITY 3305 3306 | If either number is zero return zero, unless the other is +/-INFINITY, 3307 | or NaN, in which case we return NaN. 3308 Lmulsf$b$0: 3309 | Here d1 (==b) is zero. 3310 movel a6@(8),d1 | get a again to check for non-finiteness 3311 bra 1f 3312 Lmulsf$a$0: 3313 movel a6@(12),d1 | get b again to check for non-finiteness 3314 1: bclr IMM (31),d1 | clear sign bit 3315 cmpl IMM (INFINITY),d1 | and check for a large exponent 3316 bge Lf$inop | if b is +/-INFINITY or NaN return NaN 3317 movel d7,d0 | else return signed zero 3318 PICLEA SYM (_fpCCR),a0 | 3319 movew IMM (0),a0@ | 3320 #ifndef __mcoldfire__ 3321 moveml sp@+,d2-d7 | 3322 #else 3323 moveml sp@,d2-d7 3324 | XXX if frame pointer is ever removed, stack pointer must 3325 | be adjusted here. 3326 #endif 3327 unlk a6 | 3328 rts | 3329 3330 | If a number is denormalized we put an exponent of 1 but do not put the 3331 | hidden bit back into the fraction; instead we shift left until bit 23 3332 | (the hidden bit) is set, adjusting the exponent accordingly. We do this 3333 | to ensure that the product of the fractions is close to 1. 3334 Lmulsf$a$den: 3335 movel IMM (1),d2 3336 andl d5,d0 3337 1: addl d0,d0 | shift a left (until bit 23 is set) 3338 #ifndef __mcoldfire__ 3339 subw IMM (1),d2 | and adjust exponent 3340 #else 3341 subql IMM (1),d2 | and adjust exponent 3342 #endif 3343 btst IMM (FLT_MANT_DIG-1),d0 3344 bne Lmulsf$1 | 3345 bra 1b | else loop back 3346 3347 Lmulsf$b$den: 3348 movel IMM (1),d3 3349 andl d5,d1 3350 1: addl d1,d1 | shift b left until bit 23 is set 3351 #ifndef __mcoldfire__ 3352 subw IMM (1),d3 | and adjust exponent 3353 #else 3354 subql IMM (1),d3 | and adjust exponent 3355 #endif 3356 btst IMM (FLT_MANT_DIG-1),d1 3357 bne Lmulsf$2 | 3358 bra 1b | else loop back 3359 3360 |============================================================================= 3361 | __divsf3 3362 |============================================================================= 3363 3364 | float __divsf3(float, float); 3365 FUNC(__divsf3) 3366 SYM (__divsf3): 3367 #ifndef __mcoldfire__ 3368 link a6,IMM (0) 3369 moveml d2-d7,sp@- 3370 #else 3371 link a6,IMM (-24) 3372 moveml d2-d7,sp@ 3373 #endif 3374 movel a6@(8),d0 | get a into d0 3375 movel a6@(12),d1 | and b into d1 3376 movel d0,d7 | d7 will hold the sign of the result 3377 eorl d1,d7 | 3378 andl IMM (0x80000000),d7 | 3379 movel IMM (INFINITY),d6 | useful constant (+INFINITY) 3380 movel d6,d5 | another (mask for fraction) 3381 notl d5 | 3382 movel IMM (0x00800000),d4 | this is to put hidden bit back 3383 bclr IMM (31),d0 | get rid of a's sign bit ' 3384 movel d0,d2 | 3385 beq Ldivsf$a$0 | branch if a is zero 3386 bclr IMM (31),d1 | get rid of b's sign bit ' 3387 movel d1,d3 | 3388 beq Ldivsf$b$0 | branch if b is zero 3389 cmpl d6,d0 | is a big? 3390 bhi Ldivsf$inop | if a is NaN return NaN 3391 beq Ldivsf$inf | if a is INFINITY we have to check b 3392 cmpl d6,d1 | now compare b with INFINITY 3393 bhi Ldivsf$inop | if b is NaN return NaN 3394 beq Ldivsf$underflow 3395 | Here we have both numbers finite and nonzero (and with no sign bit). 3396 | Now we get the exponents into d2 and d3 and normalize the numbers to 3397 | ensure that the ratio of the fractions is close to 1. We do this by 3398 | making sure that bit #FLT_MANT_DIG-1 (hidden bit) is set. 3399 andl d6,d2 | and isolate exponent in d2 3400 beq Ldivsf$a$den | if exponent is zero we have a denormalized 3401 andl d5,d0 | and isolate fraction 3402 orl d4,d0 | and put hidden bit back 3403 swap d2 | I like exponents in the first byte 3404 #ifndef __mcoldfire__ 3405 lsrw IMM (7),d2 | 3406 #else 3407 lsrl IMM (7),d2 | 3408 #endif 3409 Ldivsf$1: | 3410 andl d6,d3 | 3411 beq Ldivsf$b$den | 3412 andl d5,d1 | 3413 orl d4,d1 | 3414 swap d3 | 3415 #ifndef __mcoldfire__ 3416 lsrw IMM (7),d3 | 3417 #else 3418 lsrl IMM (7),d3 | 3419 #endif 3420 Ldivsf$2: | 3421 #ifndef __mcoldfire__ 3422 subw d3,d2 | subtract exponents 3423 addw IMM (F_BIAS),d2 | and add bias 3424 #else 3425 subl d3,d2 | subtract exponents 3426 addl IMM (F_BIAS),d2 | and add bias 3427 #endif 3428 3429 | We are now ready to do the division. We have prepared things in such a way 3430 | that the ratio of the fractions will be less than 2 but greater than 1/2. 3431 | At this point the registers in use are: 3432 | d0 holds a (first operand, bit FLT_MANT_DIG=0, bit FLT_MANT_DIG-1=1) 3433 | d1 holds b (second operand, bit FLT_MANT_DIG=1) 3434 | d2 holds the difference of the exponents, corrected by the bias 3435 | d7 holds the sign of the ratio 3436 | d4, d5, d6 hold some constants 3437 movel d7,a0 | d6-d7 will hold the ratio of the fractions 3438 movel IMM (0),d6 | 3439 movel d6,d7 3440 3441 moveq IMM (FLT_MANT_DIG+1),d3 3442 1: cmpl d0,d1 | is a < b? 3443 bhi 2f | 3444 bset d3,d6 | set a bit in d6 3445 subl d1,d0 | if a >= b a <-- a-b 3446 beq 3f | if a is zero, exit 3447 2: addl d0,d0 | multiply a by 2 3448 #ifndef __mcoldfire__ 3449 dbra d3,1b 3450 #else 3451 subql IMM (1),d3 3452 bpl 1b 3453 #endif 3454 3455 | Now we keep going to set the sticky bit ... 3456 moveq IMM (FLT_MANT_DIG),d3 3457 1: cmpl d0,d1 3458 ble 2f 3459 addl d0,d0 3460 #ifndef __mcoldfire__ 3461 dbra d3,1b 3462 #else 3463 subql IMM(1),d3 3464 bpl 1b 3465 #endif 3466 movel IMM (0),d1 3467 bra 3f 3468 2: movel IMM (0),d1 3469 #ifndef __mcoldfire__ 3470 subw IMM (FLT_MANT_DIG),d3 3471 addw IMM (31),d3 3472 #else 3473 subl IMM (FLT_MANT_DIG),d3 3474 addl IMM (31),d3 3475 #endif 3476 bset d3,d1 3477 3: 3478 movel d6,d0 | put the ratio in d0-d1 3479 movel a0,d7 | get sign back 3480 3481 | Because of the normalization we did before we are guaranteed that 3482 | d0 is smaller than 2^26 but larger than 2^24. Thus bit 26 is not set, 3483 | bit 25 could be set, and if it is not set then bit 24 is necessarily set. 3484 btst IMM (FLT_MANT_DIG+1),d0 3485 beq 1f | if it is not set, then bit 24 is set 3486 lsrl IMM (1),d0 | 3487 #ifndef __mcoldfire__ 3488 addw IMM (1),d2 | 3489 #else 3490 addl IMM (1),d2 | 3491 #endif 3492 1: 3493 | Now round, check for over- and underflow, and exit. 3494 moveq IMM (DIVIDE),d5 3495 bra Lround$exit 3496 3497 Ldivsf$inop: 3498 moveq IMM (DIVIDE),d5 3499 bra Lf$inop 3500 3501 Ldivsf$overflow: 3502 moveq IMM (DIVIDE),d5 3503 bra Lf$overflow 3504 3505 Ldivsf$underflow: 3506 moveq IMM (DIVIDE),d5 3507 bra Lf$underflow 3508 3509 Ldivsf$a$0: 3510 moveq IMM (DIVIDE),d5 3511 | If a is zero check to see whether b is zero also. In that case return 3512 | NaN; then check if b is NaN, and return NaN also in that case. Else 3513 | return a properly signed zero. 3514 andl IMM (0x7fffffff),d1 | clear sign bit and test b 3515 beq Lf$inop | if b is also zero return NaN 3516 cmpl IMM (INFINITY),d1 | check for NaN 3517 bhi Lf$inop | 3518 movel d7,d0 | else return signed zero 3519 PICLEA SYM (_fpCCR),a0 | 3520 movew IMM (0),a0@ | 3521 #ifndef __mcoldfire__ 3522 moveml sp@+,d2-d7 | 3523 #else 3524 moveml sp@,d2-d7 | 3525 | XXX if frame pointer is ever removed, stack pointer must 3526 | be adjusted here. 3527 #endif 3528 unlk a6 | 3529 rts | 3530 3531 Ldivsf$b$0: 3532 moveq IMM (DIVIDE),d5 3533 | If we got here a is not zero. Check if a is NaN; in that case return NaN, 3534 | else return +/-INFINITY. Remember that a is in d0 with the sign bit 3535 | cleared already. 3536 cmpl IMM (INFINITY),d0 | compare d0 with INFINITY 3537 bhi Lf$inop | if larger it is NaN 3538 bra Lf$div$0 | else signal DIVIDE_BY_ZERO 3539 3540 Ldivsf$inf: 3541 moveq IMM (DIVIDE),d5 3542 | If a is INFINITY we have to check b 3543 cmpl IMM (INFINITY),d1 | compare b with INFINITY 3544 bge Lf$inop | if b is NaN or INFINITY return NaN 3545 bra Lf$overflow | else return overflow 3546 3547 | If a number is denormalized we put an exponent of 1 but do not put the 3548 | bit back into the fraction. 3549 Ldivsf$a$den: 3550 movel IMM (1),d2 3551 andl d5,d0 3552 1: addl d0,d0 | shift a left until bit FLT_MANT_DIG-1 is set 3553 #ifndef __mcoldfire__ 3554 subw IMM (1),d2 | and adjust exponent 3555 #else 3556 subl IMM (1),d2 | and adjust exponent 3557 #endif 3558 btst IMM (FLT_MANT_DIG-1),d0 3559 bne Ldivsf$1 3560 bra 1b 3561 3562 Ldivsf$b$den: 3563 movel IMM (1),d3 3564 andl d5,d1 3565 1: addl d1,d1 | shift b left until bit FLT_MANT_DIG is set 3566 #ifndef __mcoldfire__ 3567 subw IMM (1),d3 | and adjust exponent 3568 #else 3569 subl IMM (1),d3 | and adjust exponent 3570 #endif 3571 btst IMM (FLT_MANT_DIG-1),d1 3572 bne Ldivsf$2 3573 bra 1b 3574 3575 Lround$exit: 3576 | This is a common exit point for __mulsf3 and __divsf3. 3577 3578 | First check for underlow in the exponent: 3579 #ifndef __mcoldfire__ 3580 cmpw IMM (-FLT_MANT_DIG-1),d2 3581 #else 3582 cmpl IMM (-FLT_MANT_DIG-1),d2 3583 #endif 3584 blt Lf$underflow 3585 | It could happen that the exponent is less than 1, in which case the 3586 | number is denormalized. In this case we shift right and adjust the 3587 | exponent until it becomes 1 or the fraction is zero (in the latter case 3588 | we signal underflow and return zero). 3589 movel IMM (0),d6 | d6 is used temporarily 3590 #ifndef __mcoldfire__ 3591 cmpw IMM (1),d2 | if the exponent is less than 1 we 3592 #else 3593 cmpl IMM (1),d2 | if the exponent is less than 1 we 3594 #endif 3595 bge 2f | have to shift right (denormalize) 3596 1: 3597 #ifndef __mcoldfire__ 3598 addw IMM (1),d2 | adjust the exponent 3599 lsrl IMM (1),d0 | shift right once 3600 roxrl IMM (1),d1 | 3601 roxrl IMM (1),d6 | d6 collect bits we would lose otherwise 3602 cmpw IMM (1),d2 | is the exponent 1 already? 3603 #else 3604 addql IMM (1),d2 | adjust the exponent 3605 lsrl IMM (1),d6 3606 btst IMM (0),d1 3607 beq 11f 3608 bset IMM (31),d6 3609 11: lsrl IMM (1),d1 3610 btst IMM (0),d0 3611 beq 10f 3612 bset IMM (31),d1 3613 10: lsrl IMM (1),d0 3614 cmpl IMM (1),d2 | is the exponent 1 already? 3615 #endif 3616 beq 2f | if not loop back 3617 bra 1b | 3618 bra Lf$underflow | safety check, shouldn't execute ' 3619 2: orl d6,d1 | this is a trick so we don't lose ' 3620 | the extra bits which were flushed right 3621 | Now call the rounding routine (which takes care of denormalized numbers): 3622 lea pc@(Lround$0),a0 | to return from rounding routine 3623 PICLEA SYM (_fpCCR),a1 | check the rounding mode 3624 #ifdef __mcoldfire__ 3625 clrl d6 3626 #endif 3627 movew a1@(6),d6 | rounding mode in d6 3628 beq Lround$to$nearest 3629 #ifndef __mcoldfire__ 3630 cmpw IMM (ROUND_TO_PLUS),d6 3631 #else 3632 cmpl IMM (ROUND_TO_PLUS),d6 3633 #endif 3634 bhi Lround$to$minus 3635 blt Lround$to$zero 3636 bra Lround$to$plus 3637 Lround$0: 3638 | Here we have a correctly rounded result (either normalized or denormalized). 3639 3640 | Here we should have either a normalized number or a denormalized one, and 3641 | the exponent is necessarily larger or equal to 1 (so we don't have to ' 3642 | check again for underflow!). We have to check for overflow or for a 3643 | denormalized number (which also signals underflow). 3644 | Check for overflow (i.e., exponent >= 255). 3645 #ifndef __mcoldfire__ 3646 cmpw IMM (0x00ff),d2 3647 #else 3648 cmpl IMM (0x00ff),d2 3649 #endif 3650 bge Lf$overflow 3651 | Now check for a denormalized number (exponent==0). 3652 movew d2,d2 3653 beq Lf$den 3654 1: 3655 | Put back the exponents and sign and return. 3656 #ifndef __mcoldfire__ 3657 lslw IMM (7),d2 | exponent back to fourth byte 3658 #else 3659 lsll IMM (7),d2 | exponent back to fourth byte 3660 #endif 3661 bclr IMM (FLT_MANT_DIG-1),d0 3662 swap d0 | and put back exponent 3663 #ifndef __mcoldfire__ 3664 orw d2,d0 | 3665 #else 3666 orl d2,d0 3667 #endif 3668 swap d0 | 3669 orl d7,d0 | and sign also 3670 3671 PICLEA SYM (_fpCCR),a0 3672 movew IMM (0),a0@ 3673 #ifndef __mcoldfire__ 3674 moveml sp@+,d2-d7 3675 #else 3676 moveml sp@,d2-d7 3677 | XXX if frame pointer is ever removed, stack pointer must 3678 | be adjusted here. 3679 #endif 3680 unlk a6 3681 rts 3682 3683 |============================================================================= 3684 | __negsf2 3685 |============================================================================= 3686 3687 | This is trivial and could be shorter if we didn't bother checking for NaN ' 3688 | and +/-INFINITY. 3689 3690 | float __negsf2(float); 3691 FUNC(__negsf2) 3692 SYM (__negsf2): 3693 #ifndef __mcoldfire__ 3694 link a6,IMM (0) 3695 moveml d2-d7,sp@- 3696 #else 3697 link a6,IMM (-24) 3698 moveml d2-d7,sp@ 3699 #endif 3700 moveq IMM (NEGATE),d5 3701 movel a6@(8),d0 | get number to negate in d0 3702 bchg IMM (31),d0 | negate 3703 movel d0,d1 | make a positive copy 3704 bclr IMM (31),d1 | 3705 tstl d1 | check for zero 3706 beq 2f | if zero (either sign) return +zero 3707 cmpl IMM (INFINITY),d1 | compare to +INFINITY 3708 blt 1f | 3709 bhi Lf$inop | if larger (fraction not zero) is NaN 3710 movel d0,d7 | else get sign and return INFINITY 3711 andl IMM (0x80000000),d7 3712 bra Lf$infty 3713 1: PICLEA SYM (_fpCCR),a0 3714 movew IMM (0),a0@ 3715 #ifndef __mcoldfire__ 3716 moveml sp@+,d2-d7 3717 #else 3718 moveml sp@,d2-d7 3719 | XXX if frame pointer is ever removed, stack pointer must 3720 | be adjusted here. 3721 #endif 3722 unlk a6 3723 rts 3724 2: bclr IMM (31),d0 3725 bra 1b 3726 3727 |============================================================================= 3728 | __cmpsf2 3729 |============================================================================= 3730 3731 GREATER = 1 3732 LESS = -1 3733 EQUAL = 0 3734 3735 | int __cmpsf2_internal(float, float, int); 3736 SYM (__cmpsf2_internal): 3737 #ifndef __mcoldfire__ 3738 link a6,IMM (0) 3739 moveml d2-d7,sp@- | save registers 3740 #else 3741 link a6,IMM (-24) 3742 moveml d2-d7,sp@ 3743 #endif 3744 moveq IMM (COMPARE),d5 3745 movel a6@(8),d0 | get first operand 3746 movel a6@(12),d1 | get second operand 3747 | Check if either is NaN, and in that case return garbage and signal 3748 | INVALID_OPERATION. Check also if either is zero, and clear the signs 3749 | if necessary. 3750 movel d0,d6 3751 andl IMM (0x7fffffff),d0 3752 beq Lcmpsf$a$0 3753 cmpl IMM (0x7f800000),d0 3754 bhi Lcmpf$inop 3755 Lcmpsf$1: 3756 movel d1,d7 3757 andl IMM (0x7fffffff),d1 3758 beq Lcmpsf$b$0 3759 cmpl IMM (0x7f800000),d1 3760 bhi Lcmpf$inop 3761 Lcmpsf$2: 3762 | Check the signs 3763 eorl d6,d7 3764 bpl 1f 3765 | If the signs are not equal check if a >= 0 3766 tstl d6 3767 bpl Lcmpsf$a$gt$b | if (a >= 0 && b < 0) => a > b 3768 bmi Lcmpsf$b$gt$a | if (a < 0 && b >= 0) => a < b 3769 1: 3770 | If the signs are equal check for < 0 3771 tstl d6 3772 bpl 1f 3773 | If both are negative exchange them 3774 #ifndef __mcoldfire__ 3775 exg d0,d1 3776 #else 3777 movel d0,d7 3778 movel d1,d0 3779 movel d7,d1 3780 #endif 3781 1: 3782 | Now that they are positive we just compare them as longs (does this also 3783 | work for denormalized numbers?). 3784 cmpl d0,d1 3785 bhi Lcmpsf$b$gt$a | |b| > |a| 3786 bne Lcmpsf$a$gt$b | |b| < |a| 3787 | If we got here a == b. 3788 movel IMM (EQUAL),d0 3789 #ifndef __mcoldfire__ 3790 moveml sp@+,d2-d7 | put back the registers 3791 #else 3792 moveml sp@,d2-d7 3793 #endif 3794 unlk a6 3795 rts 3796 Lcmpsf$a$gt$b: 3797 movel IMM (GREATER),d0 3798 #ifndef __mcoldfire__ 3799 moveml sp@+,d2-d7 | put back the registers 3800 #else 3801 moveml sp@,d2-d7 3802 | XXX if frame pointer is ever removed, stack pointer must 3803 | be adjusted here. 3804 #endif 3805 unlk a6 3806 rts 3807 Lcmpsf$b$gt$a: 3808 movel IMM (LESS),d0 3809 #ifndef __mcoldfire__ 3810 moveml sp@+,d2-d7 | put back the registers 3811 #else 3812 moveml sp@,d2-d7 3813 | XXX if frame pointer is ever removed, stack pointer must 3814 | be adjusted here. 3815 #endif 3816 unlk a6 3817 rts 3818 3819 Lcmpsf$a$0: 3820 bclr IMM (31),d6 3821 bra Lcmpsf$1 3822 Lcmpsf$b$0: 3823 bclr IMM (31),d7 3824 bra Lcmpsf$2 3825 3826 Lcmpf$inop: 3827 movl a6@(16),d0 3828 moveq IMM (INEXACT_RESULT+INVALID_OPERATION),d7 3829 moveq IMM (SINGLE_FLOAT),d6 3830 PICJUMP $_exception_handler 3831 3832 | int __cmpsf2(float, float); 3833 FUNC(__cmpsf2) 3834 SYM (__cmpsf2): 3835 link a6,IMM (0) 3836 pea 1 3837 movl a6@(12),sp@- 3838 movl a6@(8),sp@- 3839 PICCALL SYM (__cmpsf2_internal) 3840 unlk a6 3841 rts 3842 3843 |============================================================================= 3844 | rounding routines 3845 |============================================================================= 3846 3847 | The rounding routines expect the number to be normalized in registers 3848 | d0-d1, with the exponent in register d2. They assume that the 3849 | exponent is larger or equal to 1. They return a properly normalized number 3850 | if possible, and a denormalized number otherwise. The exponent is returned 3851 | in d2. 3852 3853 Lround$to$nearest: 3854 | We now normalize as suggested by D. Knuth ("Seminumerical Algorithms"): 3855 | Here we assume that the exponent is not too small (this should be checked 3856 | before entering the rounding routine), but the number could be denormalized. 3857 3858 | Check for denormalized numbers: 3859 1: btst IMM (FLT_MANT_DIG),d0 3860 bne 2f | if set the number is normalized 3861 | Normalize shifting left until bit #FLT_MANT_DIG is set or the exponent 3862 | is one (remember that a denormalized number corresponds to an 3863 | exponent of -F_BIAS+1). 3864 #ifndef __mcoldfire__ 3865 cmpw IMM (1),d2 | remember that the exponent is at least one 3866 #else 3867 cmpl IMM (1),d2 | remember that the exponent is at least one 3868 #endif 3869 beq 2f | an exponent of one means denormalized 3870 addl d1,d1 | else shift and adjust the exponent 3871 addxl d0,d0 | 3872 #ifndef __mcoldfire__ 3873 dbra d2,1b | 3874 #else 3875 subql IMM (1),d2 3876 bpl 1b 3877 #endif 3878 2: 3879 | Now round: we do it as follows: after the shifting we can write the 3880 | fraction part as f + delta, where 1 < f < 2^25, and 0 <= delta <= 2. 3881 | If delta < 1, do nothing. If delta > 1, add 1 to f. 3882 | If delta == 1, we make sure the rounded number will be even (odd?) 3883 | (after shifting). 3884 btst IMM (0),d0 | is delta < 1? 3885 beq 2f | if so, do not do anything 3886 tstl d1 | is delta == 1? 3887 bne 1f | if so round to even 3888 movel d0,d1 | 3889 andl IMM (2),d1 | bit 1 is the last significant bit 3890 addl d1,d0 | 3891 bra 2f | 3892 1: movel IMM (1),d1 | else add 1 3893 addl d1,d0 | 3894 | Shift right once (because we used bit #FLT_MANT_DIG!). 3895 2: lsrl IMM (1),d0 3896 | Now check again bit #FLT_MANT_DIG (rounding could have produced a 3897 | 'fraction overflow' ...). 3898 btst IMM (FLT_MANT_DIG),d0 3899 beq 1f 3900 lsrl IMM (1),d0 3901 #ifndef __mcoldfire__ 3902 addw IMM (1),d2 3903 #else 3904 addql IMM (1),d2 3905 #endif 3906 1: 3907 | If bit #FLT_MANT_DIG-1 is clear we have a denormalized number, so we 3908 | have to put the exponent to zero and return a denormalized number. 3909 btst IMM (FLT_MANT_DIG-1),d0 3910 beq 1f 3911 jmp a0@ 3912 1: movel IMM (0),d2 3913 jmp a0@ 3914 3915 Lround$to$zero: 3916 Lround$to$plus: 3917 Lround$to$minus: 3918 jmp a0@ 3919 #endif /* L_float */ 3920 3921 | gcc expects the routines __eqdf2, __nedf2, __gtdf2, __gedf2, 3922 | __ledf2, __ltdf2 to all return the same value as a direct call to 3923 | __cmpdf2 would. In this implementation, each of these routines 3924 | simply calls __cmpdf2. It would be more efficient to give the 3925 | __cmpdf2 routine several names, but separating them out will make it 3926 | easier to write efficient versions of these routines someday. 3927 | If the operands recompare unordered unordered __gtdf2 and __gedf2 return -1. 3928 | The other routines return 1. 3929 3930 #ifdef L_eqdf2 3931 .text 3932 FUNC(__eqdf2) 3933 .globl SYM (__eqdf2) 3934 SYM (__eqdf2): 3935 link a6,IMM (0) 3936 pea 1 3937 movl a6@(20),sp@- 3938 movl a6@(16),sp@- 3939 movl a6@(12),sp@- 3940 movl a6@(8),sp@- 3941 PICCALL SYM (__cmpdf2_internal) 3942 unlk a6 3943 rts 3944 #endif /* L_eqdf2 */ 3945 3946 #ifdef L_nedf2 3947 .text 3948 FUNC(__nedf2) 3949 .globl SYM (__nedf2) 3950 SYM (__nedf2): 3951 link a6,IMM (0) 3952 pea 1 3953 movl a6@(20),sp@- 3954 movl a6@(16),sp@- 3955 movl a6@(12),sp@- 3956 movl a6@(8),sp@- 3957 PICCALL SYM (__cmpdf2_internal) 3958 unlk a6 3959 rts 3960 #endif /* L_nedf2 */ 3961 3962 #ifdef L_gtdf2 3963 .text 3964 FUNC(__gtdf2) 3965 .globl SYM (__gtdf2) 3966 SYM (__gtdf2): 3967 link a6,IMM (0) 3968 pea -1 3969 movl a6@(20),sp@- 3970 movl a6@(16),sp@- 3971 movl a6@(12),sp@- 3972 movl a6@(8),sp@- 3973 PICCALL SYM (__cmpdf2_internal) 3974 unlk a6 3975 rts 3976 #endif /* L_gtdf2 */ 3977 3978 #ifdef L_gedf2 3979 .text 3980 FUNC(__gedf2) 3981 .globl SYM (__gedf2) 3982 SYM (__gedf2): 3983 link a6,IMM (0) 3984 pea -1 3985 movl a6@(20),sp@- 3986 movl a6@(16),sp@- 3987 movl a6@(12),sp@- 3988 movl a6@(8),sp@- 3989 PICCALL SYM (__cmpdf2_internal) 3990 unlk a6 3991 rts 3992 #endif /* L_gedf2 */ 3993 3994 #ifdef L_ltdf2 3995 .text 3996 FUNC(__ltdf2) 3997 .globl SYM (__ltdf2) 3998 SYM (__ltdf2): 3999 link a6,IMM (0) 4000 pea 1 4001 movl a6@(20),sp@- 4002 movl a6@(16),sp@- 4003 movl a6@(12),sp@- 4004 movl a6@(8),sp@- 4005 PICCALL SYM (__cmpdf2_internal) 4006 unlk a6 4007 rts 4008 #endif /* L_ltdf2 */ 4009 4010 #ifdef L_ledf2 4011 .text 4012 FUNC(__ledf2) 4013 .globl SYM (__ledf2) 4014 SYM (__ledf2): 4015 link a6,IMM (0) 4016 pea 1 4017 movl a6@(20),sp@- 4018 movl a6@(16),sp@- 4019 movl a6@(12),sp@- 4020 movl a6@(8),sp@- 4021 PICCALL SYM (__cmpdf2_internal) 4022 unlk a6 4023 rts 4024 #endif /* L_ledf2 */ 4025 4026 | The comments above about __eqdf2, et. al., also apply to __eqsf2, 4027 | et. al., except that the latter call __cmpsf2 rather than __cmpdf2. 4028 4029 #ifdef L_eqsf2 4030 .text 4031 FUNC(__eqsf2) 4032 .globl SYM (__eqsf2) 4033 SYM (__eqsf2): 4034 link a6,IMM (0) 4035 pea 1 4036 movl a6@(12),sp@- 4037 movl a6@(8),sp@- 4038 PICCALL SYM (__cmpsf2_internal) 4039 unlk a6 4040 rts 4041 #endif /* L_eqsf2 */ 4042 4043 #ifdef L_nesf2 4044 .text 4045 FUNC(__nesf2) 4046 .globl SYM (__nesf2) 4047 SYM (__nesf2): 4048 link a6,IMM (0) 4049 pea 1 4050 movl a6@(12),sp@- 4051 movl a6@(8),sp@- 4052 PICCALL SYM (__cmpsf2_internal) 4053 unlk a6 4054 rts 4055 #endif /* L_nesf2 */ 4056 4057 #ifdef L_gtsf2 4058 .text 4059 FUNC(__gtsf2) 4060 .globl SYM (__gtsf2) 4061 SYM (__gtsf2): 4062 link a6,IMM (0) 4063 pea -1 4064 movl a6@(12),sp@- 4065 movl a6@(8),sp@- 4066 PICCALL SYM (__cmpsf2_internal) 4067 unlk a6 4068 rts 4069 #endif /* L_gtsf2 */ 4070 4071 #ifdef L_gesf2 4072 .text 4073 FUNC(__gesf2) 4074 .globl SYM (__gesf2) 4075 SYM (__gesf2): 4076 link a6,IMM (0) 4077 pea -1 4078 movl a6@(12),sp@- 4079 movl a6@(8),sp@- 4080 PICCALL SYM (__cmpsf2_internal) 4081 unlk a6 4082 rts 4083 #endif /* L_gesf2 */ 4084 4085 #ifdef L_ltsf2 4086 .text 4087 FUNC(__ltsf2) 4088 .globl SYM (__ltsf2) 4089 SYM (__ltsf2): 4090 link a6,IMM (0) 4091 pea 1 4092 movl a6@(12),sp@- 4093 movl a6@(8),sp@- 4094 PICCALL SYM (__cmpsf2_internal) 4095 unlk a6 4096 rts 4097 #endif /* L_ltsf2 */ 4098 4099 #ifdef L_lesf2 4100 .text 4101 FUNC(__lesf2) 4102 .globl SYM (__lesf2) 4103 SYM (__lesf2): 4104 link a6,IMM (0) 4105 pea 1 4106 movl a6@(12),sp@- 4107 movl a6@(8),sp@- 4108 PICCALL SYM (__cmpsf2_internal) 4109 unlk a6 4110 rts 4111 #endif /* L_lesf2 */ 4112 4113 #if defined (__ELF__) && defined (__linux__) 4114 /* Make stack non-executable for ELF linux targets. */ 4115 .section .note.GNU-stack,"",@progbits 4116 #endif 4117
Indexes created Sat Feb 28 05:31:39 UTC 2026