summaryrefslogtreecommitdiff
path: root/test/math/ieetst.c
diff options
context:
space:
mode:
Diffstat (limited to 'test/math/ieetst.c')
-rw-r--r--test/math/ieetst.c1700
1 files changed, 850 insertions, 850 deletions
diff --git a/test/math/ieetst.c b/test/math/ieetst.c
index 3c89145e3..8c2453898 100644
--- a/test/math/ieetst.c
+++ b/test/math/ieetst.c
@@ -1,850 +1,850 @@
-/* Floating point to ASCII input and output string test program.
- *
- * Numbers in the native machine data structure are converted
- * to e type, then to and from decimal ASCII strings. Native
- * printf() and scanf() functions are also used to produce
- * and read strings. The resulting e type binary values
- * are compared, with diagnostic printouts of any discrepancies.
- *
- * Steve Moshier, 16 Dec 88
- * last revision: 16 May 92
- */
-
-#include "ehead.h"
-#include "mconf.h"
-
-/* Include tests of 80-bit long double precision: */
-#define LDOUBLE 0
-/* Abort subtest after getting this many errors: */
-#define MAXERR 5
-/* Number of random arguments to try (set as large as you have
- * patience for): */
-#define NRAND 100
-/* Perform internal consistency test: */
-#define CHKINTERNAL 0
-
-static unsigned short fullp[NE], rounded[NE];
-float prec24, sprec24, ssprec24;
-double prec53, sprec53, ssprec53;
-#if LDOUBLE
-long double prec64, sprec64, ssprec64;
-#endif
-
-static unsigned short rprint[NE], rscan[NE];
-static unsigned short q1[NE], q2[NE], q5[NE];
-static unsigned short e1[NE], e2[NE], e3[NE];
-static double d1, d2;
-static int errprint = 0;
-static int errscan = 0;
-static int identerr = 0;
-static int errtot = 0;
-static int count = 0;
-static char str0[80], str1[80], str2[80], str3[80];
-static unsigned short eten[NE], maxm[NE];
-
-int m, n, k2, mprec, SPREC;
-
-char *Ten = "10.0";
-char tformat[10];
-char *format24 = "%.8e";
-#ifdef DEC
-char *format53 = "%.17e";
-#else
-char *format53 = "%.16e";
-#endif
-char *fformat24 = "%e";
-char *fformat53 = "%le";
-char *pct = "%";
-char *quo = "\042";
-#if LDOUBLE
-char *format64 = "%.20Le";
-char *fformat64 = "%Le";
-#endif
-char *format;
-char *fformat;
-char *toomany = "Too many errors; aborting this test.\n";
-
-static int mnrflag;
-static int etrflag;
-void chkit(), printerr(), mnrand(), etrand(), shownoncrit();
-void chkid(), pvec();
-
-main()
-{
-int i, iprec;
-
-printf( "Steve Moshier's printf/scanf tester, version 0.2.\n\n" );
-#ifdef DEC
- /* DEC PDP-11/VAX single precision not yet implemented */
-for( iprec = 1; iprec<2; iprec++ )
-#else
-for( iprec = 0; iprec<3; iprec++ )
-#endif
- {
- errscan = 0;
- identerr = 0;
- errprint = 0;
- eclear( rprint );
- eclear( rscan );
-
-switch( iprec )
- {
- case 0:
- SPREC = 8; /* # digits after the decimal point */
- mprec = 24; /* # bits in the significand */
- m = 9; /* max # decimal digits for correct rounding */
- n = 13; /* max power of ten for correct rounding */
- k2 = -125; /* underflow beyond 2^-k2 */
- format = format24; /* printf format string */
- fformat = fformat24; /* scanf format string */
- mnrflag = 1; /* sets interval for random numbers */
- etrflag = 1;
- printf( "Testing FLOAT precision.\n" );
- break;
-
- case 1:
-#ifdef DEC
- SPREC = 17;
- mprec = 56;
- m = 17;
- n = 27;
- k2 = -125;
- format = format53;
- fformat = fformat53;
- mnrflag = 2;
- etrflag = 1;
- printf( "Testing DEC DOUBLE precision.\n" );
- break;
-#else
- SPREC = 16;
- mprec = 53;
- m = 17;
- n = 27;
- k2 = -1021;
- format = format53;
- fformat = fformat53;
- mnrflag = 2;
- etrflag = 2;
- printf( "Testing DOUBLE precision.\n" );
- break;
-#endif
- case 2:
-#if LDOUBLE
- SPREC = 20;
- mprec = 64;
- m = 20;
- n = 34;
- k2 = -16382;
- format = format64;
- fformat = fformat64;
- mnrflag = 3;
- etrflag = 3;
- printf( "Testing LONG DOUBLE precision.\n" );
- break;
-#else
- goto nodenorm;
-#endif
- }
-
- asctoe( Ten, eten );
-/* 10^m - 1 */
- d2 = m;
- e53toe( &d2, e1 );
- epow( eten, e1, maxm );
- esub( eone, maxm, maxm );
-
-/* test 1 */
- printf( "1. Checking 10^n - 1 for n = %d to %d.\n", -m, m );
- emov( eone, q5 );
- for( count=0; count<=m; count++ )
- {
- esub( eone, q5, fullp );
- chkit( 1 );
- ediv( q5, eone, q2 );
- esub( eone, q2, fullp );
- chkit( 1 );
- emul( eten, q5, q5 );
- if( errtot >= MAXERR )
- {
- printf( "%s", toomany );
- goto end1;
- }
- }
-end1:
- printerr();
-
-
-/* test 2 */
- printf( "2. Checking powers of 10 from 10^-%d to 10^%d.\n", n, n );
- emov( eone, q5 );
- for( count=0; count<=n; count++ )
- {
- emov( q5, fullp );
- chkit( 2 );
- ediv( q5, eone, fullp );
- chkit( 2 );
- emul( eten, q5, q5 );
- if( errtot >= MAXERR )
- {
- printf( "%s", toomany );
- goto end2;
- }
- }
-end2:
- printerr();
-
-/* test 3 */
- printf( "3. Checking (10^%d-1)*10^n from n = -%d to %d.\n", m, n, n );
- emov( eone, q5 );
- for( count= -n; count<=n; count++ )
- {
- emul( maxm, q5, fullp );
- chkit( 3 );
- emov( q5, fullp );
- ediv( fullp, eone, fullp );
- emul( maxm, fullp, fullp );
- chkit( 3 );
- emul( eten, q5, q5 );
- if( errtot >= MAXERR )
- {
- printf( "%s", toomany );
- goto end3;
- }
- }
-end3:
- printerr();
-
-
-
-/* test 4 */
- printf( "4. Checking powers of 2 from 2^-24 to 2^+56.\n" );
- d1 = -24.0;
- e53toe( &d1, q1 );
- epow( etwo, q1, q5 );
-
- for( count = -24; count <= 56; count++ )
- {
- emov( q5, fullp );
- chkit( 4 );
- emul( etwo, q5, q5 );
- if( errtot >= MAXERR )
- {
- printf( "%s", toomany );
- goto end4;
- }
- }
-end4:
- printerr();
-
-
-/* test 5 */
- printf( "5. Checking 2^n - 1 for n = 0 to %d.\n", mprec );
- emov( eone, q5 );
- for( count=0; count<=mprec; count++ )
- {
- esub( eone, q5, fullp );
- chkit( 5 );
- emul( etwo, q5, q5 );
- if( errtot >= MAXERR )
- {
- printf( "%s", toomany );
- goto end5;
- }
- }
-end5:
- printerr();
-
-/* test 6 */
- printf( "6. Checking 2^n + 1 for n = 0 to %d.\n", mprec );
- emov( eone, q5 );
- for( count=0; count<=mprec; count++ )
- {
- eadd( eone, q5, fullp );
- chkit( 6 );
- emul( etwo, q5, q5 );
- if( errtot >= MAXERR )
- {
- printf( "%s", toomany );
- goto end6;
- }
- }
-end6:
- printerr();
-
-/* test 7 */
- printf(
- "7. Checking %d values M * 10^N with random integer M and N,\n",
- NRAND );
- printf(" 1 <= M <= 10^%d - 1 and -%d <= N <= +%d.\n", m, n, n );
- for( i=0; i<NRAND; i++ )
- {
- mnrand( fullp );
- chkit( 7 );
- if( errtot >= MAXERR )
- {
- printf( "%s", toomany );
- goto end7;
- }
- }
-end7:
- printerr();
-
-/* test 8 */
- printf("8. Checking critical rounding cases.\n" );
- for( i=0; i<20; i++ )
- {
- mnrand( fullp );
- eabs( fullp );
- if( ecmp( fullp, eone ) < 0 )
- ediv( fullp, eone, fullp );
- efloor( fullp, fullp );
- eadd( ehalf, fullp, fullp );
- chkit( 8 );
- if( errtot >= MAXERR )
- {
- printf( "%s", toomany );
- goto end8;
- }
- }
-end8:
- printerr();
-
-
-
-/* test 9 */
- printf("9. Testing on %d random non-denormal values.\n", NRAND );
- for( i=0; i<NRAND; i++ )
- {
- etrand( fullp );
- chkit( 9 );
- }
- printerr();
- shownoncrit();
-
-/* test 10 */
- printf(
- "Do you want to check denormal numbers in this precision ? (y/n) " );
- gets( str0 );
- if( str0[0] != 'y' )
- goto nodenorm;
-
- printf( "10. Checking denormal numbers.\n" );
-
-/* Form 2^-starting power */
- d1 = k2;
- e53toe( &d1, q1 );
- epow( etwo, q1, e1 );
-
-/* Find 2^-mprec less than starting power */
- d1 = -mprec + 4;
- e53toe( &d1, q1 );
- epow( etwo, q1, e3 );
- emul( e1, e3, e3 );
- emov( e3, e2 );
- ediv( etwo, e2, e2 );
-
- while( ecmp(e1,e2) != 0 )
- {
- eadd( e1, e2, fullp );
- switch( mprec )
- {
-#if LDOUBLE
- case 64:
- etoe64( e1, &sprec64 );
- e64toe( &sprec64, q1 );
- etoe64( fullp, &prec64 );
- e64toe( &prec64, q2 );
- break;
-#endif
-#ifdef DEC
- case 56:
-#endif
- case 53:
- etoe53( e1, &sprec53 );
- e53toe( &sprec53, q1 );
- etoe53( fullp, &prec53 );
- e53toe( &prec53, q2 );
- break;
-
- case 24:
- etoe24( e1, &sprec24 );
- e24toe( &sprec24, q1 );
- etoe24( fullp, &prec24 );
- e24toe( &prec24, q2 );
- break;
- }
- if( ecmp( q2, ezero ) == 0 )
- goto maxden;
- chkit(10);
- if( ecmp(q1,q2) == 0 )
- {
- ediv( etwo, e1, e1 );
- emov( e3, e2 );
- }
- if( errtot >= MAXERR )
- {
- printf( "%s", toomany );
- goto maxden;
- }
- ediv( etwo, e2, e2 );
- }
-maxden:
- printerr();
-nodenorm:
- printf( "\n" );
- } /* loop on precision */
-printf( "End of test.\n" );
-}
-
-#if CHKINTERNAL
-long double xprec64;
-double xprec53;
-float xprec24;
-
-/* Check binary -> printf -> scanf -> binary identity
- * of internal routines
- */
-void chkinternal( ref, tst, string )
-unsigned short ref[], tst[];
-char *string;
-{
-
-if( ecmp(ref,tst) != 0 )
- {
- printf( "internal identity compare error!\n" );
- chkid( ref, tst, string );
- }
-}
-#endif
-
-
-/* Check binary -> printf -> scanf -> binary identity
- */
-void chkid( print, scan, string )
-unsigned short print[], scan[];
-char *string;
-{
-/* Test printf-scanf identity */
-if( ecmp( print, scan ) != 0 )
- {
- pvec( print, NE );
- printf( " ->printf-> %s ->scanf->\n", string );
- pvec( scan, NE );
- printf( " is not an identity.\n" );
- ++identerr;
- }
-}
-
-
-/* Check scanf result
- */
-void chkscan( ref, tst, string )
-unsigned short ref[], tst[];
-char *string;
-{
-/* Test scanf() */
-if( ecmp( ref, tst ) != 0 )
- {
- printf( "scanf(%s) -> ", string );
- pvec( tst, NE );
- printf( "\n should be " );
- pvec( ref, NE );
- printf( ".\n" );
- ++errscan;
- ++errtot;
- }
-}
-
-
-/* Test printf() result
- */
-void chkprint( ref, tst, string )
-unsigned short ref[], tst[];
-char *string;
-{
-if( ecmp(ref, tst) != 0 )
- {
- printf( "printf( ");
- pvec( ref, NE );
- printf( ") -> %s\n", string );
- printf( " = " );
- pvec( tst, NE );
- printf( ".\n" );
- ++errprint;
- ++errtot;
- }
-}
-
-
-/* Print array of n 16-bit shorts
- */
-void pvec( x, n )
-unsigned short x[];
-int n;
-{
-int i;
-
-for( i=0; i<n; i++ )
- {
- printf( "%04x ", x[i] );
- }
-}
-
-/* Measure worst case printf rounding error
- */
-void cmpprint( ref, tst )
-unsigned short ref[], tst[];
-{
-unsigned short e[NE];
-
-if( ecmp( ref, ezero ) != 0 )
- {
- esub( ref, tst, e );
- ediv( ref, e, e );
- eabs( e );
- if( ecmp( e, rprint ) > 0 )
- emov( e, rprint );
- }
-}
-
-/* Measure worst case scanf rounding error
- */
-void cmpscan( ref, tst )
-unsigned short ref[], tst[];
-{
-unsigned short er[NE];
-
-if( ecmp( ref, ezero ) != 0 )
- {
- esub( ref, tst, er );
- ediv( ref, er, er );
- eabs( er );
- if( ecmp( er, rscan ) > 0 )
- emov( er, rscan );
- if( ecmp( er, ehalf ) > 0 )
- {
- etoasc( tst, str1, 21 );
- printf( "Bad error: scanf(%s) = %s !\n", str0, str1 );
- }
- }
-}
-
-/* Check rounded-down decimal string output of printf
- */
-void cmptrunc( ref, tst )
-unsigned short ref[], tst[];
-{
-if( ecmp( ref, tst ) != 0 )
- {
- printf( "printf(%s%s%s, %s) -> %s\n", quo, tformat, quo, str1, str2 );
- printf( "should be %s .\n", str3 );
- errprint += 1;
- }
-}
-
-
-void shownoncrit()
-{
-
-etoasc( rprint, str0, 3 );
-printf( "Maximum relative printf error found = %s .\n", str0 );
-etoasc( rscan, str0, 3 );
-printf( "Maximum relative scanf error found = %s .\n", str0 );
-}
-
-
-
-/* Produce arguments and call comparison subroutines.
- */
-void chkit( testno )
-int testno;
-{
-unsigned short t[NE], u[NE], v[NE];
-int j;
-
-switch( mprec )
- {
-#if LDOUBLE
- case 64:
- etoe64( fullp, &prec64 );
- e64toe( &prec64, rounded );
-#if CHKINTERNAL
- e64toasc( &prec64, str1, SPREC );
- asctoe64( str1, &xprec64 );
- e64toe( &xprec64, t );
- chkinternal( rounded, t, str1 );
-#endif
-/* check printf and scanf */
- sprintf( str2, format, prec64 );
- sscanf( str2, fformat, &sprec64 );
- e64toe( &sprec64, u );
- chkid( rounded, u, str2 );
- asctoe64( str2, &ssprec64 );
- e64toe( &ssprec64, v );
- chkscan( v, u, str2 );
- chkprint( rounded, v, str2 );
- if( testno < 8 )
- break;
-/* rounding error measurement */
- etoasc( fullp, str0, 24 );
- etoe64( fullp, &ssprec64 );
- e64toe( &ssprec64, u );
- sprintf( str2, format, ssprec64 );
- asctoe( str2, t );
- cmpprint( u, t );
- sscanf( str0, fformat, &sprec64 );
- e64toe( &sprec64, t );
- cmpscan( fullp, t );
- if( testno < 8 )
- break;
-/* strings rounded to less than maximum precision */
- e64toasc( &ssprec64, str1, 24 );
- for( j=SPREC-1; j>0; j-- )
- {
- e64toasc( &ssprec64, str3, j );
- asctoe( str3, v );
- sprintf( tformat, "%s.%dLe", pct, j );
- sprintf( str2, tformat, ssprec64 );
- asctoe( str2, t );
- cmptrunc( v, t );
- }
- break;
-#endif
-#ifdef DEC
- case 56:
-#endif
- case 53:
- etoe53( fullp, &prec53 );
- e53toe( &prec53, rounded );
-#if CHKINTERNAL
- e53toasc( &prec53, str1, SPREC );
- asctoe53( str1, &xprec53 );
- e53toe( &xprec53, t );
- chkinternal( rounded, t, str1 );
-#endif
- sprintf( str2, format, prec53 );
- sscanf( str2, fformat, &sprec53 );
- e53toe( &sprec53, u );
- chkid( rounded, u, str2 );
- asctoe53( str2, &ssprec53 );
- e53toe( &ssprec53, v );
- chkscan( v, u, str2 );
- chkprint( rounded, v, str2 );
- if( testno < 8 )
- break;
-/* rounding error measurement */
- etoasc( fullp, str0, 24 );
- etoe53( fullp, &ssprec53 );
- e53toe( &ssprec53, u );
- sprintf( str2, format, ssprec53 );
- asctoe( str2, t );
- cmpprint( u, t );
- sscanf( str0, fformat, &sprec53 );
- e53toe( &sprec53, t );
- cmpscan( fullp, t );
- if( testno < 8 )
- break;
- e53toasc( &ssprec53, str1, 24 );
- for( j=SPREC-1; j>0; j-- )
- {
- e53toasc( &ssprec53, str3, j );
- asctoe( str3, v );
- sprintf( tformat, "%s.%de", pct, j );
- sprintf( str2, tformat, ssprec53 );
- asctoe( str2, t );
- cmptrunc( v, t );
- }
- break;
-
- case 24:
- etoe24( fullp, &prec24 );
- e24toe( &prec24, rounded );
-#if CHKINTERNAL
- e24toasc( &prec24, str1, SPREC );
- asctoe24( str1, &xprec24 );
- e24toe( &xprec24, t );
- chkinternal( rounded, t, str1 );
-#endif
- sprintf( str2, format, prec24 );
- sscanf( str2, fformat, &sprec24 );
- e24toe( &sprec24, u );
- chkid( rounded, u, str2 );
- asctoe24( str2, &ssprec24 );
- e24toe( &ssprec24, v );
- chkscan( v, u, str2 );
- chkprint( rounded, v, str2 );
- if( testno < 8 )
- break;
-/* rounding error measurement */
- etoasc( fullp, str0, 24 );
- etoe24( fullp, &ssprec24 );
- e24toe( &ssprec24, u );
- sprintf( str2, format, ssprec24 );
- asctoe( str2, t );
- cmpprint( u, t );
- sscanf( str0, fformat, &sprec24 );
- e24toe( &sprec24, t );
- cmpscan( fullp, t );
-/*
- if( testno < 8 )
- break;
-*/
- e24toasc( &ssprec24, str1, 24 );
- for( j=SPREC-1; j>0; j-- )
- {
- e24toasc( &ssprec24, str3, j );
- asctoe( str3, v );
- sprintf( tformat, "%s.%de", pct, j );
- sprintf( str2, tformat, ssprec24 );
- asctoe( str2, t );
- cmptrunc( v, t );
- }
- break;
- }
-}
-
-
-void printerr()
-{
-if( (errscan == 0) && (identerr == 0) && (errprint == 0) )
- printf( "No errors found.\n" );
-else
- {
- printf( "%d binary -> decimal errors found.\n", errprint );
- printf( "%d decimal -> binary errors found.\n", errscan );
- }
-errscan = 0; /* reset for next test */
-identerr = 0;
-errprint = 0;
-errtot = 0;
-}
-
-
-/* Random number generator
- * in the range M * 10^N, where 1 <= M <= 10^17 - 1
- * and -27 <= N <= +27. Test values of M are logarithmically distributed
- * random integers; test values of N are uniformly distributed random integers.
- */
-
-static char *fwidth = "1.036163291797320557783096e1"; /* log(sqrt(10^9-1)) */
-static char *dwidth = "1.957197329044938830915E1"; /* log(sqrt(10^17-1)) */
-static char *ldwidth = "2.302585092994045684017491e1"; /* log(sqrt(10^20-1)) */
-
-static char *a13 = "13.0";
-static char *a27 = "27.0";
-static char *a34 = "34.0";
-static char *a10m13 = "1.0e-13";
-static unsigned short LOW[ NE ], WIDTH[NE], e27[NE], e10m13[NE];
-
-
-void mnrand( erand )
-unsigned short erand[];
-{
-unsigned short ea[NE], em[NE], en[NE], ex[NE];
-double x, a;
-
-if( mnrflag )
- {
- if( mnrflag == 3 )
- {
- asctoe( ldwidth, WIDTH );
- asctoe( a34, e27 );
- }
- if( mnrflag == 2 )
- {
- asctoe( dwidth, WIDTH );
- asctoe( a27, e27 );
- }
- if( mnrflag == 1 )
- {
- asctoe( fwidth, WIDTH );
- asctoe( a13, e27 );
- }
- asctoe( a10m13, e10m13 );
- mnrflag = 0;
- }
-drand( &x );
-e53toe( &x, ex ); /* x = WIDTH * ( x - 1.0 ) + LOW; */
-esub( eone, ex, ex );
-emul( WIDTH, ex, ex );
-eexp( ex, ex ); /* x = exp(x); */
-
-drand( &a );
-e53toe( &a, ea );
-emul( ea, ex, ea ); /* a = 1.0e-13 * x * a; */
-emul( e10m13, ea, ea );
-eabs( ea );
-eadd( ea, ex, ex ); /* add fuzz */
-emul( ex, ex, ex ); /* square it, to get range to 10^17 - 1 */
-efloor( ex, em ); /* this is M */
-
-/* Random power of 10 */
-drand( &a );
-e53toe( &a, ex );
-esub( eone, ex, ex ); /* y3 = 54.0 * ( y3 - 1.0 ) + 0.5; */
-emul( e27, ex, ex );
-eadd( ex, ex, ex );
-eadd( ehalf, ex, ex );
-efloor( ex, ex ); /* y3 = floor( y3 ) - 27.0; */
-esub( e27, ex, en ); /* this is N */
-epow( eten, en, ex );
-emul( ex, em, erand );
-}
-
-/* -ln 2^16382 */
-char *ldemin = "-1.1355137111933024058873097E4";
-char *ldewid = "2.2710274223866048117746193E4";
-/* -ln 2^1022 */
-char *demin = "-7.0839641853226410622441123E2";
-char *dewid = "1.4167928370645282124488225E3";
-/* -ln 2^125 */
-char *femin = "-8.6643397569993163677154015E1";
-char *fewid = "1.7328679513998632735430803E2";
-
-void etrand( erand )
-unsigned short erand[];
-{
-unsigned short ea[NE], ex[NE];
-double x, a;
-
-if( etrflag )
- {
- if( etrflag == 3 )
- {
- asctoe( ldemin, LOW );
- asctoe( ldewid, WIDTH );
- asctoe( a34, e27 );
- }
- if( etrflag == 2 )
- {
- asctoe( demin, LOW );
- asctoe( dewid, WIDTH );
- asctoe( a27, e27 );
- }
- if( etrflag == 1 )
- {
- asctoe( femin, LOW );
- asctoe( fewid, WIDTH );
- asctoe( a13, e27 );
- }
- asctoe( a10m13, e10m13 );
- etrflag = 0;
- }
-drand( &x );
-e53toe( &x, ex ); /* x = WIDTH * ( x - 1.0 ) + LOW; */
-esub( eone, ex, ex );
-emul( WIDTH, ex, ex );
-eadd( LOW, ex, ex );
-eexp( ex, ex ); /* x = exp(x); */
-
-/* add fuzz
- */
-drand( &a );
-e53toe( &a, ea );
-emul( ea, ex, ea ); /* a = 1.0e-13 * x * a; */
-emul( e10m13, ea, ea );
-if( ecmp( ex, ezero ) > 0 )
- eneg( ea );
-eadd( ea, ex, erand );
-}
-
+/* Floating point to ASCII input and output string test program.
+ *
+ * Numbers in the native machine data structure are converted
+ * to e type, then to and from decimal ASCII strings. Native
+ * printf() and scanf() functions are also used to produce
+ * and read strings. The resulting e type binary values
+ * are compared, with diagnostic printouts of any discrepancies.
+ *
+ * Steve Moshier, 16 Dec 88
+ * last revision: 16 May 92
+ */
+
+#include "ehead.h"
+#include "mconf.h"
+
+/* Include tests of 80-bit long double precision: */
+#define LDOUBLE 0
+/* Abort subtest after getting this many errors: */
+#define MAXERR 5
+/* Number of random arguments to try (set as large as you have
+ * patience for): */
+#define NRAND 100
+/* Perform internal consistency test: */
+#define CHKINTERNAL 0
+
+static unsigned short fullp[NE], rounded[NE];
+float prec24, sprec24, ssprec24;
+double prec53, sprec53, ssprec53;
+#if LDOUBLE
+long double prec64, sprec64, ssprec64;
+#endif
+
+static unsigned short rprint[NE], rscan[NE];
+static unsigned short q1[NE], q2[NE], q5[NE];
+static unsigned short e1[NE], e2[NE], e3[NE];
+static double d1, d2;
+static int errprint = 0;
+static int errscan = 0;
+static int identerr = 0;
+static int errtot = 0;
+static int count = 0;
+static char str0[80], str1[80], str2[80], str3[80];
+static unsigned short eten[NE], maxm[NE];
+
+int m, n, k2, mprec, SPREC;
+
+char *Ten = "10.0";
+char tformat[10];
+char *format24 = "%.8e";
+#ifdef DEC
+char *format53 = "%.17e";
+#else
+char *format53 = "%.16e";
+#endif
+char *fformat24 = "%e";
+char *fformat53 = "%le";
+char *pct = "%";
+char *quo = "\042";
+#if LDOUBLE
+char *format64 = "%.20Le";
+char *fformat64 = "%Le";
+#endif
+char *format;
+char *fformat;
+char *toomany = "Too many errors; aborting this test.\n";
+
+static int mnrflag;
+static int etrflag;
+void chkit(), printerr(), mnrand(), etrand(), shownoncrit();
+void chkid(), pvec();
+
+main()
+{
+int i, iprec;
+
+printf( "Steve Moshier's printf/scanf tester, version 0.2.\n\n" );
+#ifdef DEC
+ /* DEC PDP-11/VAX single precision not yet implemented */
+for( iprec = 1; iprec<2; iprec++ )
+#else
+for( iprec = 0; iprec<3; iprec++ )
+#endif
+ {
+ errscan = 0;
+ identerr = 0;
+ errprint = 0;
+ eclear( rprint );
+ eclear( rscan );
+
+switch( iprec )
+ {
+ case 0:
+ SPREC = 8; /* # digits after the decimal point */
+ mprec = 24; /* # bits in the significand */
+ m = 9; /* max # decimal digits for correct rounding */
+ n = 13; /* max power of ten for correct rounding */
+ k2 = -125; /* underflow beyond 2^-k2 */
+ format = format24; /* printf format string */
+ fformat = fformat24; /* scanf format string */
+ mnrflag = 1; /* sets interval for random numbers */
+ etrflag = 1;
+ printf( "Testing FLOAT precision.\n" );
+ break;
+
+ case 1:
+#ifdef DEC
+ SPREC = 17;
+ mprec = 56;
+ m = 17;
+ n = 27;
+ k2 = -125;
+ format = format53;
+ fformat = fformat53;
+ mnrflag = 2;
+ etrflag = 1;
+ printf( "Testing DEC DOUBLE precision.\n" );
+ break;
+#else
+ SPREC = 16;
+ mprec = 53;
+ m = 17;
+ n = 27;
+ k2 = -1021;
+ format = format53;
+ fformat = fformat53;
+ mnrflag = 2;
+ etrflag = 2;
+ printf( "Testing DOUBLE precision.\n" );
+ break;
+#endif
+ case 2:
+#if LDOUBLE
+ SPREC = 20;
+ mprec = 64;
+ m = 20;
+ n = 34;
+ k2 = -16382;
+ format = format64;
+ fformat = fformat64;
+ mnrflag = 3;
+ etrflag = 3;
+ printf( "Testing LONG DOUBLE precision.\n" );
+ break;
+#else
+ goto nodenorm;
+#endif
+ }
+
+ asctoe( Ten, eten );
+/* 10^m - 1 */
+ d2 = m;
+ e53toe( &d2, e1 );
+ epow( eten, e1, maxm );
+ esub( eone, maxm, maxm );
+
+/* test 1 */
+ printf( "1. Checking 10^n - 1 for n = %d to %d.\n", -m, m );
+ emov( eone, q5 );
+ for( count=0; count<=m; count++ )
+ {
+ esub( eone, q5, fullp );
+ chkit( 1 );
+ ediv( q5, eone, q2 );
+ esub( eone, q2, fullp );
+ chkit( 1 );
+ emul( eten, q5, q5 );
+ if( errtot >= MAXERR )
+ {
+ printf( "%s", toomany );
+ goto end1;
+ }
+ }
+end1:
+ printerr();
+
+
+/* test 2 */
+ printf( "2. Checking powers of 10 from 10^-%d to 10^%d.\n", n, n );
+ emov( eone, q5 );
+ for( count=0; count<=n; count++ )
+ {
+ emov( q5, fullp );
+ chkit( 2 );
+ ediv( q5, eone, fullp );
+ chkit( 2 );
+ emul( eten, q5, q5 );
+ if( errtot >= MAXERR )
+ {
+ printf( "%s", toomany );
+ goto end2;
+ }
+ }
+end2:
+ printerr();
+
+/* test 3 */
+ printf( "3. Checking (10^%d-1)*10^n from n = -%d to %d.\n", m, n, n );
+ emov( eone, q5 );
+ for( count= -n; count<=n; count++ )
+ {
+ emul( maxm, q5, fullp );
+ chkit( 3 );
+ emov( q5, fullp );
+ ediv( fullp, eone, fullp );
+ emul( maxm, fullp, fullp );
+ chkit( 3 );
+ emul( eten, q5, q5 );
+ if( errtot >= MAXERR )
+ {
+ printf( "%s", toomany );
+ goto end3;
+ }
+ }
+end3:
+ printerr();
+
+
+
+/* test 4 */
+ printf( "4. Checking powers of 2 from 2^-24 to 2^+56.\n" );
+ d1 = -24.0;
+ e53toe( &d1, q1 );
+ epow( etwo, q1, q5 );
+
+ for( count = -24; count <= 56; count++ )
+ {
+ emov( q5, fullp );
+ chkit( 4 );
+ emul( etwo, q5, q5 );
+ if( errtot >= MAXERR )
+ {
+ printf( "%s", toomany );
+ goto end4;
+ }
+ }
+end4:
+ printerr();
+
+
+/* test 5 */
+ printf( "5. Checking 2^n - 1 for n = 0 to %d.\n", mprec );
+ emov( eone, q5 );
+ for( count=0; count<=mprec; count++ )
+ {
+ esub( eone, q5, fullp );
+ chkit( 5 );
+ emul( etwo, q5, q5 );
+ if( errtot >= MAXERR )
+ {
+ printf( "%s", toomany );
+ goto end5;
+ }
+ }
+end5:
+ printerr();
+
+/* test 6 */
+ printf( "6. Checking 2^n + 1 for n = 0 to %d.\n", mprec );
+ emov( eone, q5 );
+ for( count=0; count<=mprec; count++ )
+ {
+ eadd( eone, q5, fullp );
+ chkit( 6 );
+ emul( etwo, q5, q5 );
+ if( errtot >= MAXERR )
+ {
+ printf( "%s", toomany );
+ goto end6;
+ }
+ }
+end6:
+ printerr();
+
+/* test 7 */
+ printf(
+ "7. Checking %d values M * 10^N with random integer M and N,\n",
+ NRAND );
+ printf(" 1 <= M <= 10^%d - 1 and -%d <= N <= +%d.\n", m, n, n );
+ for( i=0; i<NRAND; i++ )
+ {
+ mnrand( fullp );
+ chkit( 7 );
+ if( errtot >= MAXERR )
+ {
+ printf( "%s", toomany );
+ goto end7;
+ }
+ }
+end7:
+ printerr();
+
+/* test 8 */
+ printf("8. Checking critical rounding cases.\n" );
+ for( i=0; i<20; i++ )
+ {
+ mnrand( fullp );
+ eabs( fullp );
+ if( ecmp( fullp, eone ) < 0 )
+ ediv( fullp, eone, fullp );
+ efloor( fullp, fullp );
+ eadd( ehalf, fullp, fullp );
+ chkit( 8 );
+ if( errtot >= MAXERR )
+ {
+ printf( "%s", toomany );
+ goto end8;
+ }
+ }
+end8:
+ printerr();
+
+
+
+/* test 9 */
+ printf("9. Testing on %d random non-denormal values.\n", NRAND );
+ for( i=0; i<NRAND; i++ )
+ {
+ etrand( fullp );
+ chkit( 9 );
+ }
+ printerr();
+ shownoncrit();
+
+/* test 10 */
+ printf(
+ "Do you want to check denormal numbers in this precision ? (y/n) " );
+ gets( str0 );
+ if( str0[0] != 'y' )
+ goto nodenorm;
+
+ printf( "10. Checking denormal numbers.\n" );
+
+/* Form 2^-starting power */
+ d1 = k2;
+ e53toe( &d1, q1 );
+ epow( etwo, q1, e1 );
+
+/* Find 2^-mprec less than starting power */
+ d1 = -mprec + 4;
+ e53toe( &d1, q1 );
+ epow( etwo, q1, e3 );
+ emul( e1, e3, e3 );
+ emov( e3, e2 );
+ ediv( etwo, e2, e2 );
+
+ while( ecmp(e1,e2) != 0 )
+ {
+ eadd( e1, e2, fullp );
+ switch( mprec )
+ {
+#if LDOUBLE
+ case 64:
+ etoe64( e1, &sprec64 );
+ e64toe( &sprec64, q1 );
+ etoe64( fullp, &prec64 );
+ e64toe( &prec64, q2 );
+ break;
+#endif
+#ifdef DEC
+ case 56:
+#endif
+ case 53:
+ etoe53( e1, &sprec53 );
+ e53toe( &sprec53, q1 );
+ etoe53( fullp, &prec53 );
+ e53toe( &prec53, q2 );
+ break;
+
+ case 24:
+ etoe24( e1, &sprec24 );
+ e24toe( &sprec24, q1 );
+ etoe24( fullp, &prec24 );
+ e24toe( &prec24, q2 );
+ break;
+ }
+ if( ecmp( q2, ezero ) == 0 )
+ goto maxden;
+ chkit(10);
+ if( ecmp(q1,q2) == 0 )
+ {
+ ediv( etwo, e1, e1 );
+ emov( e3, e2 );
+ }
+ if( errtot >= MAXERR )
+ {
+ printf( "%s", toomany );
+ goto maxden;
+ }
+ ediv( etwo, e2, e2 );
+ }
+maxden:
+ printerr();
+nodenorm:
+ printf( "\n" );
+ } /* loop on precision */
+printf( "End of test.\n" );
+}
+
+#if CHKINTERNAL
+long double xprec64;
+double xprec53;
+float xprec24;
+
+/* Check binary -> printf -> scanf -> binary identity
+ * of internal routines
+ */
+void chkinternal( ref, tst, string )
+unsigned short ref[], tst[];
+char *string;
+{
+
+if( ecmp(ref,tst) != 0 )
+ {
+ printf( "internal identity compare error!\n" );
+ chkid( ref, tst, string );
+ }
+}
+#endif
+
+
+/* Check binary -> printf -> scanf -> binary identity
+ */
+void chkid( print, scan, string )
+unsigned short print[], scan[];
+char *string;
+{
+/* Test printf-scanf identity */
+if( ecmp( print, scan ) != 0 )
+ {
+ pvec( print, NE );
+ printf( " ->printf-> %s ->scanf->\n", string );
+ pvec( scan, NE );
+ printf( " is not an identity.\n" );
+ ++identerr;
+ }
+}
+
+
+/* Check scanf result
+ */
+void chkscan( ref, tst, string )
+unsigned short ref[], tst[];
+char *string;
+{
+/* Test scanf() */
+if( ecmp( ref, tst ) != 0 )
+ {
+ printf( "scanf(%s) -> ", string );
+ pvec( tst, NE );
+ printf( "\n should be " );
+ pvec( ref, NE );
+ printf( ".\n" );
+ ++errscan;
+ ++errtot;
+ }
+}
+
+
+/* Test printf() result
+ */
+void chkprint( ref, tst, string )
+unsigned short ref[], tst[];
+char *string;
+{
+if( ecmp(ref, tst) != 0 )
+ {
+ printf( "printf( ");
+ pvec( ref, NE );
+ printf( ") -> %s\n", string );
+ printf( " = " );
+ pvec( tst, NE );
+ printf( ".\n" );
+ ++errprint;
+ ++errtot;
+ }
+}
+
+
+/* Print array of n 16-bit shorts
+ */
+void pvec( x, n )
+unsigned short x[];
+int n;
+{
+int i;
+
+for( i=0; i<n; i++ )
+ {
+ printf( "%04x ", x[i] );
+ }
+}
+
+/* Measure worst case printf rounding error
+ */
+void cmpprint( ref, tst )
+unsigned short ref[], tst[];
+{
+unsigned short e[NE];
+
+if( ecmp( ref, ezero ) != 0 )
+ {
+ esub( ref, tst, e );
+ ediv( ref, e, e );
+ eabs( e );
+ if( ecmp( e, rprint ) > 0 )
+ emov( e, rprint );
+ }
+}
+
+/* Measure worst case scanf rounding error
+ */
+void cmpscan( ref, tst )
+unsigned short ref[], tst[];
+{
+unsigned short er[NE];
+
+if( ecmp( ref, ezero ) != 0 )
+ {
+ esub( ref, tst, er );
+ ediv( ref, er, er );
+ eabs( er );
+ if( ecmp( er, rscan ) > 0 )
+ emov( er, rscan );
+ if( ecmp( er, ehalf ) > 0 )
+ {
+ etoasc( tst, str1, 21 );
+ printf( "Bad error: scanf(%s) = %s !\n", str0, str1 );
+ }
+ }
+}
+
+/* Check rounded-down decimal string output of printf
+ */
+void cmptrunc( ref, tst )
+unsigned short ref[], tst[];
+{
+if( ecmp( ref, tst ) != 0 )
+ {
+ printf( "printf(%s%s%s, %s) -> %s\n", quo, tformat, quo, str1, str2 );
+ printf( "should be %s .\n", str3 );
+ errprint += 1;
+ }
+}
+
+
+void shownoncrit()
+{
+
+etoasc( rprint, str0, 3 );
+printf( "Maximum relative printf error found = %s .\n", str0 );
+etoasc( rscan, str0, 3 );
+printf( "Maximum relative scanf error found = %s .\n", str0 );
+}
+
+
+
+/* Produce arguments and call comparison subroutines.
+ */
+void chkit( testno )
+int testno;
+{
+unsigned short t[NE], u[NE], v[NE];
+int j;
+
+switch( mprec )
+ {
+#if LDOUBLE
+ case 64:
+ etoe64( fullp, &prec64 );
+ e64toe( &prec64, rounded );
+#if CHKINTERNAL
+ e64toasc( &prec64, str1, SPREC );
+ asctoe64( str1, &xprec64 );
+ e64toe( &xprec64, t );
+ chkinternal( rounded, t, str1 );
+#endif
+/* check printf and scanf */
+ sprintf( str2, format, prec64 );
+ sscanf( str2, fformat, &sprec64 );
+ e64toe( &sprec64, u );
+ chkid( rounded, u, str2 );
+ asctoe64( str2, &ssprec64 );
+ e64toe( &ssprec64, v );
+ chkscan( v, u, str2 );
+ chkprint( rounded, v, str2 );
+ if( testno < 8 )
+ break;
+/* rounding error measurement */
+ etoasc( fullp, str0, 24 );
+ etoe64( fullp, &ssprec64 );
+ e64toe( &ssprec64, u );
+ sprintf( str2, format, ssprec64 );
+ asctoe( str2, t );
+ cmpprint( u, t );
+ sscanf( str0, fformat, &sprec64 );
+ e64toe( &sprec64, t );
+ cmpscan( fullp, t );
+ if( testno < 8 )
+ break;
+/* strings rounded to less than maximum precision */
+ e64toasc( &ssprec64, str1, 24 );
+ for( j=SPREC-1; j>0; j-- )
+ {
+ e64toasc( &ssprec64, str3, j );
+ asctoe( str3, v );
+ sprintf( tformat, "%s.%dLe", pct, j );
+ sprintf( str2, tformat, ssprec64 );
+ asctoe( str2, t );
+ cmptrunc( v, t );
+ }
+ break;
+#endif
+#ifdef DEC
+ case 56:
+#endif
+ case 53:
+ etoe53( fullp, &prec53 );
+ e53toe( &prec53, rounded );
+#if CHKINTERNAL
+ e53toasc( &prec53, str1, SPREC );
+ asctoe53( str1, &xprec53 );
+ e53toe( &xprec53, t );
+ chkinternal( rounded, t, str1 );
+#endif
+ sprintf( str2, format, prec53 );
+ sscanf( str2, fformat, &sprec53 );
+ e53toe( &sprec53, u );
+ chkid( rounded, u, str2 );
+ asctoe53( str2, &ssprec53 );
+ e53toe( &ssprec53, v );
+ chkscan( v, u, str2 );
+ chkprint( rounded, v, str2 );
+ if( testno < 8 )
+ break;
+/* rounding error measurement */
+ etoasc( fullp, str0, 24 );
+ etoe53( fullp, &ssprec53 );
+ e53toe( &ssprec53, u );
+ sprintf( str2, format, ssprec53 );
+ asctoe( str2, t );
+ cmpprint( u, t );
+ sscanf( str0, fformat, &sprec53 );
+ e53toe( &sprec53, t );
+ cmpscan( fullp, t );
+ if( testno < 8 )
+ break;
+ e53toasc( &ssprec53, str1, 24 );
+ for( j=SPREC-1; j>0; j-- )
+ {
+ e53toasc( &ssprec53, str3, j );
+ asctoe( str3, v );
+ sprintf( tformat, "%s.%de", pct, j );
+ sprintf( str2, tformat, ssprec53 );
+ asctoe( str2, t );
+ cmptrunc( v, t );
+ }
+ break;
+
+ case 24:
+ etoe24( fullp, &prec24 );
+ e24toe( &prec24, rounded );
+#if CHKINTERNAL
+ e24toasc( &prec24, str1, SPREC );
+ asctoe24( str1, &xprec24 );
+ e24toe( &xprec24, t );
+ chkinternal( rounded, t, str1 );
+#endif
+ sprintf( str2, format, prec24 );
+ sscanf( str2, fformat, &sprec24 );
+ e24toe( &sprec24, u );
+ chkid( rounded, u, str2 );
+ asctoe24( str2, &ssprec24 );
+ e24toe( &ssprec24, v );
+ chkscan( v, u, str2 );
+ chkprint( rounded, v, str2 );
+ if( testno < 8 )
+ break;
+/* rounding error measurement */
+ etoasc( fullp, str0, 24 );
+ etoe24( fullp, &ssprec24 );
+ e24toe( &ssprec24, u );
+ sprintf( str2, format, ssprec24 );
+ asctoe( str2, t );
+ cmpprint( u, t );
+ sscanf( str0, fformat, &sprec24 );
+ e24toe( &sprec24, t );
+ cmpscan( fullp, t );
+/*
+ if( testno < 8 )
+ break;
+*/
+ e24toasc( &ssprec24, str1, 24 );
+ for( j=SPREC-1; j>0; j-- )
+ {
+ e24toasc( &ssprec24, str3, j );
+ asctoe( str3, v );
+ sprintf( tformat, "%s.%de", pct, j );
+ sprintf( str2, tformat, ssprec24 );
+ asctoe( str2, t );
+ cmptrunc( v, t );
+ }
+ break;
+ }
+}
+
+
+void printerr()
+{
+if( (errscan == 0) && (identerr == 0) && (errprint == 0) )
+ printf( "No errors found.\n" );
+else
+ {
+ printf( "%d binary -> decimal errors found.\n", errprint );
+ printf( "%d decimal -> binary errors found.\n", errscan );
+ }
+errscan = 0; /* reset for next test */
+identerr = 0;
+errprint = 0;
+errtot = 0;
+}
+
+
+/* Random number generator
+ * in the range M * 10^N, where 1 <= M <= 10^17 - 1
+ * and -27 <= N <= +27. Test values of M are logarithmically distributed
+ * random integers; test values of N are uniformly distributed random integers.
+ */
+
+static char *fwidth = "1.036163291797320557783096e1"; /* log(sqrt(10^9-1)) */
+static char *dwidth = "1.957197329044938830915E1"; /* log(sqrt(10^17-1)) */
+static char *ldwidth = "2.302585092994045684017491e1"; /* log(sqrt(10^20-1)) */
+
+static char *a13 = "13.0";
+static char *a27 = "27.0";
+static char *a34 = "34.0";
+static char *a10m13 = "1.0e-13";
+static unsigned short LOW[ NE ], WIDTH[NE], e27[NE], e10m13[NE];
+
+
+void mnrand( erand )
+unsigned short erand[];
+{
+unsigned short ea[NE], em[NE], en[NE], ex[NE];
+double x, a;
+
+if( mnrflag )
+ {
+ if( mnrflag == 3 )
+ {
+ asctoe( ldwidth, WIDTH );
+ asctoe( a34, e27 );
+ }
+ if( mnrflag == 2 )
+ {
+ asctoe( dwidth, WIDTH );
+ asctoe( a27, e27 );
+ }
+ if( mnrflag == 1 )
+ {
+ asctoe( fwidth, WIDTH );
+ asctoe( a13, e27 );
+ }
+ asctoe( a10m13, e10m13 );
+ mnrflag = 0;
+ }
+drand( &x );
+e53toe( &x, ex ); /* x = WIDTH * ( x - 1.0 ) + LOW; */
+esub( eone, ex, ex );
+emul( WIDTH, ex, ex );
+eexp( ex, ex ); /* x = exp(x); */
+
+drand( &a );
+e53toe( &a, ea );
+emul( ea, ex, ea ); /* a = 1.0e-13 * x * a; */
+emul( e10m13, ea, ea );
+eabs( ea );
+eadd( ea, ex, ex ); /* add fuzz */
+emul( ex, ex, ex ); /* square it, to get range to 10^17 - 1 */
+efloor( ex, em ); /* this is M */
+
+/* Random power of 10 */
+drand( &a );
+e53toe( &a, ex );
+esub( eone, ex, ex ); /* y3 = 54.0 * ( y3 - 1.0 ) + 0.5; */
+emul( e27, ex, ex );
+eadd( ex, ex, ex );
+eadd( ehalf, ex, ex );
+efloor( ex, ex ); /* y3 = floor( y3 ) - 27.0; */
+esub( e27, ex, en ); /* this is N */
+epow( eten, en, ex );
+emul( ex, em, erand );
+}
+
+/* -ln 2^16382 */
+char *ldemin = "-1.1355137111933024058873097E4";
+char *ldewid = "2.2710274223866048117746193E4";
+/* -ln 2^1022 */
+char *demin = "-7.0839641853226410622441123E2";
+char *dewid = "1.4167928370645282124488225E3";
+/* -ln 2^125 */
+char *femin = "-8.6643397569993163677154015E1";
+char *fewid = "1.7328679513998632735430803E2";
+
+void etrand( erand )
+unsigned short erand[];
+{
+unsigned short ea[NE], ex[NE];
+double x, a;
+
+if( etrflag )
+ {
+ if( etrflag == 3 )
+ {
+ asctoe( ldemin, LOW );
+ asctoe( ldewid, WIDTH );
+ asctoe( a34, e27 );
+ }
+ if( etrflag == 2 )
+ {
+ asctoe( demin, LOW );
+ asctoe( dewid, WIDTH );
+ asctoe( a27, e27 );
+ }
+ if( etrflag == 1 )
+ {
+ asctoe( femin, LOW );
+ asctoe( fewid, WIDTH );
+ asctoe( a13, e27 );
+ }
+ asctoe( a10m13, e10m13 );
+ etrflag = 0;
+ }
+drand( &x );
+e53toe( &x, ex ); /* x = WIDTH * ( x - 1.0 ) + LOW; */
+esub( eone, ex, ex );
+emul( WIDTH, ex, ex );
+eadd( LOW, ex, ex );
+eexp( ex, ex ); /* x = exp(x); */
+
+/* add fuzz
+ */
+drand( &a );
+e53toe( &a, ea );
+emul( ea, ex, ea ); /* a = 1.0e-13 * x * a; */
+emul( e10m13, ea, ea );
+if( ecmp( ex, ezero ) > 0 )
+ eneg( ea );
+eadd( ea, ex, erand );
+}
+
generated by cgit v1.2.3 (git 2.25.1) at 2024-05-16 07:14:41 +0000