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