Home | History | Annotate | Line # | Download | only in softfloat
softfloat-specialize revision 1.10
      1 /*	$NetBSD: softfloat-specialize,v 1.10 2025/04/27 15:23:27 riastradh Exp $	*/
      2 
      3 /* This is a derivative work. */
      4 
      5 /*
      6 ===============================================================================
      7 
      8 This C source fragment is part of the SoftFloat IEC/IEEE Floating-point
      9 Arithmetic Package, Release 2a.
     10 
     11 Written by John R. Hauser.  This work was made possible in part by the
     12 International Computer Science Institute, located at Suite 600, 1947 Center
     13 Street, Berkeley, California 94704.  Funding was partially provided by the
     14 National Science Foundation under grant MIP-9311980.  The original version
     15 of this code was written as part of a project to build a fixed-point vector
     16 processor in collaboration with the University of California at Berkeley,
     17 overseen by Profs. Nelson Morgan and John Wawrzynek.  More information
     18 is available through the Web page `http://HTTP.CS.Berkeley.EDU/~jhauser/
     19 arithmetic/SoftFloat.html'.
     20 
     21 THIS SOFTWARE IS DISTRIBUTED AS IS, FOR FREE.  Although reasonable effort
     22 has been made to avoid it, THIS SOFTWARE MAY CONTAIN FAULTS THAT WILL AT
     23 TIMES RESULT IN INCORRECT BEHAVIOR.  USE OF THIS SOFTWARE IS RESTRICTED TO
     24 PERSONS AND ORGANIZATIONS WHO CAN AND WILL TAKE FULL RESPONSIBILITY FOR ANY
     25 AND ALL LOSSES, COSTS, OR OTHER PROBLEMS ARISING FROM ITS USE.
     26 
     27 Derivative works are acceptable, even for commercial purposes, so long as
     28 (1) they include prominent notice that the work is derivative, and (2) they
     29 include prominent notice akin to these four paragraphs for those parts of
     30 this code that are retained.
     31 
     32 ===============================================================================
     33 */
     34 
     35 #include <signal.h>
     36 #include <string.h>
     37 #include <unistd.h>
     38 
     39 #include "reentrant.h"
     40 
     41 /*
     42 -------------------------------------------------------------------------------
     43 Underflow tininess-detection mode, statically initialized to default value.
     44 (The declaration in `softfloat.h' must match the `int8' type here.)
     45 -------------------------------------------------------------------------------
     46 */
     47 #ifdef SOFTFLOAT_FOR_GCC
     48 static
     49 #endif
     50 int8 float_detect_tininess = float_tininess_after_rounding;
     51 
     52 /*
     53 -------------------------------------------------------------------------------
     54 Raises the exceptions specified by `flags'.  Floating-point traps can be
     55 defined here if desired.  It is currently not possible for such a trap to
     56 substitute a result value.  If traps are not implemented, this routine
     57 should be simply `float_exception_flags |= flags;'.
     58 -------------------------------------------------------------------------------
     59 */
     60 #ifdef SOFTFLOAT_FOR_GCC
     61 #ifndef set_float_exception_mask
     62 #define float_exception_mask	_softfloat_float_exception_mask
     63 #endif
     64 #endif
     65 #ifndef set_float_exception_mask
     66 fp_except float_exception_mask = 0;
     67 #endif
     68 void
     69 float_raise( fp_except newflags )
     70 {
     71     siginfo_t info;
     72     struct sigaction sa;
     73     sigset_t sigmask, osigmask;
     74     fp_except flags;
     75 
     76     for (;;) {
     77 #ifdef set_float_exception_mask
     78 	flags = newflags | set_float_exception_flags(newflags, 0);
     79 #else
     80 	float_exception_flags |= newflags;
     81 	flags = float_exception_flags;
     82 #endif
     83 
     84 	/*
     85 	 * If none of the sticky flags are trapped (i.e., enabled in
     86 	 * float_exception_mask), we're done.  Trapping is unusual and
     87 	 * costly anyway, so take the non-trapping path as the fast
     88 	 * path.
     89 	 */
     90 	flags &= float_exception_mask;
     91 	if (__predict_true(flags == 0))
     92 	    break;
     93 
     94 	/*
     95 	 * Block all signals while we figure out how to deliver a
     96 	 * non-maskable (as a signal), non-ignorable SIGFPE, and obtain
     97 	 * the current signal mask.
     98 	 */
     99 	sigfillset(&sigmask);
    100 	thr_sigsetmask(SIG_BLOCK, &sigmask, &osigmask);
    101 
    102 	/*
    103 	 * Find the current signal disposition of SIGFPE.
    104 	 */
    105 	sigaction(SIGFPE, NULL, &sa);
    106 
    107 	/*
    108 	 * If SIGFPE is masked or ignored, unmask it and reset it to
    109 	 * the default disposition to deliver the signal.
    110 	 */
    111 	if (sigismember(&osigmask, SIGFPE) ||
    112 	    ((sa.sa_flags & SA_SIGINFO) == 0 &&
    113 		sa.sa_handler == SIG_IGN)) {
    114 		/*
    115 		 * Prepare to unmask SIGFPE.  This will take effect
    116 		 * when we use thr_sigsetmask(SIG_SETMASK, ...) below,
    117 		 * once the signal has been queued, so that it happens
    118 		 * atomically with respect to other signal delivery.
    119 		 */
    120 		sigdelset(&osigmask, SIGFPE);
    121 
    122 		/*
    123 		 * Reset SIGFPE to the default disposition, which is to
    124 		 * terminate the process.
    125 		 */
    126 		memset(&sa, 0, sizeof(sa));
    127 		sa.sa_handler = SIG_DFL;
    128 		sigemptyset(&sa.sa_mask);
    129 		sa.sa_flags = 0;
    130 		sigaction(SIGFPE, &sa, NULL);
    131 	}
    132 
    133 	/*
    134 	 * Queue the signal for delivery.  It won't trigger the signal
    135 	 * handler yet, because it's still masked, but as soon as we
    136 	 * unmask it either the process will terminate or the signal
    137 	 * handler will be called.
    138 	 */
    139 	memset(&info, 0, sizeof info);
    140 	info.si_signo = SIGFPE;
    141 	info.si_pid = getpid();
    142 	info.si_uid = geteuid();
    143 	if (flags & float_flag_underflow)
    144 	    info.si_code = FPE_FLTUND;
    145 	else if (flags & float_flag_overflow)
    146 	    info.si_code = FPE_FLTOVF;
    147 	else if (flags & float_flag_divbyzero)
    148 	    info.si_code = FPE_FLTDIV;
    149 	else if (flags & float_flag_invalid)
    150 	    info.si_code = FPE_FLTINV;
    151 	else if (flags & float_flag_inexact)
    152 	    info.si_code = FPE_FLTRES;
    153 	sigqueueinfo(getpid(), &info);
    154 
    155 	/*
    156 	 * Restore the old signal mask, except with SIGFPE unmasked
    157 	 * even if it was masked before.
    158 	 *
    159 	 * At this point, either the process will terminate (if SIGFPE
    160 	 * had or now has the default disposition) or the signal
    161 	 * handler will be called (if SIGFPE had a non-default,
    162 	 * non-ignored disposition).
    163 	 *
    164 	 * If the signal handler returns, it can't change the set of
    165 	 * exceptions raised by this floating-point operation -- but it
    166 	 * can change the sticky set from previous operations, and it
    167 	 * can change the set of exceptions that are trapped, so loop
    168 	 * around; next time we might make progress instead of calling
    169 	 * the signal handler again.
    170 	 */
    171 	thr_sigsetmask(SIG_SETMASK, &osigmask, NULL);
    172     }
    173 }
    174 #undef float_exception_mask
    175 
    176 /*
    177 -------------------------------------------------------------------------------
    178 Internal canonical NaN format.
    179 -------------------------------------------------------------------------------
    180 */
    181 typedef struct {
    182     flag sign;
    183     bits64 high, low;
    184 } commonNaNT;
    185 
    186 /*
    187 -------------------------------------------------------------------------------
    188 The pattern for a default generated single-precision NaN.
    189 -------------------------------------------------------------------------------
    190 */
    191 #define float32_default_nan 0xFFFFFFFF
    192 
    193 /*
    194 -------------------------------------------------------------------------------
    195 Returns 1 if the single-precision floating-point value `a' is a NaN;
    196 otherwise returns 0.
    197 -------------------------------------------------------------------------------
    198 */
    199 #ifdef SOFTFLOAT_FOR_GCC
    200 static
    201 #endif
    202 flag float32_is_nan( float32 a )
    203 {
    204 
    205     return ( (bits32)0xFF000000 < (bits32) ( a<<1 ) );
    206 
    207 }
    208 
    209 /*
    210 -------------------------------------------------------------------------------
    211 Returns 1 if the single-precision floating-point value `a' is a signaling
    212 NaN; otherwise returns 0.
    213 -------------------------------------------------------------------------------
    214 */
    215 #if defined(SOFTFLOAT_FOR_GCC) \
    216     && !defined(SOFTFLOATAARCH64_FOR_GCC) \
    217     && !defined(SOFTFLOATSPARC64_FOR_GCC) \
    218     && !defined(SOFTFLOATM68K_FOR_GCC)
    219 static
    220 #endif
    221 flag float32_is_signaling_nan( float32 a )
    222 {
    223 
    224     return ( ( ( a>>22 ) & 0x1FF ) == 0x1FE ) && ( a & 0x003FFFFF );
    225 
    226 }
    227 
    228 /*
    229 -------------------------------------------------------------------------------
    230 Returns the result of converting the single-precision floating-point NaN
    231 `a' to the canonical NaN format.  If `a' is a signaling NaN, the invalid
    232 exception is raised.
    233 -------------------------------------------------------------------------------
    234 */
    235 static commonNaNT float32ToCommonNaN( float32 a )
    236 {
    237     commonNaNT z;
    238 
    239     if ( float32_is_signaling_nan( a ) ) float_raise( float_flag_invalid );
    240     z.sign = a>>31;
    241     z.low = 0;
    242     z.high = ( (bits64) a )<<41;
    243     return z;
    244 
    245 }
    246 
    247 /*
    248 -------------------------------------------------------------------------------
    249 Returns the result of converting the canonical NaN `a' to the single-
    250 precision floating-point format.
    251 -------------------------------------------------------------------------------
    252 */
    253 static float32 commonNaNToFloat32( commonNaNT a )
    254 {
    255 
    256     return ( ( (bits32) a.sign )<<31 ) | 0x7FC00000 | (bits32)( a.high>>41 );
    257 
    258 }
    259 
    260 /*
    261 -------------------------------------------------------------------------------
    262 Takes two single-precision floating-point values `a' and `b', one of which
    263 is a NaN, and returns the appropriate NaN result.  If either `a' or `b' is a
    264 signaling NaN, the invalid exception is raised.
    265 -------------------------------------------------------------------------------
    266 */
    267 static float32 propagateFloat32NaN( float32 a, float32 b )
    268 {
    269     flag aIsNaN, aIsSignalingNaN, bIsNaN, bIsSignalingNaN;
    270 
    271     aIsNaN = float32_is_nan( a );
    272     aIsSignalingNaN = float32_is_signaling_nan( a );
    273     bIsNaN = float32_is_nan( b );
    274     bIsSignalingNaN = float32_is_signaling_nan( b );
    275     a |= 0x00400000;
    276     b |= 0x00400000;
    277     if ( aIsSignalingNaN | bIsSignalingNaN ) float_raise( float_flag_invalid );
    278     if ( aIsNaN ) {
    279         return ( aIsSignalingNaN & bIsNaN ) ? b : a;
    280     }
    281     else {
    282         return b;
    283     }
    284 
    285 }
    286 
    287 /*
    288 -------------------------------------------------------------------------------
    289 The pattern for a default generated double-precision NaN.
    290 -------------------------------------------------------------------------------
    291 */
    292 #define float64_default_nan LIT64( 0xFFFFFFFFFFFFFFFF )
    293 
    294 /*
    295 -------------------------------------------------------------------------------
    296 Returns 1 if the double-precision floating-point value `a' is a NaN;
    297 otherwise returns 0.
    298 -------------------------------------------------------------------------------
    299 */
    300 #ifdef SOFTFLOAT_FOR_GCC
    301 static
    302 #endif
    303 flag float64_is_nan( float64 a )
    304 {
    305 
    306     return ( (bits64)LIT64( 0xFFE0000000000000 ) <
    307 	     (bits64) ( FLOAT64_DEMANGLE(a)<<1 ) );
    308 
    309 }
    310 
    311 /*
    312 -------------------------------------------------------------------------------
    313 Returns 1 if the double-precision floating-point value `a' is a signaling
    314 NaN; otherwise returns 0.
    315 -------------------------------------------------------------------------------
    316 */
    317 #if defined(SOFTFLOAT_FOR_GCC) \
    318     && !defined(SOFTFLOATAARCH64_FOR_GCC) \
    319     && !defined(SOFTFLOATSPARC64_FOR_GCC) \
    320     && !defined(SOFTFLOATM68K_FOR_GCC)
    321 static
    322 #endif
    323 flag float64_is_signaling_nan( float64 a )
    324 {
    325 
    326     return
    327            ( ( ( FLOAT64_DEMANGLE(a)>>51 ) & 0xFFF ) == 0xFFE )
    328         && ( FLOAT64_DEMANGLE(a) & LIT64( 0x0007FFFFFFFFFFFF ) );
    329 
    330 }
    331 
    332 /*
    333 -------------------------------------------------------------------------------
    334 Returns the result of converting the double-precision floating-point NaN
    335 `a' to the canonical NaN format.  If `a' is a signaling NaN, the invalid
    336 exception is raised.
    337 -------------------------------------------------------------------------------
    338 */
    339 static commonNaNT float64ToCommonNaN( float64 a )
    340 {
    341     commonNaNT z;
    342 
    343     if ( float64_is_signaling_nan( a ) ) float_raise( float_flag_invalid );
    344     z.sign = (flag)(FLOAT64_DEMANGLE(a)>>63);
    345     z.low = 0;
    346     z.high = FLOAT64_DEMANGLE(a)<<12;
    347     return z;
    348 
    349 }
    350 
    351 /*
    352 -------------------------------------------------------------------------------
    353 Returns the result of converting the canonical NaN `a' to the double-
    354 precision floating-point format.
    355 -------------------------------------------------------------------------------
    356 */
    357 static float64 commonNaNToFloat64( commonNaNT a )
    358 {
    359 
    360     return FLOAT64_MANGLE(
    361 	( ( (bits64) a.sign )<<63 )
    362         | LIT64( 0x7FF8000000000000 )
    363         | ( a.high>>12 ) );
    364 
    365 }
    366 
    367 /*
    368 -------------------------------------------------------------------------------
    369 Takes two double-precision floating-point values `a' and `b', one of which
    370 is a NaN, and returns the appropriate NaN result.  If either `a' or `b' is a
    371 signaling NaN, the invalid exception is raised.
    372 -------------------------------------------------------------------------------
    373 */
    374 static float64 propagateFloat64NaN( float64 a, float64 b )
    375 {
    376     flag aIsNaN, aIsSignalingNaN, bIsNaN, bIsSignalingNaN;
    377 
    378     aIsNaN = float64_is_nan( a );
    379     aIsSignalingNaN = float64_is_signaling_nan( a );
    380     bIsNaN = float64_is_nan( b );
    381     bIsSignalingNaN = float64_is_signaling_nan( b );
    382     a |= FLOAT64_MANGLE(LIT64( 0x0008000000000000 ));
    383     b |= FLOAT64_MANGLE(LIT64( 0x0008000000000000 ));
    384     if ( aIsSignalingNaN | bIsSignalingNaN ) float_raise( float_flag_invalid );
    385     if ( aIsNaN ) {
    386         return ( aIsSignalingNaN & bIsNaN ) ? b : a;
    387     }
    388     else {
    389         return b;
    390     }
    391 
    392 }
    393 
    394 #ifdef FLOATX80
    395 
    396 /*
    397 -------------------------------------------------------------------------------
    398 The pattern for a default generated extended double-precision NaN.  The
    399 `high' and `low' values hold the most- and least-significant bits,
    400 respectively.
    401 -------------------------------------------------------------------------------
    402 */
    403 #define floatx80_default_nan_high 0xFFFF
    404 #define floatx80_default_nan_low  LIT64( 0xFFFFFFFFFFFFFFFF )
    405 
    406 /*
    407 -------------------------------------------------------------------------------
    408 Returns 1 if the extended double-precision floating-point value `a' is a
    409 NaN; otherwise returns 0.
    410 -------------------------------------------------------------------------------
    411 */
    412 flag floatx80_is_nan( floatx80 a )
    413 {
    414 
    415     return ( ( a.high & 0x7FFF ) == 0x7FFF ) && (bits64) ( a.low<<1 );
    416 
    417 }
    418 
    419 /*
    420 -------------------------------------------------------------------------------
    421 Returns 1 if the extended double-precision floating-point value `a' is a
    422 signaling NaN; otherwise returns 0.
    423 -------------------------------------------------------------------------------
    424 */
    425 flag floatx80_is_signaling_nan( floatx80 a )
    426 {
    427     bits64 aLow;
    428 
    429     aLow = a.low & ~ LIT64( 0x4000000000000000 );
    430     return
    431            ( ( a.high & 0x7FFF ) == 0x7FFF )
    432         && (bits64) ( aLow<<1 )
    433         && ( a.low == aLow );
    434 
    435 }
    436 
    437 /*
    438 -------------------------------------------------------------------------------
    439 Returns the result of converting the extended double-precision floating-
    440 point NaN `a' to the canonical NaN format.  If `a' is a signaling NaN, the
    441 invalid exception is raised.
    442 -------------------------------------------------------------------------------
    443 */
    444 static commonNaNT floatx80ToCommonNaN( floatx80 a )
    445 {
    446     commonNaNT z;
    447 
    448     if ( floatx80_is_signaling_nan( a ) ) float_raise( float_flag_invalid );
    449     z.sign = a.high>>15;
    450     z.low = 0;
    451     z.high = a.low<<1;
    452     return z;
    453 
    454 }
    455 
    456 /*
    457 -------------------------------------------------------------------------------
    458 Returns the result of converting the canonical NaN `a' to the extended
    459 double-precision floating-point format.
    460 -------------------------------------------------------------------------------
    461 */
    462 static floatx80 commonNaNToFloatx80( commonNaNT a )
    463 {
    464     floatx80 z;
    465 
    466     z.low = LIT64( 0xC000000000000000 ) | ( a.high>>1 );
    467     z.high = ( ( (bits16) a.sign )<<15 ) | 0x7FFF;
    468     return z;
    469 
    470 }
    471 
    472 /*
    473 -------------------------------------------------------------------------------
    474 Takes two extended double-precision floating-point values `a' and `b', one
    475 of which is a NaN, and returns the appropriate NaN result.  If either `a' or
    476 `b' is a signaling NaN, the invalid exception is raised.
    477 -------------------------------------------------------------------------------
    478 */
    479 static floatx80 propagateFloatx80NaN( floatx80 a, floatx80 b )
    480 {
    481     flag aIsNaN, aIsSignalingNaN, bIsNaN, bIsSignalingNaN;
    482 
    483     aIsNaN = floatx80_is_nan( a );
    484     aIsSignalingNaN = floatx80_is_signaling_nan( a );
    485     bIsNaN = floatx80_is_nan( b );
    486     bIsSignalingNaN = floatx80_is_signaling_nan( b );
    487     a.low |= LIT64( 0xC000000000000000 );
    488     b.low |= LIT64( 0xC000000000000000 );
    489     if ( aIsSignalingNaN | bIsSignalingNaN ) float_raise( float_flag_invalid );
    490     if ( aIsNaN ) {
    491         return ( aIsSignalingNaN & bIsNaN ) ? b : a;
    492     }
    493     else {
    494         return b;
    495     }
    496 
    497 }
    498 
    499 #endif
    500 
    501 #ifdef FLOAT128
    502 
    503 /*
    504 -------------------------------------------------------------------------------
    505 The pattern for a default generated quadruple-precision NaN.  The `high' and
    506 `low' values hold the most- and least-significant bits, respectively.
    507 -------------------------------------------------------------------------------
    508 */
    509 #define float128_default_nan_high LIT64( 0xFFFFFFFFFFFFFFFF )
    510 #define float128_default_nan_low  LIT64( 0xFFFFFFFFFFFFFFFF )
    511 
    512 /*
    513 -------------------------------------------------------------------------------
    514 Returns 1 if the quadruple-precision floating-point value `a' is a NaN;
    515 otherwise returns 0.
    516 -------------------------------------------------------------------------------
    517 */
    518 flag float128_is_nan( float128 a )
    519 {
    520 
    521     return
    522            ( (bits64)LIT64( 0xFFFE000000000000 ) <= (bits64) ( a.high<<1 ) )
    523         && ( a.low || ( a.high & LIT64( 0x0000FFFFFFFFFFFF ) ) );
    524 
    525 }
    526 
    527 /*
    528 -------------------------------------------------------------------------------
    529 Returns 1 if the quadruple-precision floating-point value `a' is a
    530 signaling NaN; otherwise returns 0.
    531 -------------------------------------------------------------------------------
    532 */
    533 flag float128_is_signaling_nan( float128 a )
    534 {
    535 
    536     return
    537            ( ( ( a.high>>47 ) & 0xFFFF ) == 0xFFFE )
    538         && ( a.low || ( a.high & LIT64( 0x00007FFFFFFFFFFF ) ) );
    539 
    540 }
    541 
    542 /*
    543 -------------------------------------------------------------------------------
    544 Returns the result of converting the quadruple-precision floating-point NaN
    545 `a' to the canonical NaN format.  If `a' is a signaling NaN, the invalid
    546 exception is raised.
    547 -------------------------------------------------------------------------------
    548 */
    549 static commonNaNT float128ToCommonNaN( float128 a )
    550 {
    551     commonNaNT z;
    552 
    553     if ( float128_is_signaling_nan( a ) ) float_raise( float_flag_invalid );
    554     z.sign = (flag)(a.high>>63);
    555     shortShift128Left( a.high, a.low, 16, &z.high, &z.low );
    556     return z;
    557 
    558 }
    559 
    560 /*
    561 -------------------------------------------------------------------------------
    562 Returns the result of converting the canonical NaN `a' to the quadruple-
    563 precision floating-point format.
    564 -------------------------------------------------------------------------------
    565 */
    566 static float128 commonNaNToFloat128( commonNaNT a )
    567 {
    568     float128 z;
    569 
    570     shift128Right( a.high, a.low, 16, &z.high, &z.low );
    571     z.high |= ( ( (bits64) a.sign )<<63 ) | LIT64( 0x7FFF800000000000 );
    572     return z;
    573 
    574 }
    575 
    576 /*
    577 -------------------------------------------------------------------------------
    578 Takes two quadruple-precision floating-point values `a' and `b', one of
    579 which is a NaN, and returns the appropriate NaN result.  If either `a' or
    580 `b' is a signaling NaN, the invalid exception is raised.
    581 -------------------------------------------------------------------------------
    582 */
    583 static float128 propagateFloat128NaN( float128 a, float128 b )
    584 {
    585     flag aIsNaN, aIsSignalingNaN, bIsNaN, bIsSignalingNaN;
    586 
    587     aIsNaN = float128_is_nan( a );
    588     aIsSignalingNaN = float128_is_signaling_nan( a );
    589     bIsNaN = float128_is_nan( b );
    590     bIsSignalingNaN = float128_is_signaling_nan( b );
    591     a.high |= LIT64( 0x0000800000000000 );
    592     b.high |= LIT64( 0x0000800000000000 );
    593     if ( aIsSignalingNaN | bIsSignalingNaN ) float_raise( float_flag_invalid );
    594     if ( aIsNaN ) {
    595         return ( aIsSignalingNaN & bIsNaN ) ? b : a;
    596     }
    597     else {
    598         return b;
    599     }
    600 
    601 }
    602 
    603 #endif
    604 
    605