run dos2unix on these files

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 );
+}
+