trionan.c 23 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914
  1. /*************************************************************************
  2. *
  3. * $Id$
  4. *
  5. * Copyright (C) 2001 Bjorn Reese <breese@users.sourceforge.net>
  6. *
  7. * Permission to use, copy, modify, and distribute this software for any
  8. * purpose with or without fee is hereby granted, provided that the above
  9. * copyright notice and this permission notice appear in all copies.
  10. *
  11. * THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR IMPLIED
  12. * WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED WARRANTIES OF
  13. * MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE AUTHORS AND
  14. * CONTRIBUTORS ACCEPT NO RESPONSIBILITY IN ANY CONCEIVABLE MANNER.
  15. *
  16. ************************************************************************
  17. *
  18. * Functions to handle special quantities in floating-point numbers
  19. * (that is, NaNs and infinity). They provide the capability to detect
  20. * and fabricate special quantities.
  21. *
  22. * Although written to be as portable as possible, it can never be
  23. * guaranteed to work on all platforms, as not all hardware supports
  24. * special quantities.
  25. *
  26. * The approach used here (approximately) is to:
  27. *
  28. * 1. Use C99 functionality when available.
  29. * 2. Use IEEE 754 bit-patterns if possible.
  30. * 3. Use platform-specific techniques.
  31. *
  32. ************************************************************************/
  33. /*
  34. * TODO:
  35. * o Put all the magic into trio_fpclassify_and_signbit(), and use this from
  36. * trio_isnan() etc.
  37. */
  38. /*************************************************************************
  39. * Include files
  40. */
  41. #include "triodef.h"
  42. #include "trionan.h"
  43. #include <math.h>
  44. #include <string.h>
  45. #include <limits.h>
  46. #include <float.h>
  47. #if defined(TRIO_PLATFORM_UNIX)
  48. # include <signal.h>
  49. #endif
  50. #if defined(TRIO_COMPILER_DECC)
  51. # if defined(__linux__)
  52. # include <cpml.h>
  53. # else
  54. # include <fp_class.h>
  55. # endif
  56. #endif
  57. #include <assert.h>
  58. #if defined(TRIO_DOCUMENTATION)
  59. # include "doc/doc_nan.h"
  60. #endif
  61. /** @addtogroup SpecialQuantities
  62. @{
  63. */
  64. /*************************************************************************
  65. * Definitions
  66. */
  67. #define TRIO_TRUE (1 == 1)
  68. #define TRIO_FALSE (0 == 1)
  69. /*
  70. * We must enable IEEE floating-point on Alpha
  71. */
  72. #if defined(__alpha) && !defined(_IEEE_FP)
  73. # if defined(TRIO_COMPILER_DECC)
  74. # if defined(TRIO_PLATFORM_VMS)
  75. # error "Must be compiled with option /IEEE_MODE=UNDERFLOW_TO_ZERO/FLOAT=IEEE"
  76. # else
  77. # if !defined(_CFE)
  78. # error "Must be compiled with option -ieee"
  79. # endif
  80. # endif
  81. # elif defined(TRIO_COMPILER_GCC) && (defined(__osf__) || defined(__linux__))
  82. # error "Must be compiled with option -mieee"
  83. # endif
  84. #endif /* __alpha && ! _IEEE_FP */
  85. /*
  86. * In ANSI/IEEE 754-1985 64-bits double format numbers have the
  87. * following properties (amongst others)
  88. *
  89. * o FLT_RADIX == 2: binary encoding
  90. * o DBL_MAX_EXP == 1024: 11 bits exponent, where one bit is used
  91. * to indicate special numbers (e.g. NaN and Infinity), so the
  92. * maximum exponent is 10 bits wide (2^10 == 1024).
  93. * o DBL_MANT_DIG == 53: The mantissa is 52 bits wide, but because
  94. * numbers are normalized the initial binary 1 is represented
  95. * implicitly (the so-called "hidden bit"), which leaves us with
  96. * the ability to represent 53 bits wide mantissa.
  97. */
  98. #if (FLT_RADIX == 2) && (DBL_MAX_EXP == 1024) && (DBL_MANT_DIG == 53)
  99. # define USE_IEEE_754
  100. #endif
  101. /*************************************************************************
  102. * Constants
  103. */
  104. static TRIO_CONST char rcsid[] = "@(#)$Id$";
  105. #if defined(USE_IEEE_754)
  106. /*
  107. * Endian-agnostic indexing macro.
  108. *
  109. * The value of internalEndianMagic, when converted into a 64-bit
  110. * integer, becomes 0x0706050403020100 (we could have used a 64-bit
  111. * integer value instead of a double, but not all platforms supports
  112. * that type). The value is automatically encoded with the correct
  113. * endianness by the compiler, which means that we can support any
  114. * kind of endianness. The individual bytes are then used as an index
  115. * for the IEEE 754 bit-patterns and masks.
  116. */
  117. #define TRIO_DOUBLE_INDEX(x) (((unsigned char *)&internalEndianMagic)[7-(x)])
  118. #if (defined(__BORLANDC__) && __BORLANDC__ >= 0x0590)
  119. static TRIO_CONST double internalEndianMagic = 7.949928895127362e-275;
  120. #else
  121. static TRIO_CONST double internalEndianMagic = 7.949928895127363e-275;
  122. #endif
  123. /* Mask for the exponent */
  124. static TRIO_CONST unsigned char ieee_754_exponent_mask[] = {
  125. 0x7F, 0xF0, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00
  126. };
  127. /* Mask for the mantissa */
  128. static TRIO_CONST unsigned char ieee_754_mantissa_mask[] = {
  129. 0x00, 0x0F, 0xFF, 0xFF, 0xFF, 0xFF, 0xFF, 0xFF
  130. };
  131. /* Mask for the sign bit */
  132. static TRIO_CONST unsigned char ieee_754_sign_mask[] = {
  133. 0x80, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00
  134. };
  135. /* Bit-pattern for negative zero */
  136. static TRIO_CONST unsigned char ieee_754_negzero_array[] = {
  137. 0x80, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00
  138. };
  139. /* Bit-pattern for infinity */
  140. static TRIO_CONST unsigned char ieee_754_infinity_array[] = {
  141. 0x7F, 0xF0, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00
  142. };
  143. /* Bit-pattern for quiet NaN */
  144. static TRIO_CONST unsigned char ieee_754_qnan_array[] = {
  145. 0x7F, 0xF8, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00
  146. };
  147. /*************************************************************************
  148. * Functions
  149. */
  150. /*
  151. * trio_make_double
  152. */
  153. TRIO_PRIVATE double
  154. trio_make_double
  155. TRIO_ARGS1((values),
  156. TRIO_CONST unsigned char *values)
  157. {
  158. TRIO_VOLATILE double result;
  159. int i;
  160. for (i = 0; i < (int)sizeof(double); i++) {
  161. ((TRIO_VOLATILE unsigned char *)&result)[TRIO_DOUBLE_INDEX(i)] = values[i];
  162. }
  163. return result;
  164. }
  165. /*
  166. * trio_is_special_quantity
  167. */
  168. TRIO_PRIVATE int
  169. trio_is_special_quantity
  170. TRIO_ARGS2((number, has_mantissa),
  171. double number,
  172. int *has_mantissa)
  173. {
  174. unsigned int i;
  175. unsigned char current;
  176. int is_special_quantity = TRIO_TRUE;
  177. *has_mantissa = 0;
  178. for (i = 0; i < (unsigned int)sizeof(double); i++) {
  179. current = ((unsigned char *)&number)[TRIO_DOUBLE_INDEX(i)];
  180. is_special_quantity
  181. &= ((current & ieee_754_exponent_mask[i]) == ieee_754_exponent_mask[i]);
  182. *has_mantissa |= (current & ieee_754_mantissa_mask[i]);
  183. }
  184. return is_special_quantity;
  185. }
  186. /*
  187. * trio_is_negative
  188. */
  189. TRIO_PRIVATE int
  190. trio_is_negative
  191. TRIO_ARGS1((number),
  192. double number)
  193. {
  194. unsigned int i;
  195. int is_negative = TRIO_FALSE;
  196. for (i = 0; i < (unsigned int)sizeof(double); i++) {
  197. is_negative |= (((unsigned char *)&number)[TRIO_DOUBLE_INDEX(i)]
  198. & ieee_754_sign_mask[i]);
  199. }
  200. return is_negative;
  201. }
  202. #endif /* USE_IEEE_754 */
  203. /**
  204. Generate negative zero.
  205. @return Floating-point representation of negative zero.
  206. */
  207. TRIO_PUBLIC double
  208. trio_nzero(TRIO_NOARGS)
  209. {
  210. #if defined(USE_IEEE_754)
  211. return trio_make_double(ieee_754_negzero_array);
  212. #else
  213. TRIO_VOLATILE double zero = 0.0;
  214. return -zero;
  215. #endif
  216. }
  217. /**
  218. Generate positive infinity.
  219. @return Floating-point representation of positive infinity.
  220. */
  221. TRIO_PUBLIC double
  222. trio_pinf(TRIO_NOARGS)
  223. {
  224. /* Cache the result */
  225. static double result = 0.0;
  226. if (result == 0.0) {
  227. #if defined(INFINITY) && defined(__STDC_IEC_559__)
  228. result = (double)INFINITY;
  229. #elif defined(USE_IEEE_754)
  230. result = trio_make_double(ieee_754_infinity_array);
  231. #else
  232. /*
  233. * If HUGE_VAL is different from DBL_MAX, then HUGE_VAL is used
  234. * as infinity. Otherwise we have to resort to an overflow
  235. * operation to generate infinity.
  236. */
  237. # if defined(TRIO_PLATFORM_UNIX)
  238. void (*signal_handler)(int) = signal(SIGFPE, SIG_IGN);
  239. # endif
  240. result = HUGE_VAL;
  241. if (HUGE_VAL == DBL_MAX) {
  242. /* Force overflow */
  243. result += HUGE_VAL;
  244. }
  245. # if defined(TRIO_PLATFORM_UNIX)
  246. signal(SIGFPE, signal_handler);
  247. # endif
  248. #endif
  249. }
  250. return result;
  251. }
  252. /**
  253. Generate negative infinity.
  254. @return Floating-point value of negative infinity.
  255. */
  256. TRIO_PUBLIC double
  257. trio_ninf(TRIO_NOARGS)
  258. {
  259. static double result = 0.0;
  260. if (result == 0.0) {
  261. /*
  262. * Negative infinity is calculated by negating positive infinity,
  263. * which can be done because it is legal to do calculations on
  264. * infinity (for example, 1 / infinity == 0).
  265. */
  266. result = -trio_pinf();
  267. }
  268. return result;
  269. }
  270. /**
  271. Generate NaN.
  272. @return Floating-point representation of NaN.
  273. */
  274. TRIO_PUBLIC double
  275. trio_nan(TRIO_NOARGS)
  276. {
  277. /* Cache the result */
  278. static double result = 0.0;
  279. if (result == 0.0) {
  280. #if defined(TRIO_COMPILER_SUPPORTS_C99)
  281. result = nan("");
  282. #elif defined(NAN) && defined(__STDC_IEC_559__)
  283. result = (double)NAN;
  284. #elif defined(USE_IEEE_754)
  285. result = trio_make_double(ieee_754_qnan_array);
  286. #else
  287. /*
  288. * There are several ways to generate NaN. The one used here is
  289. * to divide infinity by infinity. I would have preferred to add
  290. * negative infinity to positive infinity, but that yields wrong
  291. * result (infinity) on FreeBSD.
  292. *
  293. * This may fail if the hardware does not support NaN, or if
  294. * the Invalid Operation floating-point exception is unmasked.
  295. */
  296. # if defined(TRIO_PLATFORM_UNIX)
  297. void (*signal_handler)(int) = signal(SIGFPE, SIG_IGN);
  298. # endif
  299. result = trio_pinf() / trio_pinf();
  300. # if defined(TRIO_PLATFORM_UNIX)
  301. signal(SIGFPE, signal_handler);
  302. # endif
  303. #endif
  304. }
  305. return result;
  306. }
  307. /**
  308. Check for NaN.
  309. @param number An arbitrary floating-point number.
  310. @return Boolean value indicating whether or not the number is a NaN.
  311. */
  312. TRIO_PUBLIC int
  313. trio_isnan
  314. TRIO_ARGS1((number),
  315. double number)
  316. {
  317. #if (defined(TRIO_COMPILER_SUPPORTS_C99) && defined(isnan)) \
  318. || defined(TRIO_COMPILER_SUPPORTS_UNIX95)
  319. /*
  320. * C99 defines isnan() as a macro. UNIX95 defines isnan() as a
  321. * function. This function was already present in XPG4, but this
  322. * is a bit tricky to detect with compiler defines, so we choose
  323. * the conservative approach and only use it for UNIX95.
  324. */
  325. return isnan(number);
  326. #elif defined(TRIO_COMPILER_MSVC) || defined(TRIO_COMPILER_BCB)
  327. /*
  328. * Microsoft Visual C++ and Borland C++ Builder have an _isnan()
  329. * function.
  330. */
  331. return _isnan(number) ? TRIO_TRUE : TRIO_FALSE;
  332. #elif defined(USE_IEEE_754)
  333. /*
  334. * Examine IEEE 754 bit-pattern. A NaN must have a special exponent
  335. * pattern, and a non-empty mantissa.
  336. */
  337. int has_mantissa;
  338. int is_special_quantity;
  339. is_special_quantity = trio_is_special_quantity(number, &has_mantissa);
  340. return (is_special_quantity && has_mantissa);
  341. #else
  342. /*
  343. * Fallback solution
  344. */
  345. int status;
  346. double integral, fraction;
  347. # if defined(TRIO_PLATFORM_UNIX)
  348. void (*signal_handler)(int) = signal(SIGFPE, SIG_IGN);
  349. # endif
  350. status = (/*
  351. * NaN is the only number which does not compare to itself
  352. */
  353. ((TRIO_VOLATILE double)number != (TRIO_VOLATILE double)number) ||
  354. /*
  355. * Fallback solution if NaN compares to NaN
  356. */
  357. ((number != 0.0) &&
  358. (fraction = modf(number, &integral),
  359. integral == fraction)));
  360. # if defined(TRIO_PLATFORM_UNIX)
  361. signal(SIGFPE, signal_handler);
  362. # endif
  363. return status;
  364. #endif
  365. }
  366. /**
  367. Check for infinity.
  368. @param number An arbitrary floating-point number.
  369. @return 1 if positive infinity, -1 if negative infinity, 0 otherwise.
  370. */
  371. TRIO_PUBLIC int
  372. trio_isinf
  373. TRIO_ARGS1((number),
  374. double number)
  375. {
  376. #if defined(TRIO_COMPILER_DECC) && !defined(__linux__)
  377. /*
  378. * DECC has an isinf() macro, but it works differently than that
  379. * of C99, so we use the fp_class() function instead.
  380. */
  381. return ((fp_class(number) == FP_POS_INF)
  382. ? 1
  383. : ((fp_class(number) == FP_NEG_INF) ? -1 : 0));
  384. #elif defined(isinf)
  385. /*
  386. * C99 defines isinf() as a macro.
  387. */
  388. return isinf(number)
  389. ? ((number > 0.0) ? 1 : -1)
  390. : 0;
  391. #elif defined(TRIO_COMPILER_MSVC) || defined(TRIO_COMPILER_BCB)
  392. /*
  393. * Microsoft Visual C++ and Borland C++ Builder have an _fpclass()
  394. * function that can be used to detect infinity.
  395. */
  396. return ((_fpclass(number) == _FPCLASS_PINF)
  397. ? 1
  398. : ((_fpclass(number) == _FPCLASS_NINF) ? -1 : 0));
  399. #elif defined(USE_IEEE_754)
  400. /*
  401. * Examine IEEE 754 bit-pattern. Infinity must have a special exponent
  402. * pattern, and an empty mantissa.
  403. */
  404. int has_mantissa;
  405. int is_special_quantity;
  406. is_special_quantity = trio_is_special_quantity(number, &has_mantissa);
  407. return (is_special_quantity && !has_mantissa)
  408. ? ((number < 0.0) ? -1 : 1)
  409. : 0;
  410. #else
  411. /*
  412. * Fallback solution.
  413. */
  414. int status;
  415. # if defined(TRIO_PLATFORM_UNIX)
  416. void (*signal_handler)(int) = signal(SIGFPE, SIG_IGN);
  417. # endif
  418. double infinity = trio_pinf();
  419. status = ((number == infinity)
  420. ? 1
  421. : ((number == -infinity) ? -1 : 0));
  422. # if defined(TRIO_PLATFORM_UNIX)
  423. signal(SIGFPE, signal_handler);
  424. # endif
  425. return status;
  426. #endif
  427. }
  428. #if 0
  429. /* Temporary fix - this routine is not used anywhere */
  430. /**
  431. Check for finity.
  432. @param number An arbitrary floating-point number.
  433. @return Boolean value indicating whether or not the number is a finite.
  434. */
  435. TRIO_PUBLIC int
  436. trio_isfinite
  437. TRIO_ARGS1((number),
  438. double number)
  439. {
  440. #if defined(TRIO_COMPILER_SUPPORTS_C99) && defined(isfinite)
  441. /*
  442. * C99 defines isfinite() as a macro.
  443. */
  444. return isfinite(number);
  445. #elif defined(TRIO_COMPILER_MSVC) || defined(TRIO_COMPILER_BCB)
  446. /*
  447. * Microsoft Visual C++ and Borland C++ Builder use _finite().
  448. */
  449. return _finite(number);
  450. #elif defined(USE_IEEE_754)
  451. /*
  452. * Examine IEEE 754 bit-pattern. For finity we do not care about the
  453. * mantissa.
  454. */
  455. int dummy;
  456. return (! trio_is_special_quantity(number, &dummy));
  457. #else
  458. /*
  459. * Fallback solution.
  460. */
  461. return ((trio_isinf(number) == 0) && (trio_isnan(number) == 0));
  462. #endif
  463. }
  464. #endif
  465. /*
  466. * The sign of NaN is always false
  467. */
  468. TRIO_PUBLIC int
  469. trio_fpclassify_and_signbit
  470. TRIO_ARGS2((number, is_negative),
  471. double number,
  472. int *is_negative)
  473. {
  474. #if defined(fpclassify) && defined(signbit)
  475. /*
  476. * C99 defines fpclassify() and signbit() as a macros
  477. */
  478. *is_negative = signbit(number);
  479. switch (fpclassify(number)) {
  480. case FP_NAN:
  481. return TRIO_FP_NAN;
  482. case FP_INFINITE:
  483. return TRIO_FP_INFINITE;
  484. case FP_SUBNORMAL:
  485. return TRIO_FP_SUBNORMAL;
  486. case FP_ZERO:
  487. return TRIO_FP_ZERO;
  488. default:
  489. return TRIO_FP_NORMAL;
  490. }
  491. #else
  492. # if defined(TRIO_COMPILER_DECC)
  493. /*
  494. * DECC has an fp_class() function.
  495. */
  496. # define TRIO_FPCLASSIFY(n) fp_class(n)
  497. # define TRIO_QUIET_NAN FP_QNAN
  498. # define TRIO_SIGNALLING_NAN FP_SNAN
  499. # define TRIO_POSITIVE_INFINITY FP_POS_INF
  500. # define TRIO_NEGATIVE_INFINITY FP_NEG_INF
  501. # define TRIO_POSITIVE_SUBNORMAL FP_POS_DENORM
  502. # define TRIO_NEGATIVE_SUBNORMAL FP_NEG_DENORM
  503. # define TRIO_POSITIVE_ZERO FP_POS_ZERO
  504. # define TRIO_NEGATIVE_ZERO FP_NEG_ZERO
  505. # define TRIO_POSITIVE_NORMAL FP_POS_NORM
  506. # define TRIO_NEGATIVE_NORMAL FP_NEG_NORM
  507. # elif defined(TRIO_COMPILER_MSVC) || defined(TRIO_COMPILER_BCB)
  508. /*
  509. * Microsoft Visual C++ and Borland C++ Builder have an _fpclass()
  510. * function.
  511. */
  512. # define TRIO_FPCLASSIFY(n) _fpclass(n)
  513. # define TRIO_QUIET_NAN _FPCLASS_QNAN
  514. # define TRIO_SIGNALLING_NAN _FPCLASS_SNAN
  515. # define TRIO_POSITIVE_INFINITY _FPCLASS_PINF
  516. # define TRIO_NEGATIVE_INFINITY _FPCLASS_NINF
  517. # define TRIO_POSITIVE_SUBNORMAL _FPCLASS_PD
  518. # define TRIO_NEGATIVE_SUBNORMAL _FPCLASS_ND
  519. # define TRIO_POSITIVE_ZERO _FPCLASS_PZ
  520. # define TRIO_NEGATIVE_ZERO _FPCLASS_NZ
  521. # define TRIO_POSITIVE_NORMAL _FPCLASS_PN
  522. # define TRIO_NEGATIVE_NORMAL _FPCLASS_NN
  523. # elif defined(FP_PLUS_NORM)
  524. /*
  525. * HP-UX 9.x and 10.x have an fpclassify() function, that is different
  526. * from the C99 fpclassify() macro supported on HP-UX 11.x.
  527. *
  528. * AIX has class() for C, and _class() for C++, which returns the
  529. * same values as the HP-UX fpclassify() function.
  530. */
  531. # if defined(TRIO_PLATFORM_AIX)
  532. # if defined(__cplusplus)
  533. # define TRIO_FPCLASSIFY(n) _class(n)
  534. # else
  535. # define TRIO_FPCLASSIFY(n) class(n)
  536. # endif
  537. # else
  538. # define TRIO_FPCLASSIFY(n) fpclassify(n)
  539. # endif
  540. # define TRIO_QUIET_NAN FP_QNAN
  541. # define TRIO_SIGNALLING_NAN FP_SNAN
  542. # define TRIO_POSITIVE_INFINITY FP_PLUS_INF
  543. # define TRIO_NEGATIVE_INFINITY FP_MINUS_INF
  544. # define TRIO_POSITIVE_SUBNORMAL FP_PLUS_DENORM
  545. # define TRIO_NEGATIVE_SUBNORMAL FP_MINUS_DENORM
  546. # define TRIO_POSITIVE_ZERO FP_PLUS_ZERO
  547. # define TRIO_NEGATIVE_ZERO FP_MINUS_ZERO
  548. # define TRIO_POSITIVE_NORMAL FP_PLUS_NORM
  549. # define TRIO_NEGATIVE_NORMAL FP_MINUS_NORM
  550. # endif
  551. # if defined(TRIO_FPCLASSIFY)
  552. switch (TRIO_FPCLASSIFY(number)) {
  553. case TRIO_QUIET_NAN:
  554. case TRIO_SIGNALLING_NAN:
  555. *is_negative = TRIO_FALSE; /* NaN has no sign */
  556. return TRIO_FP_NAN;
  557. case TRIO_POSITIVE_INFINITY:
  558. *is_negative = TRIO_FALSE;
  559. return TRIO_FP_INFINITE;
  560. case TRIO_NEGATIVE_INFINITY:
  561. *is_negative = TRIO_TRUE;
  562. return TRIO_FP_INFINITE;
  563. case TRIO_POSITIVE_SUBNORMAL:
  564. *is_negative = TRIO_FALSE;
  565. return TRIO_FP_SUBNORMAL;
  566. case TRIO_NEGATIVE_SUBNORMAL:
  567. *is_negative = TRIO_TRUE;
  568. return TRIO_FP_SUBNORMAL;
  569. case TRIO_POSITIVE_ZERO:
  570. *is_negative = TRIO_FALSE;
  571. return TRIO_FP_ZERO;
  572. case TRIO_NEGATIVE_ZERO:
  573. *is_negative = TRIO_TRUE;
  574. return TRIO_FP_ZERO;
  575. case TRIO_POSITIVE_NORMAL:
  576. *is_negative = TRIO_FALSE;
  577. return TRIO_FP_NORMAL;
  578. case TRIO_NEGATIVE_NORMAL:
  579. *is_negative = TRIO_TRUE;
  580. return TRIO_FP_NORMAL;
  581. default:
  582. /* Just in case... */
  583. *is_negative = (number < 0.0);
  584. return TRIO_FP_NORMAL;
  585. }
  586. # else
  587. /*
  588. * Fallback solution.
  589. */
  590. int rc;
  591. if (number == 0.0) {
  592. /*
  593. * In IEEE 754 the sign of zero is ignored in comparisons, so we
  594. * have to handle this as a special case by examining the sign bit
  595. * directly.
  596. */
  597. # if defined(USE_IEEE_754)
  598. *is_negative = trio_is_negative(number);
  599. # else
  600. *is_negative = TRIO_FALSE; /* FIXME */
  601. # endif
  602. return TRIO_FP_ZERO;
  603. }
  604. if (trio_isnan(number)) {
  605. *is_negative = TRIO_FALSE;
  606. return TRIO_FP_NAN;
  607. }
  608. if ((rc = trio_isinf(number))) {
  609. *is_negative = (rc == -1);
  610. return TRIO_FP_INFINITE;
  611. }
  612. if ((number > 0.0) && (number < DBL_MIN)) {
  613. *is_negative = TRIO_FALSE;
  614. return TRIO_FP_SUBNORMAL;
  615. }
  616. if ((number < 0.0) && (number > -DBL_MIN)) {
  617. *is_negative = TRIO_TRUE;
  618. return TRIO_FP_SUBNORMAL;
  619. }
  620. *is_negative = (number < 0.0);
  621. return TRIO_FP_NORMAL;
  622. # endif
  623. #endif
  624. }
  625. /**
  626. Examine the sign of a number.
  627. @param number An arbitrary floating-point number.
  628. @return Boolean value indicating whether or not the number has the
  629. sign bit set (i.e. is negative).
  630. */
  631. TRIO_PUBLIC int
  632. trio_signbit
  633. TRIO_ARGS1((number),
  634. double number)
  635. {
  636. int is_negative;
  637. (void)trio_fpclassify_and_signbit(number, &is_negative);
  638. return is_negative;
  639. }
  640. #if 0
  641. /* Temporary fix - this routine is not used in libxml */
  642. /**
  643. Examine the class of a number.
  644. @param number An arbitrary floating-point number.
  645. @return Enumerable value indicating the class of @p number
  646. */
  647. TRIO_PUBLIC int
  648. trio_fpclassify
  649. TRIO_ARGS1((number),
  650. double number)
  651. {
  652. int dummy;
  653. return trio_fpclassify_and_signbit(number, &dummy);
  654. }
  655. #endif
  656. /** @} SpecialQuantities */
  657. /*************************************************************************
  658. * For test purposes.
  659. *
  660. * Add the following compiler option to include this test code.
  661. *
  662. * Unix : -DSTANDALONE
  663. * VMS : /DEFINE=(STANDALONE)
  664. */
  665. #if defined(STANDALONE)
  666. # include <stdio.h>
  667. static TRIO_CONST char *
  668. getClassification
  669. TRIO_ARGS1((type),
  670. int type)
  671. {
  672. switch (type) {
  673. case TRIO_FP_INFINITE:
  674. return "FP_INFINITE";
  675. case TRIO_FP_NAN:
  676. return "FP_NAN";
  677. case TRIO_FP_NORMAL:
  678. return "FP_NORMAL";
  679. case TRIO_FP_SUBNORMAL:
  680. return "FP_SUBNORMAL";
  681. case TRIO_FP_ZERO:
  682. return "FP_ZERO";
  683. default:
  684. return "FP_UNKNOWN";
  685. }
  686. }
  687. static void
  688. print_class
  689. TRIO_ARGS2((prefix, number),
  690. TRIO_CONST char *prefix,
  691. double number)
  692. {
  693. printf("%-6s: %s %-15s %g\n",
  694. prefix,
  695. trio_signbit(number) ? "-" : "+",
  696. getClassification(TRIO_FPCLASSIFY(number)),
  697. number);
  698. }
  699. int main(TRIO_NOARGS)
  700. {
  701. double my_nan;
  702. double my_pinf;
  703. double my_ninf;
  704. # if defined(TRIO_PLATFORM_UNIX)
  705. void (*signal_handler) TRIO_PROTO((int));
  706. # endif
  707. my_nan = trio_nan();
  708. my_pinf = trio_pinf();
  709. my_ninf = trio_ninf();
  710. print_class("Nan", my_nan);
  711. print_class("PInf", my_pinf);
  712. print_class("NInf", my_ninf);
  713. print_class("PZero", 0.0);
  714. print_class("NZero", -0.0);
  715. print_class("PNorm", 1.0);
  716. print_class("NNorm", -1.0);
  717. print_class("PSub", 1.01e-307 - 1.00e-307);
  718. print_class("NSub", 1.00e-307 - 1.01e-307);
  719. printf("NaN : %4g 0x%02x%02x%02x%02x%02x%02x%02x%02x (%2d, %2d)\n",
  720. my_nan,
  721. ((unsigned char *)&my_nan)[0],
  722. ((unsigned char *)&my_nan)[1],
  723. ((unsigned char *)&my_nan)[2],
  724. ((unsigned char *)&my_nan)[3],
  725. ((unsigned char *)&my_nan)[4],
  726. ((unsigned char *)&my_nan)[5],
  727. ((unsigned char *)&my_nan)[6],
  728. ((unsigned char *)&my_nan)[7],
  729. trio_isnan(my_nan), trio_isinf(my_nan));
  730. printf("PInf: %4g 0x%02x%02x%02x%02x%02x%02x%02x%02x (%2d, %2d)\n",
  731. my_pinf,
  732. ((unsigned char *)&my_pinf)[0],
  733. ((unsigned char *)&my_pinf)[1],
  734. ((unsigned char *)&my_pinf)[2],
  735. ((unsigned char *)&my_pinf)[3],
  736. ((unsigned char *)&my_pinf)[4],
  737. ((unsigned char *)&my_pinf)[5],
  738. ((unsigned char *)&my_pinf)[6],
  739. ((unsigned char *)&my_pinf)[7],
  740. trio_isnan(my_pinf), trio_isinf(my_pinf));
  741. printf("NInf: %4g 0x%02x%02x%02x%02x%02x%02x%02x%02x (%2d, %2d)\n",
  742. my_ninf,
  743. ((unsigned char *)&my_ninf)[0],
  744. ((unsigned char *)&my_ninf)[1],
  745. ((unsigned char *)&my_ninf)[2],
  746. ((unsigned char *)&my_ninf)[3],
  747. ((unsigned char *)&my_ninf)[4],
  748. ((unsigned char *)&my_ninf)[5],
  749. ((unsigned char *)&my_ninf)[6],
  750. ((unsigned char *)&my_ninf)[7],
  751. trio_isnan(my_ninf), trio_isinf(my_ninf));
  752. # if defined(TRIO_PLATFORM_UNIX)
  753. signal_handler = signal(SIGFPE, SIG_IGN);
  754. # endif
  755. my_pinf = DBL_MAX + DBL_MAX;
  756. my_ninf = -my_pinf;
  757. my_nan = my_pinf / my_pinf;
  758. # if defined(TRIO_PLATFORM_UNIX)
  759. signal(SIGFPE, signal_handler);
  760. # endif
  761. printf("NaN : %4g 0x%02x%02x%02x%02x%02x%02x%02x%02x (%2d, %2d)\n",
  762. my_nan,
  763. ((unsigned char *)&my_nan)[0],
  764. ((unsigned char *)&my_nan)[1],
  765. ((unsigned char *)&my_nan)[2],
  766. ((unsigned char *)&my_nan)[3],
  767. ((unsigned char *)&my_nan)[4],
  768. ((unsigned char *)&my_nan)[5],
  769. ((unsigned char *)&my_nan)[6],
  770. ((unsigned char *)&my_nan)[7],
  771. trio_isnan(my_nan), trio_isinf(my_nan));
  772. printf("PInf: %4g 0x%02x%02x%02x%02x%02x%02x%02x%02x (%2d, %2d)\n",
  773. my_pinf,
  774. ((unsigned char *)&my_pinf)[0],
  775. ((unsigned char *)&my_pinf)[1],
  776. ((unsigned char *)&my_pinf)[2],
  777. ((unsigned char *)&my_pinf)[3],
  778. ((unsigned char *)&my_pinf)[4],
  779. ((unsigned char *)&my_pinf)[5],
  780. ((unsigned char *)&my_pinf)[6],
  781. ((unsigned char *)&my_pinf)[7],
  782. trio_isnan(my_pinf), trio_isinf(my_pinf));
  783. printf("NInf: %4g 0x%02x%02x%02x%02x%02x%02x%02x%02x (%2d, %2d)\n",
  784. my_ninf,
  785. ((unsigned char *)&my_ninf)[0],
  786. ((unsigned char *)&my_ninf)[1],
  787. ((unsigned char *)&my_ninf)[2],
  788. ((unsigned char *)&my_ninf)[3],
  789. ((unsigned char *)&my_ninf)[4],
  790. ((unsigned char *)&my_ninf)[5],
  791. ((unsigned char *)&my_ninf)[6],
  792. ((unsigned char *)&my_ninf)[7],
  793. trio_isnan(my_ninf), trio_isinf(my_ninf));
  794. return 0;
  795. }
  796. #endif