Re: News: zsh-4.1.1-nonstop-fp

[Clint Adams <clint@xxxxxxx>]

> If the changes meet with zsh developer approval, I'd like to see them
> made standard in the next release of zsh.

The diff from 4.1.1, autogenerated files removed, whitespace ignored,
README removed.  Doesn't apply cleanly to HEAD (rejects on Src/init.c
and zshconfig.ac); two simple tweaks are necessary.

diff -Nurb zsh-4.1.1/ksh93-test.sh zsh-4.1.1-nonstop-fp/ksh93-test.sh
--- zsh-4.1.1/ksh93-test.sh	1969-12-31 19:00:00.000000000 -0500
+++ zsh-4.1.1-nonstop-fp/ksh93-test.sh	2004-02-27 21:07:12.000000000 -0500
@@ -0,0 +1,396 @@
+#! /usr/local/bin/ksh -
+echo
+echo "This is a test of features of ksh93 and zsh"
+echo
+
+if test -n "$ZSH_VERSION"		# zsh only: need to force load of math function features
+then
+	zmodload zsh/mathfunc || true
+fi
+
+fpmacheps()
+{
+	# set machine epsilon in global variable eps on return
+
+	typeset x
+
+	eps=$1
+	x=$1
+	while test $(( (x + eps / 2.0) != x )) -eq 1
+	do
+		eps=$(( eps / 2.0 ))
+	done
+}
+
+fpprec()
+{
+	echo ------------------------------------------------------------------------
+	echo What is floating-point precision\?
+	echo
+
+	fpmacheps 1.0
+
+	printf "Machine epsilon = %.2e\t%.2f decimal digits\n" $eps $(( -log(eps)/log(10) ))
+}
+
+fpinfinity()
+{
+	echo ------------------------------------------------------------------------
+	echo Is floating-point infinity supported\?
+	echo
+
+	typeset x y Infinity
+
+	x=2.0
+	Infinity=$(( x ** 20000))
+	printf "Infinity by 2**20000 = %.2e\n" $Infinity
+	echo
+
+	x=1.0
+	y=0.0
+	Infinity=$(( x /y ))
+	printf "Infinity by 1/0 = %.2e\n" $Infinity
+	test $((Infinity == Infinity)) -eq 1 && printf "Infinity == Infinity is true (OK)"
+	test $((Infinity == Infinity)) -ne 1 && printf "Infinity == Infinity is false (WRONG)"
+	echo
+}
+
+fpnan()
+{
+	echo ------------------------------------------------------------------------
+	echo Is floating-point NaN supported\?
+	echo
+
+	typeset x y NaN
+
+	x=0.0
+	y=0.0
+	NaN=$(( x / y))
+	echo   "echo says:   NaN =" $NaN
+	printf "printf says: NaN = %.2e\n" $NaN
+	test $((NaN != NaN)) -eq 1 && printf "NaN != NaN is true (OK)"
+	test $((NaN != NaN)) -ne 1 && printf "NaN != NaN is false (WRONG)"
+	echo
+}
+
+fpoverflow()
+{
+	echo ------------------------------------------------------------------------
+	echo What is floating-point overflow limit\?
+	echo
+
+	typeset x y
+
+	fpmacheps 1.0
+	x=$(( 1.0 - eps / 2.0 ))
+	y=1.0
+	while test $(( (y + y) > y )) -eq 1
+	do
+		# printf "DEBUG 1: x = %.15e\n" $x
+		x=$(( x + x ))
+		y=$(( y + y ))
+		# echo   "DEBUG 2: x = " $x
+		# printf "DEBUG 3: x = %.15e\n" $x
+	done
+	printf "Overflow limit near %.15e\n" $x
+}
+
+fpsignedzero()
+{
+	echo ------------------------------------------------------------------------
+	echo Are signed zeros supported\?
+	echo
+
+	typeset x y z
+
+	x=0.0
+	printf 'x=0.0:\tmust be 0.0:\tgot %.1f\n' x
+
+	x=-0.0
+	printf 'x=-0.0:\tmust be -0.0:\tgot %.1f\n' x
+
+	x=0.0
+	x=$(( -x ))
+	printf 'x=0.0\tx=$(( -x )):\tmust be -0.0:\tgot %.1f\n' x
+
+	x=0.0
+	y=$(( sqrt(x) ))
+	printf 'x=0.0\ty=$(( sqrt(x) )):\tmust be 0.0:\tgot %.1f\n' y
+
+	x=0.0
+	y=$(( sqrt(-x) ))
+	printf 'x=0.0\ty=$(( sqrt(-x) )):\tmust be -0.0:\tgot %.1f\n' y
+
+	x=1.0
+	y=$(( 1.0 / 0.0 ))
+	z=$(( 1.0 / y ))
+	printf 'x=0.0\ty=$(( 1.0 / 0.0 ))\tz=$(( 1.0 / y )):\tmust be 0.0:\tgot %.1f\n' z
+
+	x=1.0
+	y=$(( 1.0 / 0.0 ))
+	z=$(( -1.0 / y ))
+	printf 'x=0.0\ty=$(( 1.0 / 0.0 ))\tz=$(( -1.0 / y )):\tmust be -0.0:\tgot %.1f\n' z
+
+	echo
+}
+
+fpsubnormal()
+{
+	echo ------------------------------------------------------------------------
+	echo Are IEEE 754 subnormal numbers supported\?
+	echo
+
+	typeset base s x
+
+	base=2.0 			# correct for all UNIX systems (even IBM S/390 with G5 boards)
+	fpmacheps 1.0
+
+	# printf "DEBUG: eps =         %.16a\n" $eps
+
+	x=$((1 - eps / base)) 		# x has all significand bits of one: 1.fffff...p+0
+	k=0
+	for ((s = 1.0; (s > 0.0) && (((x * s) / s) == x) ; s /= base, k--))
+	do
+		# printf "DEBUG: %d\t%.16e\n" $k $s
+		:
+	done
+	if test $(( (s == 0) || (s/2 == 0) )) -eq 1
+	then
+		printf "Arithmetic on this system underflows abruptly to zero: no subnormals, sigh...\n"
+		# printf "DEBUG: x =         %.16a\n" $x
+		printf "x = (1 - eps / %d) = %.16e\n" $base $x
+		printf "s                 = %.16e (= 2.0**%d)\n" $s $k
+		printf "%d * x * s         = %.16e\n" $base          $(( base * x * s ))
+		printf "(1 + eps) * x * s = %.16e\n"                 $(( (1 + eps) * x * s ))
+		printf "%d * x * s         = %.16e\n" 1              $(( x * s ))
+	else
+		printf "This system appears to support subnormals starting at values at or below\n"
+		printf "\t%.16e\n\n" $s
+		printf "If these really are IEEE 754 64-bit subnormals, then there should be 14 values\n"
+		printf "at multiples of 1/16, before we hit zero.  Let's find out!\n\n"
+		# k=-1022
+		for (( ; s > 0.0; s /= 16, k = k - 4 ))
+		do
+			printf "%3d\t%d\t%.16e\n" $(( 1 - (k + 1022)/4 )) $k $s
+		done
+		printf "%3d\t%d\t%.16e\n" $(( 1 - (k + 1022)/4 )) $k $s
+	fi
+
+	echo
+}
+
+fpunderflow()
+{
+	echo ------------------------------------------------------------------------
+	echo What is floating-point underflow limit\?
+	echo
+
+	typeset n x
+
+	let x=1.0
+	n=0
+
+	while test $(( (x / 2.0) > 0.0 )) -gt 0
+	do
+		x=$(( x / 2.0 ))
+		n=$(( n - 1 ))
+		# printf "DEBUG: x = %.2e\t%d\n" $x $n
+	done
+	printf "Underflow limit near %.2e\n" $x
+}
+
+brace_expansion()
+{
+	echo ------------------------------------------------------------------------
+	echo Is brace expansion supported\?
+	echo
+	echo 'echo {one,two,three}.ext produces: ' {one,two,three}.ext
+	echo
+	echo Here is how other shells handle brace expansion:
+	echo
+
+	typeset sh
+
+	for sh in /bin/sh \
+		/bin/ksh \
+		/bin/csh \
+		/bin/jsh \
+		/bin/zsh \
+		/usr/local/bin/bash \
+		/usr/local/bin/ksh \
+		/usr/local/bin/pdksh \
+		/usr/local/bin/tcsh
+	do
+		if test -x $sh
+		then
+			printf "%-24s\t" $sh
+			echo 'echo {one,two,three}.ext' | $sh
+		fi
+	done
+	echo
+}
+
+indexed_arrays()
+{
+	echo ------------------------------------------------------------------------
+	echo Are indexed arrays supported\?
+	echo
+
+	typeset i name
+
+	name=(Alice Bob Carol Dave)
+
+	echo "Zeroth name:     "'${name[0]} = '${name[0]}
+	echo "First  name:     "'${name[1]} = '${name[1]}
+	echo "Second name:     "'${name[2]} = '${name[2]}
+	echo "Third  name:     "'${name[3]} = '${name[3]}
+	echo "Number of names: "'${#name[*]} = '${#name[*]}
+
+	echo
+	echo Test of array element sizes
+	echo
+
+	echo "Length of zeroth name:     "'${#name[0]} = '${#name[0]}
+	echo "Length of first  name:     "'${#name[1]} = '${#name[1]}
+	echo "Length of second name:     "'${#name[2]} = '${#name[2]}
+	echo "Length of third  name:     "'${#name[3]} = '${#name[3]}
+	echo
+}
+
+associative_arrays()
+{
+	echo ------------------------------------------------------------------------
+	echo Are associative arrays and C-style for loops supported\?
+	echo
+
+	typeset i name
+
+	name=(Alice Bob Carol Dave)
+
+	for ((i = 0; i < 4; ++i))
+	do
+		len[${name[$i]}]=${#name[$i]}
+		echo 'len[${name['$i']}] = '${len[${name[$i]}]}
+	done
+	echo
+}
+
+fparith()
+{
+	echo ------------------------------------------------------------------------
+	echo Is floating-point arithmetic supported\?
+	echo
+
+	typeset x
+
+	x=0.0
+	while test $x -lt 5.0
+	do
+		printf "printf: %s = %.3f\t" x $x
+		echo   "echo:   x = $x"
+		x=$((x + 0.25))
+	done
+	echo
+}
+
+fpfunctions()
+{
+	echo ------------------------------------------------------------------------
+	echo Are floating-point functions supported\?
+	echo
+
+	typeset x
+
+	x=0.125
+	printf "%7s\t%7s\t%7s\t%7s\t%7s\t%7s\t%7s\n" x cos exp log sin sqrt tan
+	while test $x -le 4.0
+	do
+		printf "%7.4f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\n" \
+			$x $((cos(x))) $((exp(x))) $((log(x))) $((sin(x))) $((sqrt(x))) $((tan(x)))
+		x=$((x + 0.125))
+	done
+	echo
+
+	x=0.0
+	printf "%7s\t%7s\t%7s\t%7s\t%7s\t%7s\t%7s\n" x acos asin atan cosh sinh tanh
+	while test $x -le 1.0
+	do
+		printf "%7.4f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\n" \
+			$x $((acos(x))) $((asin(x))) $((atan(x))) $((cosh(x))) $((sinh(x))) $((tanh(x)))
+		x=$((x + 0.0625))
+	done
+	echo
+
+	test -z "$ZSH_VERSION" && return
+
+	echo ------------------------------------------------------------------------
+	echo Are zsh-extension floating-point functions supported\?
+	echo
+
+	# zsh offers additional functions that we can sample (though
+	# we don't yet try the functions of two arguments)
+
+	x=0.5
+	printf "%7s\t%7s\t%7s\t%7s\t%7s\t%7s\t%7s\n" x acosh asinh atanh cbrt erf erfc
+	while test $x -le 2.0
+	do
+		printf "%7.4f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\n" \
+			$x $((acosh(x))) $((asinh(x))) $((atanh(x))) $((cbrt(x))) $((erf(x))) $((erfc(x)))
+		x=$((x + 0.0625))
+	done
+	echo
+
+	x=0.0
+	printf "%7s\t%7s\t%7s\t%7s\t%7s\t%7s\t%7s\n" x expm1 ceil fabs floor gamma ilogb
+	while test $x -le 1.0
+	do
+		printf "%7.4f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\n" \
+			$x $((expm1(x))) $((ceil(x))) $((fabs(x))) $((floor(x))) $((gamma(x))) $((ilogb(x)))
+		x=$((x + 0.0625))
+	done
+	echo
+
+	x=0.0
+	printf "%7s\t%7s\t%7s\t%7s\t%7s\t%7s\t%7s\n" x j0 j1 lgamma log10 log1p logb
+	while test $x -le 1.0
+	do
+		printf "%7.4f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\t%7.3f\n" \
+			$x $((j0(x))) $((j1(x))) $((lgamma(x))) $((log10(x))) $((log1p(x))) $((logb(x)))
+		x=$((x + 0.0625))
+	done
+	echo
+
+	x=1.0625
+	printf "%7s\t%7s\t%7s\t%7s\t%7s\t%7s\t%7s\n" x rint y0 y1
+	while test $x -le 10.0
+	do
+		printf "%7.4f\t%7.3f\t%7.3f\t%7.3f\n" \
+			$x $((rint(x))) $((y0(x))) $((y1(x)))
+		x=$((x + 0.5))
+	done
+	echo
+
+}
+
+# Run the tests: several of them show a botched implementation and
+# understanding in ksh of floating-point arithmetic, and they have to
+# be disabled, sigh...
+
+fparith			|| true
+fpprec			|| true
+fpsubnormal
+fpunderflow		|| true
+
+if test -n "$ZSH_VERSION"		# zsh only: ksh botches these
+then
+	fpoverflow		|| true
+	fpinfinity		|| true
+	fpnan		|| true
+fi
+
+fpfunctions		|| true
+fpsignedzero		|| true
+
+brace_expansion		|| true
+indexed_arrays		|| true
+associative_arrays	|| true
diff -Nurb zsh-4.1.1/Src/init.c zsh-4.1.1-nonstop-fp/Src/init.c
--- zsh-4.1.1/Src/init.c	2003-05-15 05:25:21.000000000 -0400
+++ zsh-4.1.1-nonstop-fp/Src/init.c	2004-02-27 21:37:26.000000000 -0500
@@ -1178,6 +1178,26 @@
     return 1;
 }
 
+#if defined(NONSTOP_FP)
+#if defined(HAVE_GET_FPC_CSR)
+#include <sys/fpu.h>
+#endif
+
+static void
+flush_to_zero(int on_off)		/* see "man sigfpe" on SGI IRIX 6.x for documentation */
+{
+#if defined(HAVE_GET_FPC_CSR)
+	union fpc_csr n;
+
+	n.fc_word = get_fpc_csr();
+	n.fc_struct.flush = (on_off ? 1 : 0);
+	set_fpc_csr(n.fc_word);
+#endif
+}
+#endif
+
+
+
 /*
  * This is real main entry point. This has to be mod_export'ed
  * so zsh.exe can found it on Cygwin
@@ -1193,6 +1213,10 @@
     setlocale(LC_ALL, "");
 #endif
 
+#if defined(NONSTOP_FP)
+    flush_to_zero(0);
+#endif
+
     init_hackzero(argv, environ);
 
     /*
diff -Nurb zsh-4.1.1/Src/math.c zsh-4.1.1-nonstop-fp/Src/math.c
--- zsh-4.1.1/Src/math.c	2003-04-04 08:56:27.000000000 -0500
+++ zsh-4.1.1-nonstop-fp/Src/math.c	2004-02-27 21:38:29.000000000 -0500
@@ -31,6 +31,7 @@
 #include "math.pro"
 
 #include <math.h>
+#include <string.h>
 
 /* nonzero means we are not evaluating, just parsing */
  
@@ -189,6 +190,38 @@
 /**/
 int outputradix;
 
+
+#if !defined(HAVE_ISINF)
+/**/
+int
+(isinf)(double x)
+{
+    if ((-1.0 < x) && (x < 1.0))	/* x is small, and thus finite */
+	return (0);
+    else if ((x + x) == x)		/* only true if x == Infinity */
+	return (1);
+    else				/* must be finite (normal or subnormal), or NaN */
+	return (0);
+}
+#endif
+
+#if !defined(HAVE_ISNAN)
+/**/
+static double
+(store)(double *x)
+{
+    return (*x);
+}
+
+/**/
+int
+(isnan)(double x)
+{
+    /* (x != x) should be sufficient, but some compilers incorrectly optimize it away */
+    return (store(&x) != store(&x));
+}
+#endif
+
 /**/
 static int
 zzlex(void)
@@ -391,6 +424,24 @@
 	    }
 	/* Fall through! */
 	default:
+#if defined(NONSTOP_FP)
+	    if (strcmp(ptr-1,"NaN") == 0)
+	    {
+		yyval.type = MN_FLOAT;
+		yyval.u.d = 0.0;
+		yyval.u.d /= yyval.u.d;
+		ptr += 2;
+		return NUM;
+	    }
+	    else if (strcmp(ptr-1,"Inf") == 0)
+	    {
+		yyval.type = MN_FLOAT;
+		yyval.u.d = 0.0;
+		yyval.u.d = 1.0 / yyval.u.d;
+		ptr += 2;
+		return NUM;
+	    }
+#endif
 	    if (idigit(*--ptr) || *ptr == '.') {
 		char *nptr;
 		for (nptr = ptr; idigit(*nptr); nptr++);
@@ -603,10 +654,12 @@
 static int
 notzero(mnumber a)
 {
+#if !defined(NONSTOP_FP)
     if ((a.type & MN_INTEGER) ? a.u.l == 0 : a.u.d == 0.0) {
 	zerr("division by zero", NULL, 0);
 	return 0;
     }
+#endif
     return 1;
 }
 
diff -Nurb zsh-4.1.1/Src/Modules/mathfunc.c zsh-4.1.1-nonstop-fp/Src/Modules/mathfunc.c
--- zsh-4.1.1/Src/Modules/mathfunc.c	2002-10-15 14:00:01.000000000 -0400
+++ zsh-4.1.1-nonstop-fp/Src/Modules/mathfunc.c	2004-02-27 10:19:18.000000000 -0500
@@ -211,6 +211,7 @@
   if (errflag)
     return ret;
 
+#if !defined(NONSTOP_FP)
   if (id & 0xff00) {
       int rtst = 0;
 
@@ -253,6 +254,7 @@
 	  return ret;
       }
   }
+#endif
 
   switch (id & 0xff) {
   case MF_ABS:
diff -Nurb zsh-4.1.1/Src/params.c zsh-4.1.1-nonstop-fp/Src/params.c
--- zsh-4.1.1/Src/params.c	2003-05-06 11:39:56.000000000 -0400
+++ zsh-4.1.1-nonstop-fp/Src/params.c	2004-02-27 14:44:05.000000000 -0500
@@ -32,6 +32,13 @@
 
 #include "version.h"
 
+#if defined(NONSTOP_FP)
+#include <math.h>
+int (isinf)(double);
+int (isnan)(double);
+#endif
+
+
 /* what level of localness we are at */
  
 /**/
@@ -3463,11 +3470,20 @@
 	ret = NULL;
     } else {
 	VARARR(char, buf, 512 + digits);
+#if defined(NONSTOP_FP)
+	if (isinf(dval))
+	    ret = dupstring((dval < 0.0) ? "-Inf" : "Inf");
+	else if (isnan(dval))
+	    ret = dupstring("NaN");
+	else 
+#endif
+	{
 	sprintf(buf, fmt, digits, dval);
 	if (!strchr(buf, 'e') && !strchr(buf, '.'))
 	    strcat(buf, ".");
 	ret = dupstring(buf);
     }
+    }
 #ifdef USE_LOCALE
     if (prev_locale) setlocale(LC_NUMERIC, prev_locale);
 #endif
diff -Nurb zsh-4.1.1/zshconfig.ac zsh-4.1.1-nonstop-fp/zshconfig.ac
--- zsh-4.1.1/zshconfig.ac	2003-05-06 11:39:03.000000000 -0400
+++ zsh-4.1.1-nonstop-fp/zshconfig.ac	2004-02-27 19:57:04.000000000 -0500
@@ -962,7 +962,8 @@
 	       pcre_compile pcre_study pcre_exec \
 	       nl_langinfo \
 	       erand48 open_memstream \
-	       wctomb iconv)
+	       wctomb iconv \
+	       get_fpc_csr isinf isnan)
 AC_FUNC_STRCOLL
 
 dnl  Check if tgetent accepts NULL (and will allocate its own termcap buffer)