1703 lines · c
1/*2 * mathtest.c - test rig for mathlib3 *4 * Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.5 * See https://llvm.org/LICENSE.txt for license information.6 * SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception7 */8 9#include <assert.h>10#include <stdio.h>11#include <stdlib.h>12#include <string.h>13#include <setjmp.h>14#include <ctype.h>15#include <math.h>16#include <errno.h>17#include <limits.h>18#include <fenv.h>19#include "mathlib.h"20 21#ifndef math_errhandling22# define math_errhandling 023#endif24 25#ifdef __cplusplus26 #define EXTERN_C extern "C"27#else28 #define EXTERN_C extern29#endif30 31#ifndef TRUE32#define TRUE 133#endif34#ifndef FALSE35#define FALSE 036#endif37 38#ifdef IMPORT_SYMBOL39#define STR2(x) #x40#define STR(x) STR2(x)41_Pragma(STR(import IMPORT_SYMBOL))42#endif43 44int dmsd, dlsd;45int quiet = 0;46int doround = 0;47unsigned statusmask = FE_ALL_EXCEPT;48 49#define EXTRABITS (12)50#define ULPUNIT (1<<EXTRABITS)51 52typedef int (*test) (void);53 54/*55 struct to hold info about a function (which could actually be a macro)56*/57typedef struct {58 enum {59 t_func, t_macro60 } type;61 enum {62 at_d, at_s, /* double or single precision float */63 at_d2, at_s2, /* same, but taking two args */64 at_di, at_si, /* double/single and an int */65 at_dip, at_sip, /* double/single and an int ptr */66 at_ddp, at_ssp, /* d/s and a d/s ptr */67 at_dc, at_sc, /* double or single precision complex */68 at_dc2, at_sc2 /* same, but taking two args */69 } argtype;70 enum {71 rt_d, rt_s, rt_i, /* double, single, int */72 rt_dc, rt_sc, /* double, single precision complex */73 rt_d2, rt_s2 /* also use res2 */74 } rettype;75 union {76 void* ptr;77 double (*d_d_ptr)(double);78 float (*s_s_ptr)(float);79 int (*d_i_ptr)(double);80 int (*s_i_ptr)(float);81 double (*d2_d_ptr)(double, double);82 float (*s2_s_ptr)(float, float);83 double (*di_d_ptr)(double,int);84 float (*si_s_ptr)(float,int);85 double (*dip_d_ptr)(double,int*);86 float (*sip_s_ptr)(float,int*);87 double (*ddp_d_ptr)(double,double*);88 float (*ssp_s_ptr)(float,float*);89 } func;90 enum {91 m_none,92 m_isfinite, m_isfinitef,93 m_isgreater, m_isgreaterequal,94 m_isgreaterequalf, m_isgreaterf,95 m_isinf, m_isinff,96 m_isless, m_islessequal,97 m_islessequalf, m_islessf,98 m_islessgreater, m_islessgreaterf,99 m_isnan, m_isnanf,100 m_isnormal, m_isnormalf,101 m_isunordered, m_isunorderedf,102 m_fpclassify, m_fpclassifyf,103 m_signbit, m_signbitf,104 /* not actually a macro, but makes things easier */105 m_rred, m_rredf,106 m_cadd, m_csub, m_cmul, m_cdiv,107 m_caddf, m_csubf, m_cmulf, m_cdivf108 } macro_name; /* only used if a macro/something that can't be done using func */109 long long tolerance;110 const char* name;111} test_func;112 113/* used in qsort */114int compare_tfuncs(const void* a, const void* b) {115 return strcmp(((test_func*)a)->name, ((test_func*)b)->name);116}117 118int is_double_argtype(int argtype) {119 switch(argtype) {120 case at_d:121 case at_d2:122 case at_dc:123 case at_dc2:124 return 1;125 default:126 return 0;127 }128}129 130int is_single_argtype(int argtype) {131 switch(argtype) {132 case at_s:133 case at_s2:134 case at_sc:135 case at_sc2:136 return 1;137 default:138 return 0;139 }140}141 142int is_double_rettype(int rettype) {143 switch(rettype) {144 case rt_d:145 case rt_dc:146 case rt_d2:147 return 1;148 default:149 return 0;150 }151}152 153int is_single_rettype(int rettype) {154 switch(rettype) {155 case rt_s:156 case rt_sc:157 case rt_s2:158 return 1;159 default:160 return 0;161 }162}163 164int is_complex_argtype(int argtype) {165 switch(argtype) {166 case at_dc:167 case at_sc:168 case at_dc2:169 case at_sc2:170 return 1;171 default:172 return 0;173 }174}175 176int is_complex_rettype(int rettype) {177 switch(rettype) {178 case rt_dc:179 case rt_sc:180 return 1;181 default:182 return 0;183 }184}185 186/*187 * Special-case flags indicating that some functions' error188 * tolerance handling is more complicated than a fixed relative189 * error bound.190 */191#define ABSLOWERBOUND 0x4000000000000000LL192#define PLUSMINUSPIO2 0x1000000000000000LL193 194#define ARM_PREFIX(x) x195 196#define TFUNC(arg,ret,name,tolerance) { t_func, arg, ret, (void*)&name, m_none, tolerance, #name }197#define TFUNCARM(arg,ret,name,tolerance) { t_func, arg, ret, (void*)& ARM_PREFIX(name), m_none, tolerance, #name }198#define MFUNC(arg,ret,name,tolerance) { t_macro, arg, ret, NULL, m_##name, tolerance, #name }199 200/* sincosf wrappers for easier testing. */201static float sincosf_sinf(float x) { float s,c; sincosf(x, &s, &c); return s; }202static float sincosf_cosf(float x) { float s,c; sincosf(x, &s, &c); return c; }203 204test_func tfuncs[] = {205 /* trigonometric */206 TFUNC(at_d,rt_d, acos, 4*ULPUNIT),207 TFUNC(at_d,rt_d, asin, 4*ULPUNIT),208 TFUNC(at_d,rt_d, atan, 4*ULPUNIT),209 TFUNC(at_d2,rt_d, atan2, 4*ULPUNIT),210 211 TFUNC(at_d,rt_d, tan, 2*ULPUNIT),212 TFUNC(at_d,rt_d, sin, 2*ULPUNIT),213 TFUNC(at_d,rt_d, cos, 2*ULPUNIT),214 215 TFUNC(at_s,rt_s, acosf, 4*ULPUNIT),216 TFUNC(at_s,rt_s, asinf, 4*ULPUNIT),217 TFUNC(at_s,rt_s, atanf, 4*ULPUNIT),218 TFUNC(at_s2,rt_s, atan2f, 4*ULPUNIT),219 TFUNCARM(at_s,rt_s, tanf, 4*ULPUNIT),220 TFUNCARM(at_s,rt_s, sinf, 3*ULPUNIT/4),221 TFUNCARM(at_s,rt_s, cosf, 3*ULPUNIT/4),222 TFUNCARM(at_s,rt_s, sincosf_sinf, 3*ULPUNIT/4),223 TFUNCARM(at_s,rt_s, sincosf_cosf, 3*ULPUNIT/4),224 225 /* hyperbolic */226 TFUNC(at_d, rt_d, atanh, 4*ULPUNIT),227 TFUNC(at_d, rt_d, asinh, 4*ULPUNIT),228 TFUNC(at_d, rt_d, acosh, 4*ULPUNIT),229 TFUNC(at_d,rt_d, tanh, 4*ULPUNIT),230 TFUNC(at_d,rt_d, sinh, 4*ULPUNIT),231 TFUNC(at_d,rt_d, cosh, 4*ULPUNIT),232 233 TFUNC(at_s, rt_s, atanhf, 4*ULPUNIT),234 TFUNC(at_s, rt_s, asinhf, 4*ULPUNIT),235 TFUNC(at_s, rt_s, acoshf, 4*ULPUNIT),236 TFUNC(at_s,rt_s, tanhf, 4*ULPUNIT),237 TFUNC(at_s,rt_s, sinhf, 4*ULPUNIT),238 TFUNC(at_s,rt_s, coshf, 4*ULPUNIT),239 240 /* exponential and logarithmic */241 TFUNC(at_d,rt_d, log, 3*ULPUNIT/4),242 TFUNC(at_d,rt_d, log10, 3*ULPUNIT),243 TFUNC(at_d,rt_d, log2, 3*ULPUNIT/4),244 TFUNC(at_d,rt_d, log1p, 2*ULPUNIT),245 TFUNC(at_d,rt_d, exp, 3*ULPUNIT/4),246 TFUNC(at_d,rt_d, exp2, 3*ULPUNIT/4),247 TFUNC(at_d,rt_d, expm1, ULPUNIT),248 TFUNCARM(at_s,rt_s, logf, ULPUNIT),249 TFUNC(at_s,rt_s, log10f, 3*ULPUNIT),250 TFUNCARM(at_s,rt_s, log2f, ULPUNIT),251 TFUNC(at_s,rt_s, log1pf, 2*ULPUNIT),252 TFUNCARM(at_s,rt_s, expf, 3*ULPUNIT/4),253 TFUNCARM(at_s,rt_s, exp2f, 3*ULPUNIT/4),254 TFUNC(at_s,rt_s, expm1f, ULPUNIT),255 256 /* power */257 TFUNC(at_d2,rt_d, pow, 3*ULPUNIT/4),258 TFUNC(at_d,rt_d, sqrt, ULPUNIT/2),259 TFUNC(at_d,rt_d, cbrt, 2*ULPUNIT),260 TFUNC(at_d2, rt_d, hypot, 4*ULPUNIT),261 262 TFUNCARM(at_s2,rt_s, powf, ULPUNIT),263 TFUNC(at_s,rt_s, sqrtf, ULPUNIT/2),264 TFUNC(at_s,rt_s, cbrtf, 2*ULPUNIT),265 TFUNC(at_s2, rt_s, hypotf, 4*ULPUNIT),266 267 /* error function */268 TFUNC(at_d,rt_d, erf, 16*ULPUNIT),269 TFUNC(at_s,rt_s, erff, 16*ULPUNIT),270 TFUNC(at_d,rt_d, erfc, 16*ULPUNIT),271 TFUNC(at_s,rt_s, erfcf, 16*ULPUNIT),272 273 /* gamma functions */274 TFUNC(at_d,rt_d, tgamma, 16*ULPUNIT),275 TFUNC(at_s,rt_s, tgammaf, 16*ULPUNIT),276 TFUNC(at_d,rt_d, lgamma, 16*ULPUNIT | ABSLOWERBOUND),277 TFUNC(at_s,rt_s, lgammaf, 16*ULPUNIT | ABSLOWERBOUND),278 279 TFUNC(at_d,rt_d, ceil, 0),280 TFUNC(at_s,rt_s, ceilf, 0),281 TFUNC(at_d2,rt_d, copysign, 0),282 TFUNC(at_s2,rt_s, copysignf, 0),283 TFUNC(at_d,rt_d, floor, 0),284 TFUNC(at_s,rt_s, floorf, 0),285 TFUNC(at_d2,rt_d, fmax, 0),286 TFUNC(at_s2,rt_s, fmaxf, 0),287 TFUNC(at_d2,rt_d, fmin, 0),288 TFUNC(at_s2,rt_s, fminf, 0),289 TFUNC(at_d2,rt_d, fmod, 0),290 TFUNC(at_s2,rt_s, fmodf, 0),291 MFUNC(at_d, rt_i, fpclassify, 0),292 MFUNC(at_s, rt_i, fpclassifyf, 0),293 TFUNC(at_dip,rt_d, frexp, 0),294 TFUNC(at_sip,rt_s, frexpf, 0),295 MFUNC(at_d, rt_i, isfinite, 0),296 MFUNC(at_s, rt_i, isfinitef, 0),297 MFUNC(at_d, rt_i, isgreater, 0),298 MFUNC(at_d, rt_i, isgreaterequal, 0),299 MFUNC(at_s, rt_i, isgreaterequalf, 0),300 MFUNC(at_s, rt_i, isgreaterf, 0),301 MFUNC(at_d, rt_i, isinf, 0),302 MFUNC(at_s, rt_i, isinff, 0),303 MFUNC(at_d, rt_i, isless, 0),304 MFUNC(at_d, rt_i, islessequal, 0),305 MFUNC(at_s, rt_i, islessequalf, 0),306 MFUNC(at_s, rt_i, islessf, 0),307 MFUNC(at_d, rt_i, islessgreater, 0),308 MFUNC(at_s, rt_i, islessgreaterf, 0),309 MFUNC(at_d, rt_i, isnan, 0),310 MFUNC(at_s, rt_i, isnanf, 0),311 MFUNC(at_d, rt_i, isnormal, 0),312 MFUNC(at_s, rt_i, isnormalf, 0),313 MFUNC(at_d, rt_i, isunordered, 0),314 MFUNC(at_s, rt_i, isunorderedf, 0),315 TFUNC(at_di,rt_d, ldexp, 0),316 TFUNC(at_si,rt_s, ldexpf, 0),317 TFUNC(at_ddp,rt_d2, modf, 0),318 TFUNC(at_ssp,rt_s2, modff, 0),319#ifndef BIGRANGERED320 MFUNC(at_d, rt_d, rred, 2*ULPUNIT),321#else322 MFUNC(at_d, rt_d, m_rred, ULPUNIT),323#endif324 MFUNC(at_d, rt_i, signbit, 0),325 MFUNC(at_s, rt_i, signbitf, 0),326};327 328/*329 * keywords are: func size op1 op2 result res2 errno op1r op1i op2r op2i resultr resulti330 * also we ignore: wrongresult wrongres2 wrongerrno331 * op1 equivalent to op1r, same with op2 and result332 */333 334typedef struct {335 test_func *func;336 unsigned op1r[2]; /* real part, also used for non-complex numbers */337 unsigned op1i[2]; /* imaginary part */338 unsigned op2r[2];339 unsigned op2i[2];340 unsigned resultr[3];341 unsigned resulti[3];342 enum {343 rc_none, rc_zero, rc_infinity, rc_nan, rc_finite344 } resultc; /* special complex results, rc_none means use resultr and resulti as normal */345 unsigned res2[2];346 unsigned status; /* IEEE status return, if any */347 unsigned maybestatus; /* for optional status, or allowance for spurious */348 int nresult; /* number of result words */349 int in_err, in_err_limit;350 int err;351 int maybeerr;352 int valid;353 int comment;354 int random;355} testdetail;356 357enum { /* keywords */358 k_errno, k_errno_in, k_error, k_func, k_maybeerror, k_maybestatus, k_op1, k_op1i, k_op1r, k_op2, k_op2i, k_op2r,359 k_random, k_res2, k_result, k_resultc, k_resulti, k_resultr, k_status,360 k_wrongres2, k_wrongresult, k_wrongstatus, k_wrongerrno361};362char *keywords[] = {363 "errno", "errno_in", "error", "func", "maybeerror", "maybestatus", "op1", "op1i", "op1r", "op2", "op2i", "op2r",364 "random", "res2", "result", "resultc", "resulti", "resultr", "status",365 "wrongres2", "wrongresult", "wrongstatus", "wrongerrno"366};367 368enum {369 e_0, e_EDOM, e_ERANGE,370 371 /*372 * This enum makes sure that we have the right number of errnos in the373 * errno[] array374 */375 e_number_of_errnos376};377char *errnos[] = {378 "0", "EDOM", "ERANGE"379};380 381enum {382 e_none, e_divbyzero, e_domain, e_overflow, e_underflow383};384char *errors[] = {385 "0", "divbyzero", "domain", "overflow", "underflow"386};387 388static int verbose, fo, strict;389 390/* state toggled by random=on / random=off */391static int randomstate;392 393/* Canonify a double NaN: SNaNs all become 7FF00000.00000001 and QNaNs394 * all become 7FF80000.00000001 */395void canon_dNaN(unsigned a[2]) {396 if ((a[0] & 0x7FF00000) != 0x7FF00000)397 return; /* not Inf or NaN */398 if (!(a[0] & 0xFFFFF) && !a[1])399 return; /* Inf */400 a[0] &= 0x7FF80000; /* canonify top word */401 a[1] = 0x00000001; /* canonify bottom word */402}403 404/* Canonify a single NaN: SNaNs all become 7F800001 and QNaNs405 * all become 7FC00001. Returns classification of the NaN. */406void canon_sNaN(unsigned a[1]) {407 if ((a[0] & 0x7F800000) != 0x7F800000)408 return; /* not Inf or NaN */409 if (!(a[0] & 0x7FFFFF))410 return; /* Inf */411 a[0] &= 0x7FC00000; /* canonify most bits */412 a[0] |= 0x00000001; /* canonify bottom bit */413}414 415/*416 * Detect difficult operands for FO mode.417 */418int is_dhard(unsigned a[2])419{420 if ((a[0] & 0x7FF00000) == 0x7FF00000)421 return TRUE; /* inf or NaN */422 if ((a[0] & 0x7FF00000) == 0 &&423 ((a[0] & 0x7FFFFFFF) | a[1]) != 0)424 return TRUE; /* denormal */425 return FALSE;426}427int is_shard(unsigned a[1])428{429 if ((a[0] & 0x7F800000) == 0x7F800000)430 return TRUE; /* inf or NaN */431 if ((a[0] & 0x7F800000) == 0 &&432 (a[0] & 0x7FFFFFFF) != 0)433 return TRUE; /* denormal */434 return FALSE;435}436 437/*438 * Normalise all zeroes into +0, for FO mode.439 */440void dnormzero(unsigned a[2])441{442 if (a[0] == 0x80000000 && a[1] == 0)443 a[0] = 0;444}445void snormzero(unsigned a[1])446{447 if (a[0] == 0x80000000)448 a[0] = 0;449}450 451static int find(char *word, char **array, int asize) {452 int i, j;453 454 asize /= sizeof(char *);455 456 i = -1; j = asize; /* strictly between i and j */457 while (j-i > 1) {458 int k = (i+j) / 2;459 int c = strcmp(word, array[k]);460 if (c > 0)461 i = k;462 else if (c < 0)463 j = k;464 else /* found it! */465 return k;466 }467 return -1; /* not found */468}469 470static test_func* find_testfunc(char *word) {471 int i, j, asize;472 473 asize = sizeof(tfuncs)/sizeof(test_func);474 475 i = -1; j = asize; /* strictly between i and j */476 while (j-i > 1) {477 int k = (i+j) / 2;478 int c = strcmp(word, tfuncs[k].name);479 if (c > 0)480 i = k;481 else if (c < 0)482 j = k;483 else /* found it! */484 return tfuncs + k;485 }486 return NULL; /* not found */487}488 489static long long calc_error(unsigned a[2], unsigned b[3], int shift, int rettype) {490 unsigned r0, r1, r2;491 int sign, carry;492 long long result;493 494 /*495 * If either number is infinite, require exact equality. If496 * either number is NaN, require that both are NaN. If either497 * of these requirements is broken, return INT_MAX.498 */499 if (is_double_rettype(rettype)) {500 if ((a[0] & 0x7FF00000) == 0x7FF00000 ||501 (b[0] & 0x7FF00000) == 0x7FF00000) {502 if (((a[0] & 0x800FFFFF) || a[1]) &&503 ((b[0] & 0x800FFFFF) || b[1]) &&504 (a[0] & 0x7FF00000) == 0x7FF00000 &&505 (b[0] & 0x7FF00000) == 0x7FF00000)506 return 0; /* both NaN - OK */507 if (!((a[0] & 0xFFFFF) || a[1]) &&508 !((b[0] & 0xFFFFF) || b[1]) &&509 a[0] == b[0])510 return 0; /* both same sign of Inf - OK */511 return LLONG_MAX;512 }513 } else {514 if ((a[0] & 0x7F800000) == 0x7F800000 ||515 (b[0] & 0x7F800000) == 0x7F800000) {516 if ((a[0] & 0x807FFFFF) &&517 (b[0] & 0x807FFFFF) &&518 (a[0] & 0x7F800000) == 0x7F800000 &&519 (b[0] & 0x7F800000) == 0x7F800000)520 return 0; /* both NaN - OK */521 if (!(a[0] & 0x7FFFFF) &&522 !(b[0] & 0x7FFFFF) &&523 a[0] == b[0])524 return 0; /* both same sign of Inf - OK */525 return LLONG_MAX;526 }527 }528 529 /*530 * Both finite. Return INT_MAX if the signs differ.531 */532 if ((a[0] ^ b[0]) & 0x80000000)533 return LLONG_MAX;534 535 /*536 * Now it's just straight multiple-word subtraction.537 */538 if (is_double_rettype(rettype)) {539 r2 = -b[2]; carry = (r2 == 0);540 r1 = a[1] + ~b[1] + carry; carry = (r1 < a[1] || (carry && r1 == a[1]));541 r0 = a[0] + ~b[0] + carry;542 } else {543 r2 = -b[1]; carry = (r2 == 0);544 r1 = a[0] + ~b[0] + carry; carry = (r1 < a[0] || (carry && r1 == a[0]));545 r0 = ~0 + carry;546 }547 548 /*549 * Forgive larger errors in specialised cases.550 */551 if (shift > 0) {552 if (shift > 32*3)553 return 0; /* all errors are forgiven! */554 while (shift >= 32) {555 r2 = r1;556 r1 = r0;557 r0 = -(r0 >> 31);558 shift -= 32;559 }560 561 if (shift > 0) {562 r2 = (r2 >> shift) | (r1 << (32-shift));563 r1 = (r1 >> shift) | (r0 << (32-shift));564 r0 = (r0 >> shift) | ((-(r0 >> 31)) << (32-shift));565 }566 }567 568 if (r0 & 0x80000000) {569 sign = 1;570 r2 = ~r2; carry = (r2 == 0);571 r1 = 0 + ~r1 + carry; carry = (carry && (r2 == 0));572 r0 = 0 + ~r0 + carry;573 } else {574 sign = 0;575 }576 577 if (r0 >= (1LL<<(31-EXTRABITS)))578 return LLONG_MAX; /* many ulps out */579 580 result = (r2 >> (32-EXTRABITS)) & (ULPUNIT-1);581 result |= r1 << EXTRABITS;582 result |= (long long)r0 << (32+EXTRABITS);583 if (sign)584 result = -result;585 return result;586}587 588/* special named operands */589 590typedef struct {591 unsigned op1, op2;592 char* name;593} special_op;594 595static special_op special_ops_double[] = {596 {0x00000000,0x00000000,"0"},597 {0x3FF00000,0x00000000,"1"},598 {0x7FF00000,0x00000000,"inf"},599 {0x7FF80000,0x00000001,"qnan"},600 {0x7FF00000,0x00000001,"snan"},601 {0x3ff921fb,0x54442d18,"pi2"},602 {0x400921fb,0x54442d18,"pi"},603 {0x3fe921fb,0x54442d18,"pi4"},604 {0x4002d97c,0x7f3321d2,"3pi4"},605};606 607static special_op special_ops_float[] = {608 {0x00000000,0,"0"},609 {0x3f800000,0,"1"},610 {0x7f800000,0,"inf"},611 {0x7fc00000,0,"qnan"},612 {0x7f800001,0,"snan"},613 {0x3fc90fdb,0,"pi2"},614 {0x40490fdb,0,"pi"},615 {0x3f490fdb,0,"pi4"},616 {0x4016cbe4,0,"3pi4"},617};618 619/*620 This is what is returned by the below functions.621 We need it to handle the sign of the number622*/623static special_op tmp_op = {0,0,0};624 625special_op* find_special_op_from_op(unsigned op1, unsigned op2, int is_double) {626 int i;627 special_op* sop;628 if(is_double) {629 sop = special_ops_double;630 } else {631 sop = special_ops_float;632 }633 for(i = 0; i < sizeof(special_ops_double)/sizeof(special_op); i++) {634 if(sop->op1 == (op1&0x7fffffff) && sop->op2 == op2) {635 if(tmp_op.name) free(tmp_op.name);636 tmp_op.name = malloc(strlen(sop->name)+2);637 if(op1>>31) {638 sprintf(tmp_op.name,"-%s",sop->name);639 } else {640 strcpy(tmp_op.name,sop->name);641 }642 return &tmp_op;643 }644 sop++;645 }646 return NULL;647}648 649special_op* find_special_op_from_name(const char* name, int is_double) {650 int i, neg=0;651 special_op* sop;652 if(is_double) {653 sop = special_ops_double;654 } else {655 sop = special_ops_float;656 }657 if(*name=='-') {658 neg=1;659 name++;660 } else if(*name=='+') {661 name++;662 }663 for(i = 0; i < sizeof(special_ops_double)/sizeof(special_op); i++) {664 if(0 == strcmp(name,sop->name)) {665 tmp_op.op1 = sop->op1;666 if(neg) {667 tmp_op.op1 |= 0x80000000;668 }669 tmp_op.op2 = sop->op2;670 return &tmp_op;671 }672 sop++;673 }674 return NULL;675}676 677/*678 helper function for the below679 type=0 for single, 1 for double, 2 for no sop680*/681int do_op(char* q, unsigned* op, const char* name, int num, int sop_type) {682 int i;683 int n=num;684 special_op* sop = NULL;685 for(i = 0; i < num; i++) {686 op[i] = 0;687 }688 if(sop_type<2) {689 sop = find_special_op_from_name(q,sop_type);690 }691 if(sop != NULL) {692 op[0] = sop->op1;693 op[1] = sop->op2;694 } else {695 switch(num) {696 case 1: n = sscanf(q, "%x", &op[0]); break;697 case 2: n = sscanf(q, "%x.%x", &op[0], &op[1]); break;698 case 3: n = sscanf(q, "%x.%x.%x", &op[0], &op[1], &op[2]); break;699 default: return -1;700 }701 }702 if (verbose) {703 printf("%s=",name);704 for (i = 0; (i < n); ++i) printf("%x.", op[i]);705 printf(" (n=%d)\n", n);706 }707 return n;708}709 710testdetail parsetest(char *testbuf, testdetail oldtest) {711 char *p; /* Current part of line: Option name */712 char *q; /* Current part of line: Option value */713 testdetail ret; /* What we return */714 int k; /* Function enum from k_* */715 int n; /* Used as returns for scanfs */716 int argtype=2, rettype=2; /* for do_op */717 718 /* clear ret */719 memset(&ret, 0, sizeof(ret));720 721 if (verbose) printf("Parsing line: %s\n", testbuf);722 while (*testbuf && isspace(*testbuf)) testbuf++;723 if (testbuf[0] == ';' || testbuf[0] == '#' || testbuf[0] == '!' ||724 testbuf[0] == '>' || testbuf[0] == '\0') {725 ret.comment = 1;726 if (verbose) printf("Line is a comment\n");727 return ret;728 }729 ret.comment = 0;730 731 if (*testbuf == '+') {732 if (oldtest.valid) {733 ret = oldtest; /* structure copy */734 } else {735 fprintf(stderr, "copy from invalid: ignored\n");736 }737 testbuf++;738 }739 740 ret.random = randomstate;741 742 ret.in_err = 0;743 ret.in_err_limit = e_number_of_errnos;744 745 p = strtok(testbuf, " \t");746 while (p != NULL) {747 q = strchr(p, '=');748 if (!q)749 goto balderdash;750 *q++ = '\0';751 k = find(p, keywords, sizeof(keywords));752 switch (k) {753 case k_random:754 randomstate = (!strcmp(q, "on"));755 ret.comment = 1;756 return ret; /* otherwise ignore this line */757 case k_func:758 if (verbose) printf("func=%s ", q);759 //ret.func = find(q, funcs, sizeof(funcs));760 ret.func = find_testfunc(q);761 if (ret.func == NULL)762 {763 if (verbose) printf("(id=unknown)\n");764 goto balderdash;765 }766 if(is_single_argtype(ret.func->argtype))767 argtype = 0;768 else if(is_double_argtype(ret.func->argtype))769 argtype = 1;770 if(is_single_rettype(ret.func->rettype))771 rettype = 0;772 else if(is_double_rettype(ret.func->rettype))773 rettype = 1;774 //ret.size = sizes[ret.func];775 if (verbose) printf("(name=%s) (size=%d)\n", ret.func->name, ret.func->argtype);776 break;777 case k_op1:778 case k_op1r:779 n = do_op(q,ret.op1r,"op1r",2,argtype);780 if (n < 1)781 goto balderdash;782 break;783 case k_op1i:784 n = do_op(q,ret.op1i,"op1i",2,argtype);785 if (n < 1)786 goto balderdash;787 break;788 case k_op2:789 case k_op2r:790 n = do_op(q,ret.op2r,"op2r",2,argtype);791 if (n < 1)792 goto balderdash;793 break;794 case k_op2i:795 n = do_op(q,ret.op2i,"op2i",2,argtype);796 if (n < 1)797 goto balderdash;798 break;799 case k_resultc:800 puts(q);801 if(strncmp(q,"inf",3)==0) {802 ret.resultc = rc_infinity;803 } else if(strcmp(q,"zero")==0) {804 ret.resultc = rc_zero;805 } else if(strcmp(q,"nan")==0) {806 ret.resultc = rc_nan;807 } else if(strcmp(q,"finite")==0) {808 ret.resultc = rc_finite;809 } else {810 goto balderdash;811 }812 break;813 case k_result:814 case k_resultr:815 n = (do_op)(q,ret.resultr,"resultr",3,rettype);816 if (n < 1)817 goto balderdash;818 ret.nresult = n; /* assume real and imaginary have same no. words */819 break;820 case k_resulti:821 n = do_op(q,ret.resulti,"resulti",3,rettype);822 if (n < 1)823 goto balderdash;824 break;825 case k_res2:826 n = do_op(q,ret.res2,"res2",2,rettype);827 if (n < 1)828 goto balderdash;829 break;830 case k_status:831 while (*q) {832 if (*q == 'i') ret.status |= FE_INVALID;833 if (*q == 'z') ret.status |= FE_DIVBYZERO;834 if (*q == 'o') ret.status |= FE_OVERFLOW;835 if (*q == 'u') ret.status |= FE_UNDERFLOW;836 q++;837 }838 break;839 case k_maybeerror:840 n = find(q, errors, sizeof(errors));841 if (n < 0)842 goto balderdash;843 if(math_errhandling&MATH_ERREXCEPT) {844 switch(n) {845 case e_domain: ret.maybestatus |= FE_INVALID; break;846 case e_divbyzero: ret.maybestatus |= FE_DIVBYZERO; break;847 case e_overflow: ret.maybestatus |= FE_OVERFLOW; break;848 case e_underflow: ret.maybestatus |= FE_UNDERFLOW; break;849 }850 }851 {852 switch(n) {853 case e_domain:854 ret.maybeerr = e_EDOM; break;855 case e_divbyzero:856 case e_overflow:857 case e_underflow:858 ret.maybeerr = e_ERANGE; break;859 }860 }861 case k_maybestatus:862 while (*q) {863 if (*q == 'i') ret.maybestatus |= FE_INVALID;864 if (*q == 'z') ret.maybestatus |= FE_DIVBYZERO;865 if (*q == 'o') ret.maybestatus |= FE_OVERFLOW;866 if (*q == 'u') ret.maybestatus |= FE_UNDERFLOW;867 q++;868 }869 break;870 case k_error:871 n = find(q, errors, sizeof(errors));872 if (n < 0)873 goto balderdash;874 if(math_errhandling&MATH_ERREXCEPT) {875 switch(n) {876 case e_domain: ret.status |= FE_INVALID; break;877 case e_divbyzero: ret.status |= FE_DIVBYZERO; break;878 case e_overflow: ret.status |= FE_OVERFLOW; break;879 case e_underflow: ret.status |= FE_UNDERFLOW; break;880 }881 }882 if(math_errhandling&MATH_ERRNO) {883 switch(n) {884 case e_domain:885 ret.err = e_EDOM; break;886 case e_divbyzero:887 case e_overflow:888 case e_underflow:889 ret.err = e_ERANGE; break;890 }891 }892 if(!(math_errhandling&MATH_ERRNO)) {893 switch(n) {894 case e_domain:895 ret.maybeerr = e_EDOM; break;896 case e_divbyzero:897 case e_overflow:898 case e_underflow:899 ret.maybeerr = e_ERANGE; break;900 }901 }902 break;903 case k_errno:904 ret.err = find(q, errnos, sizeof(errnos));905 if (ret.err < 0)906 goto balderdash;907 break;908 case k_errno_in:909 ret.in_err = find(q, errnos, sizeof(errnos));910 if (ret.err < 0)911 goto balderdash;912 ret.in_err_limit = ret.in_err + 1;913 break;914 case k_wrongresult:915 case k_wrongstatus:916 case k_wrongres2:917 case k_wrongerrno:918 /* quietly ignore these keys */919 break;920 default:921 goto balderdash;922 }923 p = strtok(NULL, " \t");924 }925 ret.valid = 1;926 return ret;927 928 /* come here from almost any error */929 balderdash:930 ret.valid = 0;931 return ret;932}933 934typedef enum {935 test_comment, /* deliberately not a test */936 test_invalid, /* accidentally not a test */937 test_decline, /* was a test, and wasn't run */938 test_fail, /* was a test, and failed */939 test_pass /* was a test, and passed */940} testresult;941 942char failtext[512];943 944typedef union {945 unsigned i[2];946 double f;947 double da[2];948} dbl;949 950typedef union {951 unsigned i;952 float f;953 float da[2];954} sgl;955 956/* helper function for runtest */957void print_error(int rettype, unsigned *result, char* text, char** failp) {958 special_op *sop;959 char *str;960 961 if(result) {962 *failp += sprintf(*failp," %s=",text);963 sop = find_special_op_from_op(result[0],result[1],is_double_rettype(rettype));964 if(sop) {965 *failp += sprintf(*failp,"%s",sop->name);966 } else {967 if(is_double_rettype(rettype)) {968 str="%08x.%08x";969 } else {970 str="%08x";971 }972 *failp += sprintf(*failp,str,result[0],result[1]);973 }974 }975}976 977 978void print_ulps_helper(const char *name, long long ulps, char** failp) {979 if(ulps == LLONG_MAX) {980 *failp += sprintf(*failp, " %s=HUGE", name);981 } else {982 *failp += sprintf(*failp, " %s=%.3f", name, (double)ulps / ULPUNIT);983 }984}985 986/* for complex args make ulpsr or ulpsri = 0 to not print */987void print_ulps(int rettype, long long ulpsr, long long ulpsi, char** failp) {988 if(is_complex_rettype(rettype)) {989 if (ulpsr) print_ulps_helper("ulpsr",ulpsr,failp);990 if (ulpsi) print_ulps_helper("ulpsi",ulpsi,failp);991 } else {992 if (ulpsr) print_ulps_helper("ulps",ulpsr,failp);993 }994}995 996int runtest(testdetail t) {997 int err, status;998 999 dbl d_arg1, d_arg2, d_res, d_res2;1000 sgl s_arg1, s_arg2, s_res, s_res2;1001 1002 int deferred_decline = FALSE;1003 char *failp = failtext;1004 1005 unsigned int intres=0;1006 1007 int res2_adjust = 0;1008 1009 if (t.comment)1010 return test_comment;1011 if (!t.valid)1012 return test_invalid;1013 1014 /* Set IEEE status to mathlib-normal */1015 feclearexcept(FE_ALL_EXCEPT);1016 1017 /* Deal with operands */1018#define DO_DOP(arg,op) arg.i[dmsd] = t.op[0]; arg.i[dlsd] = t.op[1]1019 DO_DOP(d_arg1,op1r);1020 DO_DOP(d_arg2,op2r);1021 s_arg1.i = t.op1r[0]; s_arg2.i = t.op2r[0];1022 1023 /*1024 * Detect NaNs, infinities and denormals on input, and set a1025 * deferred decline flag if we're in FO mode.1026 *1027 * (We defer the decline rather than doing it immediately1028 * because even in FO mode the operation is not permitted to1029 * crash or tight-loop; so we _run_ the test, and then ignore1030 * all the results.)1031 */1032 if (fo) {1033 if (is_double_argtype(t.func->argtype) && is_dhard(t.op1r))1034 deferred_decline = TRUE;1035 if (t.func->argtype==at_d2 && is_dhard(t.op2r))1036 deferred_decline = TRUE;1037 if (is_single_argtype(t.func->argtype) && is_shard(t.op1r))1038 deferred_decline = TRUE;1039 if (t.func->argtype==at_s2 && is_shard(t.op2r))1040 deferred_decline = TRUE;1041 if (is_double_rettype(t.func->rettype) && is_dhard(t.resultr))1042 deferred_decline = TRUE;1043 if (t.func->rettype==rt_d2 && is_dhard(t.res2))1044 deferred_decline = TRUE;1045 if (is_single_argtype(t.func->rettype) && is_shard(t.resultr))1046 deferred_decline = TRUE;1047 if (t.func->rettype==rt_s2 && is_shard(t.res2))1048 deferred_decline = TRUE;1049 if (t.err == e_ERANGE)1050 deferred_decline = TRUE;1051 }1052 1053 /*1054 * Perform the operation1055 */1056 1057 errno = t.in_err == e_EDOM ? EDOM : t.in_err == e_ERANGE ? ERANGE : 0;1058 if (t.err == e_0)1059 t.err = t.in_err;1060 if (t.maybeerr == e_0)1061 t.maybeerr = t.in_err;1062 1063 if(t.func->type == t_func) {1064 switch(t.func->argtype) {1065 case at_d: d_res.f = t.func->func.d_d_ptr(d_arg1.f); break;1066 case at_s: s_res.f = t.func->func.s_s_ptr(s_arg1.f); break;1067 case at_d2: d_res.f = t.func->func.d2_d_ptr(d_arg1.f, d_arg2.f); break;1068 case at_s2: s_res.f = t.func->func.s2_s_ptr(s_arg1.f, s_arg2.f); break;1069 case at_di: d_res.f = t.func->func.di_d_ptr(d_arg1.f, d_arg2.i[dmsd]); break;1070 case at_si: s_res.f = t.func->func.si_s_ptr(s_arg1.f, s_arg2.i); break;1071 case at_dip: d_res.f = t.func->func.dip_d_ptr(d_arg1.f, (int*)&intres); break;1072 case at_sip: s_res.f = t.func->func.sip_s_ptr(s_arg1.f, (int*)&intres); break;1073 case at_ddp: d_res.f = t.func->func.ddp_d_ptr(d_arg1.f, &d_res2.f); break;1074 case at_ssp: s_res.f = t.func->func.ssp_s_ptr(s_arg1.f, &s_res2.f); break;1075 default:1076 printf("unhandled function: %s\n",t.func->name);1077 return test_fail;1078 }1079 } else {1080 /* printf("macro: name=%s, num=%i, s1.i=0x%08x s1.f=%f\n",t.func->name, t.func->macro_name, s_arg1.i, (double)s_arg1.f); */1081 switch(t.func->macro_name) {1082 case m_isfinite: intres = isfinite(d_arg1.f); break;1083 case m_isinf: intres = isinf(d_arg1.f); break;1084 case m_isnan: intres = isnan(d_arg1.f); break;1085 case m_isnormal: intres = isnormal(d_arg1.f); break;1086 case m_signbit: intres = signbit(d_arg1.f); break;1087 case m_fpclassify: intres = fpclassify(d_arg1.f); break;1088 case m_isgreater: intres = isgreater(d_arg1.f, d_arg2.f); break;1089 case m_isgreaterequal: intres = isgreaterequal(d_arg1.f, d_arg2.f); break;1090 case m_isless: intres = isless(d_arg1.f, d_arg2.f); break;1091 case m_islessequal: intres = islessequal(d_arg1.f, d_arg2.f); break;1092 case m_islessgreater: intres = islessgreater(d_arg1.f, d_arg2.f); break;1093 case m_isunordered: intres = isunordered(d_arg1.f, d_arg2.f); break;1094 1095 case m_isfinitef: intres = isfinite(s_arg1.f); break;1096 case m_isinff: intres = isinf(s_arg1.f); break;1097 case m_isnanf: intres = isnan(s_arg1.f); break;1098 case m_isnormalf: intres = isnormal(s_arg1.f); break;1099 case m_signbitf: intres = signbit(s_arg1.f); break;1100 case m_fpclassifyf: intres = fpclassify(s_arg1.f); break;1101 case m_isgreaterf: intres = isgreater(s_arg1.f, s_arg2.f); break;1102 case m_isgreaterequalf: intres = isgreaterequal(s_arg1.f, s_arg2.f); break;1103 case m_islessf: intres = isless(s_arg1.f, s_arg2.f); break;1104 case m_islessequalf: intres = islessequal(s_arg1.f, s_arg2.f); break;1105 case m_islessgreaterf: intres = islessgreater(s_arg1.f, s_arg2.f); break;1106 case m_isunorderedf: intres = isunordered(s_arg1.f, s_arg2.f); break;1107 1108 default:1109 printf("unhandled macro: %s\n",t.func->name);1110 return test_fail;1111 }1112 }1113 1114 /*1115 * Decline the test if the deferred decline flag was set above.1116 */1117 if (deferred_decline)1118 return test_decline;1119 1120 /* printf("intres=%i\n",intres); */1121 1122 /* Clear the fail text (indicating a pass unless we change it) */1123 failp[0] = '\0';1124 1125 /* Check the IEEE status bits (except INX, which we disregard).1126 * We don't bother with this for complex numbers, because the1127 * complex functions are hard to get exactly right and we don't1128 * have to anyway (C99 annex G is only informative). */1129 if (!(is_complex_argtype(t.func->argtype) || is_complex_rettype(t.func->rettype))) {1130 status = fetestexcept(FE_INVALID|FE_DIVBYZERO|FE_OVERFLOW|FE_UNDERFLOW);1131 if ((status|t.maybestatus|~statusmask) != (t.status|t.maybestatus|~statusmask)) {1132 if (quiet) failtext[0]='x';1133 else {1134 failp += sprintf(failp,1135 " wrongstatus=%s%s%s%s%s",1136 (status & FE_INVALID ? "i" : ""),1137 (status & FE_DIVBYZERO ? "z" : ""),1138 (status & FE_OVERFLOW ? "o" : ""),1139 (status & FE_UNDERFLOW ? "u" : ""),1140 (status ? "" : "OK"));1141 }1142 }1143 }1144 1145 /* Check the result */1146 {1147 unsigned resultr[2], resulti[2];1148 unsigned tresultr[3], tresulti[3], wres;1149 1150 switch(t.func->rettype) {1151 case rt_d:1152 case rt_d2:1153 tresultr[0] = t.resultr[0];1154 tresultr[1] = t.resultr[1];1155 resultr[0] = d_res.i[dmsd]; resultr[1] = d_res.i[dlsd];1156 wres = 2;1157 break;1158 case rt_i:1159 tresultr[0] = t.resultr[0];1160 resultr[0] = intres;1161 wres = 1;1162 break;1163 case rt_s:1164 case rt_s2:1165 tresultr[0] = t.resultr[0];1166 resultr[0] = s_res.i;1167 wres = 1;1168 break;1169 default:1170 puts("unhandled rettype in runtest");1171 wres = 0;1172 }1173 if(t.resultc != rc_none) {1174 int err = 0;1175 switch(t.resultc) {1176 case rc_zero:1177 if(resultr[0] != 0 || resulti[0] != 0 ||1178 (wres==2 && (resultr[1] != 0 || resulti[1] != 0))) {1179 err = 1;1180 }1181 break;1182 case rc_infinity:1183 if(wres==1) {1184 if(!((resultr[0]&0x7fffffff)==0x7f800000 ||1185 (resulti[0]&0x7fffffff)==0x7f800000)) {1186 err = 1;1187 }1188 } else {1189 if(!(((resultr[0]&0x7fffffff)==0x7ff00000 && resultr[1]==0) ||1190 ((resulti[0]&0x7fffffff)==0x7ff00000 && resulti[1]==0))) {1191 err = 1;1192 }1193 }1194 break;1195 case rc_nan:1196 if(wres==1) {1197 if(!((resultr[0]&0x7fffffff)>0x7f800000 ||1198 (resulti[0]&0x7fffffff)>0x7f800000)) {1199 err = 1;1200 }1201 } else {1202 canon_dNaN(resultr);1203 canon_dNaN(resulti);1204 if(!(((resultr[0]&0x7fffffff)>0x7ff00000 && resultr[1]==1) ||1205 ((resulti[0]&0x7fffffff)>0x7ff00000 && resulti[1]==1))) {1206 err = 1;1207 }1208 }1209 break;1210 case rc_finite:1211 if(wres==1) {1212 if(!((resultr[0]&0x7fffffff)<0x7f800000 ||1213 (resulti[0]&0x7fffffff)<0x7f800000)) {1214 err = 1;1215 }1216 } else {1217 if(!((resultr[0]&0x7fffffff)<0x7ff00000 ||1218 (resulti[0]&0x7fffffff)<0x7ff00000)) {1219 err = 1;1220 }1221 }1222 break;1223 default:1224 break;1225 }1226 if(err) {1227 print_error(t.func->rettype,resultr,"wrongresultr",&failp);1228 print_error(t.func->rettype,resulti,"wrongresulti",&failp);1229 }1230 } else if (t.nresult > wres) {1231 /*1232 * The test case data has provided the result to more1233 * than double precision. Instead of testing exact1234 * equality, we test against our maximum error1235 * tolerance.1236 */1237 int rshift, ishift;1238 long long ulpsr, ulpsi, ulptolerance;1239 1240 tresultr[wres] = t.resultr[wres] << (32-EXTRABITS);1241 tresulti[wres] = t.resulti[wres] << (32-EXTRABITS);1242 if(strict) {1243 ulptolerance = 4096; /* one ulp */1244 } else {1245 ulptolerance = t.func->tolerance;1246 }1247 rshift = ishift = 0;1248 if (ulptolerance & ABSLOWERBOUND) {1249 /*1250 * Hack for the lgamma functions, which have an1251 * error behaviour that can't conveniently be1252 * characterised in pure ULPs. Really, we want to1253 * say that the error in lgamma is "at most N ULPs,1254 * or at most an absolute error of X, whichever is1255 * larger", for appropriately chosen N,X. But since1256 * these two functions are the only cases where it1257 * arises, I haven't bothered to do it in a nice way1258 * in the function table above.1259 *1260 * (The difficult cases arise with negative input1261 * values such that |gamma(x)| is very near to 1; in1262 * this situation implementations tend to separately1263 * compute lgamma(|x|) and the log of the correction1264 * term from the Euler reflection formula, and1265 * subtract - which catastrophically loses1266 * significance.)1267 *1268 * As far as I can tell, nobody cares about this:1269 * GNU libm doesn't get those cases right either,1270 * and OpenCL explicitly doesn't state a ULP error1271 * limit for lgamma. So my guess is that this is1272 * simply considered acceptable error behaviour for1273 * this particular function, and hence I feel free1274 * to allow for it here.1275 */1276 ulptolerance &= ~ABSLOWERBOUND;1277 if (t.op1r[0] & 0x80000000) {1278 if (t.func->rettype == rt_d)1279 rshift = 0x400 - ((tresultr[0] >> 20) & 0x7ff);1280 else if (t.func->rettype == rt_s)1281 rshift = 0x80 - ((tresultr[0] >> 23) & 0xff);1282 if (rshift < 0)1283 rshift = 0;1284 }1285 }1286 if (ulptolerance & PLUSMINUSPIO2) {1287 ulptolerance &= ~PLUSMINUSPIO2;1288 /*1289 * Hack for range reduction, which can reduce1290 * borderline cases in the wrong direction, i.e.1291 * return a value just outside one end of the interval1292 * [-pi/4,+pi/4] when it could have returned a value1293 * just inside the other end by subtracting an1294 * adjacent multiple of pi/2.1295 *1296 * We tolerate this, up to a point, because the1297 * trigonometric functions making use of the output of1298 * rred can cope and because making the range reducer1299 * do the exactly right thing in every case would be1300 * more expensive.1301 */1302 if (wres == 1) {1303 /* Upper bound of overshoot derived in rredf.h */1304 if ((resultr[0]&0x7FFFFFFF) <= 0x3f494b02 &&1305 (resultr[0]&0x7FFFFFFF) > 0x3f490fda &&1306 (resultr[0]&0x80000000) != (tresultr[0]&0x80000000)) {1307 unsigned long long val;1308 val = tresultr[0];1309 val = (val << 32) | tresultr[1];1310 /*1311 * Compute the alternative permitted result by1312 * subtracting from the sum of the extended1313 * single-precision bit patterns of +pi/4 and1314 * -pi/4. This is a horrible hack which only1315 * works because we can be confident that1316 * numbers in this range all have the same1317 * exponent!1318 */1319 val = 0xfe921fb54442d184ULL - val;1320 tresultr[0] = val >> 32;1321 tresultr[1] = (val >> (32-EXTRABITS)) << (32-EXTRABITS);1322 /*1323 * Also, expect a correspondingly different1324 * value of res2 as a result of this change.1325 * The adjustment depends on whether we just1326 * flipped the result from + to - or vice1327 * versa.1328 */1329 if (resultr[0] & 0x80000000) {1330 res2_adjust = +1;1331 } else {1332 res2_adjust = -1;1333 }1334 }1335 }1336 }1337 ulpsr = calc_error(resultr, tresultr, rshift, t.func->rettype);1338 if(is_complex_rettype(t.func->rettype)) {1339 ulpsi = calc_error(resulti, tresulti, ishift, t.func->rettype);1340 } else {1341 ulpsi = 0;1342 }1343 unsigned *rr = (ulpsr > ulptolerance || ulpsr < -ulptolerance) ? resultr : NULL;1344 unsigned *ri = (ulpsi > ulptolerance || ulpsi < -ulptolerance) ? resulti : NULL;1345/* printf("tolerance=%i, ulpsr=%i, ulpsi=%i, rr=%p, ri=%p\n",ulptolerance,ulpsr,ulpsi,rr,ri); */1346 if (rr || ri) {1347 if (quiet) failtext[0]='x';1348 else {1349 print_error(t.func->rettype,rr,"wrongresultr",&failp);1350 print_error(t.func->rettype,ri,"wrongresulti",&failp);1351 print_ulps(t.func->rettype,rr ? ulpsr : 0, ri ? ulpsi : 0,&failp);1352 }1353 }1354 } else {1355 if(is_complex_rettype(t.func->rettype))1356 /*1357 * Complex functions are not fully supported,1358 * this is unreachable, but prevents warnings.1359 */1360 abort();1361 /*1362 * The test case data has provided the result in1363 * exactly the output precision. Therefore we must1364 * complain about _any_ violation.1365 */1366 switch(t.func->rettype) {1367 case rt_dc:1368 canon_dNaN(tresulti);1369 canon_dNaN(resulti);1370 if (fo) {1371 dnormzero(tresulti);1372 dnormzero(resulti);1373 }1374 /* deliberate fall-through */1375 case rt_d:1376 canon_dNaN(tresultr);1377 canon_dNaN(resultr);1378 if (fo) {1379 dnormzero(tresultr);1380 dnormzero(resultr);1381 }1382 break;1383 case rt_sc:1384 canon_sNaN(tresulti);1385 canon_sNaN(resulti);1386 if (fo) {1387 snormzero(tresulti);1388 snormzero(resulti);1389 }1390 /* deliberate fall-through */1391 case rt_s:1392 canon_sNaN(tresultr);1393 canon_sNaN(resultr);1394 if (fo) {1395 snormzero(tresultr);1396 snormzero(resultr);1397 }1398 break;1399 default:1400 break;1401 }1402 if(is_complex_rettype(t.func->rettype)) {1403 unsigned *rr, *ri;1404 if(resultr[0] != tresultr[0] ||1405 (wres > 1 && resultr[1] != tresultr[1])) {1406 rr = resultr;1407 } else {1408 rr = NULL;1409 }1410 if(resulti[0] != tresulti[0] ||1411 (wres > 1 && resulti[1] != tresulti[1])) {1412 ri = resulti;1413 } else {1414 ri = NULL;1415 }1416 if(rr || ri) {1417 if (quiet) failtext[0]='x';1418 print_error(t.func->rettype,rr,"wrongresultr",&failp);1419 print_error(t.func->rettype,ri,"wrongresulti",&failp);1420 }1421 } else if (resultr[0] != tresultr[0] ||1422 (wres > 1 && resultr[1] != tresultr[1])) {1423 if (quiet) failtext[0]='x';1424 print_error(t.func->rettype,resultr,"wrongresult",&failp);1425 }1426 }1427 /*1428 * Now test res2, for those functions (frexp, modf, rred)1429 * which use it.1430 */1431 if (t.func->func.ptr == &frexp || t.func->func.ptr == &frexpf ||1432 t.func->macro_name == m_rred || t.func->macro_name == m_rredf) {1433 unsigned tres2 = t.res2[0];1434 if (res2_adjust) {1435 /* Fix for range reduction, propagated from further up */1436 tres2 = (tres2 + res2_adjust) & 3;1437 }1438 if (tres2 != intres) {1439 if (quiet) failtext[0]='x';1440 else {1441 failp += sprintf(failp,1442 " wrongres2=%08x", intres);1443 }1444 }1445 } else if (t.func->func.ptr == &modf || t.func->func.ptr == &modff) {1446 tresultr[0] = t.res2[0];1447 tresultr[1] = t.res2[1];1448 if (is_double_rettype(t.func->rettype)) {1449 canon_dNaN(tresultr);1450 resultr[0] = d_res2.i[dmsd];1451 resultr[1] = d_res2.i[dlsd];1452 canon_dNaN(resultr);1453 if (fo) {1454 dnormzero(tresultr);1455 dnormzero(resultr);1456 }1457 } else {1458 canon_sNaN(tresultr);1459 resultr[0] = s_res2.i;1460 resultr[1] = s_res2.i;1461 canon_sNaN(resultr);1462 if (fo) {1463 snormzero(tresultr);1464 snormzero(resultr);1465 }1466 }1467 if (resultr[0] != tresultr[0] ||1468 (wres > 1 && resultr[1] != tresultr[1])) {1469 if (quiet) failtext[0]='x';1470 else {1471 if (is_double_rettype(t.func->rettype))1472 failp += sprintf(failp, " wrongres2=%08x.%08x",1473 resultr[0], resultr[1]);1474 else1475 failp += sprintf(failp, " wrongres2=%08x",1476 resultr[0]);1477 }1478 }1479 }1480 }1481 1482 /* Check errno */1483 err = (errno == EDOM ? e_EDOM : errno == ERANGE ? e_ERANGE : e_0);1484 if (err != t.err && err != t.maybeerr) {1485 if (quiet) failtext[0]='x';1486 else {1487 failp += sprintf(failp, " wrongerrno=%s expecterrno=%s ", errnos[err], errnos[t.err]);1488 }1489 }1490 1491 return *failtext ? test_fail : test_pass;1492}1493 1494int passed, failed, declined;1495 1496void runtests(char *name, FILE *fp) {1497 char testbuf[512], linebuf[512];1498 int lineno = 1;1499 testdetail test;1500 1501 test.valid = 0;1502 1503 if (verbose) printf("runtests: %s\n", name);1504 while (fgets(testbuf, sizeof(testbuf), fp)) {1505 int res, print_errno;1506 testbuf[strcspn(testbuf, "\r\n")] = '\0';1507 strcpy(linebuf, testbuf);1508 test = parsetest(testbuf, test);1509 print_errno = 0;1510 while (test.in_err < test.in_err_limit) {1511 res = runtest(test);1512 if (res == test_pass) {1513 if (verbose)1514 printf("%s:%d: pass\n", name, lineno);1515 ++passed;1516 } else if (res == test_decline) {1517 if (verbose)1518 printf("%s:%d: declined\n", name, lineno);1519 ++declined;1520 } else if (res == test_fail) {1521 if (!quiet)1522 printf("%s:%d: FAIL%s: %s%s%s%s\n", name, lineno,1523 test.random ? " (random)" : "",1524 linebuf,1525 print_errno ? " errno_in=" : "",1526 print_errno ? errnos[test.in_err] : "",1527 failtext);1528 ++failed;1529 } else if (res == test_invalid) {1530 printf("%s:%d: malformed: %s\n", name, lineno, linebuf);1531 ++failed;1532 }1533 test.in_err++;1534 print_errno = 1;1535 }1536 lineno++;1537 }1538}1539 1540int main(int ac, char **av) {1541 char **files;1542 int i, nfiles = 0;1543 dbl d;1544 1545#ifdef MICROLIB1546 /*1547 * Invent argc and argv ourselves.1548 */1549 char *argv[256];1550 char args[256];1551 {1552 int sargs[2];1553 char *p;1554 1555 ac = 0;1556 1557 sargs[0]=(int)args;1558 sargs[1]=(int)sizeof(args);1559 if (!__semihost(0x15, sargs)) {1560 args[sizeof(args)-1] = '\0'; /* just in case */1561 p = args;1562 while (1) {1563 while (*p == ' ' || *p == '\t') p++;1564 if (!*p) break;1565 argv[ac++] = p;1566 while (*p && *p != ' ' && *p != '\t') p++;1567 if (*p) *p++ = '\0';1568 }1569 }1570 1571 av = argv;1572 }1573#endif1574 1575 /* Sort tfuncs */1576 qsort(tfuncs, sizeof(tfuncs)/sizeof(test_func), sizeof(test_func), &compare_tfuncs);1577 1578 /*1579 * Autodetect the `double' endianness.1580 */1581 dmsd = 0;1582 d.f = 1.0; /* 0x3ff00000 / 0x00000000 */1583 if (d.i[dmsd] == 0) {1584 dmsd = 1;1585 }1586 /*1587 * Now dmsd denotes what the compiler thinks we're at. Let's1588 * check that it agrees with what the runtime thinks.1589 */1590 d.i[0] = d.i[1] = 0x11111111;/* a random +ve number */1591 d.f /= d.f; /* must now be one */1592 if (d.i[dmsd] == 0) {1593 fprintf(stderr, "YIKES! Compiler and runtime disagree on endianness"1594 " of `double'. Bailing out\n");1595 return 1;1596 }1597 dlsd = !dmsd;1598 1599 /* default is terse */1600 verbose = 0;1601 fo = 0;1602 strict = 0;1603 1604 files = (char **)malloc((ac+1) * sizeof(char *));1605 if (!files) {1606 fprintf(stderr, "initial malloc failed!\n");1607 return 1;1608 }1609#ifdef NOCMDLINE1610 files[nfiles++] = "testfile";1611#endif1612 1613 while (--ac) {1614 char *p = *++av;1615 if (*p == '-') {1616 static char *options[] = {1617 "-fo",1618#if 01619 "-noinexact",1620 "-noround",1621#endif1622 "-nostatus",1623 "-quiet",1624 "-strict",1625 "-v",1626 "-verbose",1627 };1628 enum {1629 op_fo,1630#if 01631 op_noinexact,1632 op_noround,1633#endif1634 op_nostatus,1635 op_quiet,1636 op_strict,1637 op_v,1638 op_verbose,1639 };1640 switch (find(p, options, sizeof(options))) {1641 case op_quiet:1642 quiet = 1;1643 break;1644#if 01645 case op_noinexact:1646 statusmask &= 0x0F; /* remove bit 4 */1647 break;1648 case op_noround:1649 doround = 0;1650 break;1651#endif1652 case op_nostatus: /* no status word => noinx,noround */1653 statusmask = 0;1654 doround = 0;1655 break;1656 case op_v:1657 case op_verbose:1658 verbose = 1;1659 break;1660 case op_fo:1661 fo = 1;1662 break;1663 case op_strict: /* tolerance is 1 ulp */1664 strict = 1;1665 break;1666 default:1667 fprintf(stderr, "unrecognised option: %s\n", p);1668 break;1669 }1670 } else {1671 files[nfiles++] = p;1672 }1673 }1674 1675 passed = failed = declined = 0;1676 1677 if (nfiles) {1678 for (i = 0; i < nfiles; i++) {1679 FILE *fp = fopen(files[i], "r");1680 if (!fp) {1681 fprintf(stderr, "Couldn't open %s\n", files[i]);1682 } else1683 runtests(files[i], fp);1684 }1685 } else1686 runtests("(stdin)", stdin);1687 1688 printf("Completed. Passed %d, failed %d (total %d",1689 passed, failed, passed+failed);1690 if (declined)1691 printf(" plus %d declined", declined);1692 printf(")\n");1693 if (failed || passed == 0)1694 return 1;1695 printf("** TEST PASSED OK **\n");1696 return 0;1697}1698 1699void undef_func() {1700 failed++;1701 puts("ERROR: undefined function called");1702}1703