--- ray/src/common/calexpr.c	2010/07/25 05:50:27	2.34
+++ ray/src/common/calexpr.c	2024/09/16 17:31:14	2.51
@@ -1,5 +1,5 @@
 #ifndef lint
-static const char	RCSid[] = "$Id: calexpr.c,v 2.34 2010/07/25 05:50:27 greg Exp $";
+static const char	RCSid[] = "$Id: calexpr.c,v 2.51 2024/09/16 17:31:14 greg Exp $";
 #endif
 /*
  *  Compute data values using expression parser
@@ -19,8 +19,6 @@ static const char	RCSid[] = "$Id: calexpr.c,v 2.34 201
 
 #include "copyright.h"
 
-#include  <stdio.h>
-#include  <string.h>
 #include  <ctype.h>
 #include  <errno.h>
 #include  <math.h>
@@ -35,9 +33,11 @@ static const char	RCSid[] = "$Id: calexpr.c,v 2.34 201
 
 #define	 newnode()	(EPNODE *)ecalloc(1, sizeof(EPNODE))
 
-#define	 isdecimal(c)	(isdigit(c) || (c) == '.')
+#define	 isdecimal(c)	(isdigit(c) | ((c) == '.'))
 
-static double  euminus(EPNODE *), eargument(EPNODE *), enumber(EPNODE *);
+#define  envalue(ep)	((ep)->type==NUM ? (ep)->v.num : evalue(ep))
+
+static double  euminus(EPNODE *), enumber(EPNODE *);
 static double  echannel(EPNODE *);
 static double  eadd(EPNODE *), esubtr(EPNODE *),
                emult(EPNODE *), edivi(EPNODE *),
@@ -92,10 +92,10 @@ eparse(			/* parse an expression string */
     EPNODE  *ep;
 
     initstr(expr, NULL, 0);
-    curfunc = NULL;
+    ecurfunc = NULL;
     ep = getE1();
     if (nextc != EOF)
-	syntax("unexpected character");
+	esyntax("unexpected character");
     return(ep);
 }
 
@@ -105,20 +105,23 @@ eval(			/* evaluate an expression string */
     char  *expr
 )
 {
-    register EPNODE  *ep;
+    int  prev_support = esupport;
+    EPNODE  *ep;
     double  rval;
 
+    esupport &= ~E_RCONST;	/* don't bother reducing constant expr */
     ep = eparse(expr);
+    esupport = prev_support;	/* as you were */
     rval = evalue(ep);
-    epfree(ep);
+    epfree(ep,1);
     return(rval);
 }
 
 
 int
 epcmp(			/* compare two expressions for equivalence */
-    register EPNODE  *ep1,
-    register EPNODE  *ep2
+    EPNODE  *ep1,
+    EPNODE  *ep2
 )
 {
 	double  d;
@@ -167,10 +170,11 @@ epcmp(			/* compare two expressions for equivalence */
 
 void
 epfree(			/* free a parse tree */
-    register EPNODE	 *epar
+    EPNODE	 *epar,
+    int		frep
 )
 {
-    register EPNODE  *ep;
+    EPNODE	*ep;
 
     switch (epar->type) {
 
@@ -189,26 +193,80 @@ epfree(			/* free a parse tree */
 	    break;
 
 	default:
-	    while ((ep = epar->v.kid) != NULL) {
-		epar->v.kid = ep->sibling;
-		epfree(ep);
-	    }
+	    if (epar->nkids < 0) {
+	    	ep = epar->v.kid - epar->nkids;
+	    	while (ep > epar->v.kid)
+	    		epfree(--ep, 0);
+		efree(ep);	/* free array space */
+	    } else
+	    	while ((ep = epar->v.kid) != NULL) {
+		    epar->v.kid = ep->sibling;
+		    epfree(ep, 1);
+		}
 	    break;
 
     }
+    if (frep)
+    	efree(epar);
+    else
+    	memset(epar, 0, sizeof(EPNODE));
+}
 
-    efree((char *)epar);
+
+static void
+epflatten(			/* flatten hierarchies for '+', '*' */
+	EPNODE *epar
+)
+{
+    EPNODE	*ep;
+
+    if (epar->nkids < 0)	/* can't handle array allocations */
+    	return;
+
+    for (ep = epar->v.kid; ep != NULL; ep = ep->sibling)
+    	while (ep->type == epar->type && ep->nkids > 0) {
+	    EPNODE	*ep1 = ep->v.kid;
+	    while (ep1->sibling != NULL)
+	    	ep1 = ep1->sibling;
+	    ep1->sibling = ep->sibling;
+	    epar->nkids += ep->nkids - 1;
+	    ep1 = ep->v.kid;
+	    *ep = *ep1;
+	    efree(ep1);		/* not epfree()! */
+	}
 }
 
-				/* the following used to be a switch */
-static double
-eargument(
-    EPNODE	*ep
+
+void
+epoptimize(			/* flatten operations, lists -> arrays */
+	EPNODE	*epar
 )
 {
-    return(argument(ep->v.chan));
+    EPNODE	*ep;
+
+    if ((epar->type == '+') | (epar->type == '*'))
+    	epflatten(epar);	/* flatten associative operations */
+
+    if (epar->nkids)		/* do children if any */
+    	for (ep = epar->v.kid; ep != NULL; ep = ep->sibling)
+	    epoptimize(ep);
+
+    if (epar->nkids > 4) {	/* make list into array if > 4 kids */
+        int	n = 1;
+    	epar->v.kid = (EPNODE *)erealloc(epar->v.kid,
+    					sizeof(EPNODE)*epar->nkids);
+    	while (n < epar->nkids) {
+	    ep = epar->v.kid[n-1].sibling;
+	    epar->v.kid[n] = *ep;
+	    efree(ep);		/* not epfree()! */
+	    epar->v.kid[n-1].sibling = epar->v.kid + n;
+	    n++;
+	}
+	epar->nkids = -n;
+    }
 }
 
+				/* the following used to be a switch */
 static double
 enumber(
     EPNODE	*ep
@@ -222,7 +280,7 @@ euminus(
     EPNODE	*ep
 )
 {
-    register EPNODE  *ep1 = ep->v.kid;
+    EPNODE  *ep1 = ep->v.kid;
 
     return(-evalue(ep1));
 }
@@ -240,9 +298,14 @@ eadd(
     EPNODE	*ep
 )
 {
-    register EPNODE  *ep1 = ep->v.kid;
+    double  sum = 0;
+    EPNODE  *ep1 = ep->v.kid;
 
-    return(evalue(ep1) + evalue(ep1->sibling));
+    do
+    	sum += envalue(ep1);
+    while ((ep1 = ep1->sibling) != NULL);
+
+    return(sum);
 }
 
 static double
@@ -250,9 +313,10 @@ esubtr(
     EPNODE	*ep
 )
 {
-    register EPNODE  *ep1 = ep->v.kid;
+    EPNODE  *ep1 = ep->v.kid;
+    EPNODE  *ep2 = ep1->sibling;
 
-    return(evalue(ep1) - evalue(ep1->sibling));
+    return(envalue(ep1) - envalue(ep2));
 }
 
 static double
@@ -260,9 +324,14 @@ emult(
     EPNODE	*ep
 )
 {
-    register EPNODE  *ep1 = ep->v.kid;
+    double  prod = 1;
+    EPNODE  *ep1 = ep->v.kid;
 
-    return(evalue(ep1) * evalue(ep1->sibling));
+    do
+    	prod *= envalue(ep1);
+    while ((ep1 = ep1->sibling) != NULL);
+
+    return(prod);
 }
 
 static double
@@ -270,16 +339,15 @@ edivi(
     EPNODE	*ep
 )
 {
-    register EPNODE  *ep1 = ep->v.kid;
-    double  d;
+    EPNODE  *ep1 = ep->v.kid;
+    double  den = evalue(ep1->sibling);
 
-    d = evalue(ep1->sibling);
-    if (d == 0.0) {
+    if (den == 0.0) {
 	wputs("Division by zero\n");
 	errno = ERANGE;
 	return(0.0);
     }
-    return(evalue(ep1) / d);
+    return(envalue(ep1) / den);
 }
 
 static double
@@ -287,7 +355,7 @@ epow(
     EPNODE	*ep
 )
 {
-    register EPNODE  *ep1 = ep->v.kid;
+    EPNODE  *ep1 = ep->v.kid;
     double  d;
     int	 lasterrno;
 
@@ -295,13 +363,14 @@ epow(
     errno = 0;
     d = pow(evalue(ep1), evalue(ep1->sibling));
 #ifdef  isnan
-    if (errno == 0)
+    if (errno == 0) {
 	if (isnan(d))
 	    errno = EDOM;
 	else if (isinf(d))
 	    errno = ERANGE;
+    }
 #endif
-    if (errno == EDOM || errno == ERANGE) {
+    if ((errno == EDOM) | (errno == ERANGE)) {
 	wputs("Illegal power\n");
 	return(0.0);
     }
@@ -322,33 +391,23 @@ ebotch(
 
 EPNODE *
 ekid(			/* return pointer to a node's nth kid */
-    register EPNODE	 *ep,
-    register int  n
+    EPNODE	 *ep,
+    int  n
 )
 {
-
-    for (ep = ep->v.kid; ep != NULL; ep = ep->sibling)
-	if (--n < 0)
-	    break;
-
+    if (ep->nkids < 0) {	/* allocated array? */
+    	if (n >= -ep->nkids)
+	    return(NULL);
+    	return(ep->v.kid + n);
+    }
+    ep = ep->v.kid;		/* else get from list */
+    while (n-- > 0)
+    	if ((ep = ep->sibling) == NULL)
+		break;
     return(ep);
 }
 
 
-int
-nekids(			/* return # of kids for node ep */
-    register EPNODE	 *ep
-)
-{
-    register int  n = 0;
-
-    for (ep = ep->v.kid; ep != NULL; ep = ep->sibling)
-	n++;
-
-    return(n);
-}
-
-
 void
 initfile(		/* prepare input file */
     FILE  *fp,
@@ -364,7 +423,7 @@ initfile(		/* prepare input file */
     lineno = ln;
     linepos = 0;
     inpbuf[0] = '\0';
-    scan();
+    escan();
 }
 
 
@@ -380,7 +439,7 @@ initstr(		/* prepare input string */
     lineno = ln;
     linbuf = s;
     linepos = 0;
-    scan();
+    escan();
 }
 
 
@@ -400,9 +459,9 @@ getscanpos(	/* return current scan position */
 
 
 int
-scan(void)		/* scan next character, return literal next */
+escan(void)		/* scan next character, return literal next */
 {
-    register int  lnext = 0;
+    int  lnext = 0;
 
     do {
 	if (linbuf[linepos] == '\0')
@@ -422,13 +481,13 @@ scan(void)		/* scan next character, return literal nex
 		break;
 	}
 	if (nextc == '{') {
-	    scan();
+	    escan();
 	    while (nextc != '}')
 		if (nextc == EOF)
-		    syntax("'}' expected");
+		    esyntax("'}' expected");
 		else
-		    scan();
-	    scan();
+		    escan();
+	    escan();
 	}
     } while (isspace(nextc));
     return(lnext);
@@ -441,7 +500,7 @@ long2ascii(			      /* convert long to ascii */
 )
 {
     static char	 buf[16];
-    register char  *cp;
+    char  *cp;
     int	 neg = 0;
 
     if (l == 0)
@@ -463,13 +522,13 @@ long2ascii(			      /* convert long to ascii */
 
 
 void
-syntax(			/* report syntax error and quit */
+esyntax(			/* report syntax error and quit */
     char  *err
 )
 {
-    register int  i;
+    int  i;
 
-    if (infile != NULL || lineno != 0) {
+    if ((infile != NULL) | (lineno != 0)) {
 	if (infile != NULL) eputs(infile);
 	if (lineno != 0) {
 	    eputs(infile != NULL ? ", line " : "line ");
@@ -491,18 +550,23 @@ syntax(			/* report syntax error and quit */
 
 void
 addekid(			/* add a child to ep */
-    register EPNODE	 *ep,
-    EPNODE	*ekid
+    EPNODE	 *ep,
+    EPNODE	*ek
 )
 {
+    if (ep->nkids < 0) {
+    	eputs("Cannot add kid to EPNODE array\n");
+    	quit(1);
+    }
+    ep->nkids++;
     if (ep->v.kid == NULL)
-	ep->v.kid = ekid;
+	ep->v.kid = ek;
     else {
 	for (ep = ep->v.kid; ep->sibling != NULL; ep = ep->sibling)
 	    ;
-	ep->sibling = ekid;
+	ep->sibling = ek;
     }
-    ekid->sibling = NULL;
+    ek->sibling = NULL;		/* shouldn't be necessary */
 }
 
 
@@ -510,14 +574,14 @@ char *
 getname(void)			/* scan an identifier */
 {
     static char	 str[RMAXWORD+1];
-    register int  i, lnext;
+    int  i, lnext;
 
     lnext = nextc;
-    for (i = 0; i < RMAXWORD && isid(lnext); i++, lnext = scan())
+    for (i = 0; i < RMAXWORD && isid(lnext); i++, lnext = escan())
 	str[i] = lnext;
     str[i] = '\0';
     while (isid(lnext))		/* skip rest of name */
-	lnext = scan();
+	lnext = escan();
 
     return(str);
 }
@@ -526,13 +590,13 @@ getname(void)			/* scan an identifier */
 int
 getinum(void)			/* scan a positive integer */
 {
-    register int  n, lnext;
+    int  n, lnext;
 
     n = 0;
     lnext = nextc;
     while (isdigit(lnext)) {
 	n = n * 10 + lnext - '0';
-	lnext = scan();
+	lnext = escan();
     }
     return(n);
 }
@@ -541,37 +605,37 @@ getinum(void)			/* scan a positive integer */
 double
 getnum(void)			/* scan a positive float */
 {
-    register int  i, lnext;
+    int  i, lnext;
     char  str[RMAXWORD+1];
 
     i = 0;
     lnext = nextc;
     while (isdigit(lnext) && i < RMAXWORD) {
 	str[i++] = lnext;
-	lnext = scan();
+	lnext = escan();
     }
-    if (lnext == '.' && i < RMAXWORD) {
+    if ((lnext == '.') & (i < RMAXWORD)) {
 	str[i++] = lnext;
-	lnext = scan();
+	lnext = escan();
 	if (i == 1 && !isdigit(lnext))
-	    syntax("badly formed number");
+	    esyntax("badly formed number");
 	while (isdigit(lnext) && i < RMAXWORD) {
 	    str[i++] = lnext;
-	    lnext = scan();
+	    lnext = escan();
 	}
     }
     if ((lnext == 'e') | (lnext == 'E') && i < RMAXWORD) {
 	str[i++] = lnext;
-	lnext = scan();
+	lnext = escan();
 	if ((lnext == '-') | (lnext == '+') && i < RMAXWORD) {
 	    str[i++] = lnext;
-	    lnext = scan();
+	    lnext = escan();
 	}
 	if (!isdigit(lnext))
-	    syntax("missing exponent");
+	    esyntax("missing exponent");
 	while (isdigit(lnext) && i < RMAXWORD) {
 	    str[i++] = lnext;
-	    lnext = scan();
+	    lnext = escan();
 	}
     }
     str[i] = '\0';
@@ -584,17 +648,17 @@ EPNODE *
 getE1(void)			/* E1 -> E1 ADDOP E2 */
 				/*	 E2 */
 {
-    register EPNODE  *ep1, *ep2;
+    EPNODE  *ep1, *ep2;
 
     ep1 = getE2();
-    while (nextc == '+' || nextc == '-') {
+    while ((nextc == '+') | (nextc == '-')) {
 	ep2 = newnode();
 	ep2->type = nextc;
-	scan();
+	escan();
 	addekid(ep2, ep1);
 	addekid(ep2, getE2());
 	if (esupport&E_RCONST &&
-			ep1->type == NUM && ep1->sibling->type == NUM)
+			(ep1->type == NUM) & (ep1->sibling->type == NUM))
 		ep2 = rconst(ep2);
 	ep1 = ep2;
     }
@@ -606,29 +670,34 @@ EPNODE *
 getE2(void)			/* E2 -> E2 MULOP E3 */
 				/*	 E3 */
 {
-    register EPNODE  *ep1, *ep2;
+    EPNODE  *ep1, *ep2;
 
     ep1 = getE3();
-    while (nextc == '*' || nextc == '/') {
+    while ((nextc == '*') | (nextc == '/')) {
 	ep2 = newnode();
 	ep2->type = nextc;
-	scan();
+	escan();
 	addekid(ep2, ep1);
 	addekid(ep2, getE3());
 	if (esupport&E_RCONST) {
 		EPNODE	*ep3 = ep1->sibling;
-		if (ep1->type == NUM && ep3->type == NUM) {
+		if ((ep1->type == NUM) & (ep3->type == NUM)) {
 			ep2 = rconst(ep2);
-		} else if (ep3->type == NUM && ep3->v.num == 0) {
-			if (ep2->type == '/')
-				syntax("divide by zero constant");
-			ep1->sibling = NULL;	/* (E2 * 0) */
-			epfree(ep2);
-			ep2 = ep3;
+		} else if (ep3->type == NUM) {
+			if (ep2->type == '/') {
+				if (ep3->v.num == 0)
+					esyntax("divide by zero constant");
+				ep2->type = '*';	/* for speed */
+				ep3->v.num = 1./ep3->v.num;
+			} else if (ep3->v.num == 0) {
+				ep1->sibling = NULL;	/* (E2 * 0) */
+				epfree(ep2,1);
+				ep2 = ep3;
+			}
 		} else if (ep1->type == NUM && ep1->v.num == 0) {
-			epfree(ep3);		/* (0 * E3) or (0 / E3) */
+			epfree(ep3,1);		/* (0 * E3) or (0 / E3) */
 			ep1->sibling = NULL;
-			efree((char *)ep2);
+			efree(ep2);
 			ep2 = ep1;
 		}
 	}
@@ -642,31 +711,35 @@ EPNODE *
 getE3(void)			/* E3 -> E4 ^ E3 */
 				/*	 E4 */
 {
-	register EPNODE  *ep1, *ep2;
+	EPNODE  *ep1, *ep2;
 
 	ep1 = getE4();
 	if (nextc != '^')
 		return(ep1);
 	ep2 = newnode();
 	ep2->type = nextc;
-	scan();
+	escan();
 	addekid(ep2, ep1);
 	addekid(ep2, getE3());
 	if (esupport&E_RCONST) {
 		EPNODE	*ep3 = ep1->sibling;
-		if (ep1->type == NUM && ep3->type == NUM) {
+		if ((ep1->type == NUM) & (ep3->type == NUM)) {
 			ep2 = rconst(ep2);
 		} else if (ep1->type == NUM && ep1->v.num == 0) {
-			epfree(ep3);		/* (0 ^ E3) */
+			epfree(ep3,1);		/* (0 ^ E3) */
 			ep1->sibling = NULL;
-			efree((char *)ep2);
+			efree(ep2);
 			ep2 = ep1;
-		} else if ((ep3->type == NUM && ep3->v.num == 0) ||
+		} else if ((ep3->type == NUM && ep3->v.num == 0) |
 				(ep1->type == NUM && ep1->v.num == 1)) {
-			epfree(ep2);		/* (E4 ^ 0) or (1 ^ E3) */
-			ep2 = newnode();
+			epfree(ep2,0);		/* (E4 ^ 0) or (1 ^ E3) */
 			ep2->type = NUM;
 			ep2->v.num = 1;
+		} else if (ep3->type == NUM && ep3->v.num == 1) {
+			efree(ep3);	/* (E4 ^ 1) */
+			ep1->sibling = NULL;
+			efree(ep2);
+			ep2 = ep1;
 		}
 	}
 	return(ep2);
@@ -677,10 +750,10 @@ EPNODE *
 getE4(void)			/* E4 -> ADDOP E5 */
 				/*	 E5 */
 {
-    register EPNODE  *ep1, *ep2;
+    EPNODE  *ep1, *ep2;
 
     if (nextc == '-') {
-	scan();
+	escan();
 	ep2 = getE5();
 	if (ep2->type == NUM) {
 		ep2->v.num = -ep2->v.num;
@@ -688,7 +761,7 @@ getE4(void)			/* E4 -> ADDOP E5 */
 	}
 	if (ep2->type == UMINUS) {	/* don't generate -(-E5) */
 	    ep1 = ep2->v.kid;
-	    efree((char *)ep2);
+	    efree(ep2);
 	    return(ep1);
 	}
 	ep1 = newnode();
@@ -697,7 +770,7 @@ getE4(void)			/* E4 -> ADDOP E5 */
 	return(ep1);
     }
     if (nextc == '+')
-	scan();
+	escan();
     return(getE5());
 }
 
@@ -712,32 +785,30 @@ getE5(void)			/* E5 -> (E1) */
 {
 	int	 i;
 	char  *nam;
-	register EPNODE  *ep1, *ep2;
+	EPNODE  *ep1, *ep2;
 
 	if (nextc == '(') {
-		scan();
+		escan();
 		ep1 = getE1();
 		if (nextc != ')')
-			syntax("')' expected");
-		scan();
+			esyntax("')' expected");
+		escan();
 		return(ep1);
 	}
-
 	if (esupport&E_INCHAN && nextc == '$') {
-		scan();
+		escan();
 		ep1 = newnode();
 		ep1->type = CHAN;
 		ep1->v.chan = getinum();
 		return(ep1);
 	}
-
 	if (esupport&(E_VARIABLE|E_FUNCTION) &&
-			(isalpha(nextc) || nextc == CNTXMARK)) {
+			(isalpha(nextc) | (nextc == CNTXMARK))) {
 		nam = getname();
 		ep1 = NULL;
 		if ((esupport&(E_VARIABLE|E_FUNCTION)) == (E_VARIABLE|E_FUNCTION)
-				&& curfunc != NULL)
-			for (i = 1, ep2 = curfunc->v.kid->sibling;
+				&& ecurfunc != NULL)
+			for (i = 1, ep2 = ecurfunc->v.kid->sibling;
 					ep2 != NULL; i++, ep2 = ep2->sibling)
 				if (!strcmp(ep2->v.name, nam)) {
 					ep1 = newnode();
@@ -756,44 +827,43 @@ getE5(void)			/* E5 -> (E1) */
 			addekid(ep2, ep1);
 			ep1 = ep2;
 			do {
-				scan();
+				escan();
 				addekid(ep1, getE1());
 			} while (nextc == ',');
 			if (nextc != ')')
-				syntax("')' expected");
-			scan();
+				esyntax("')' expected");
+			escan();
 		} else if (!(esupport&E_VARIABLE))
-			syntax("'(' expected");
+			esyntax("'(' expected");
 		if (esupport&E_RCONST && isconstvar(ep1))
 			ep1 = rconst(ep1);
 		return(ep1);
 	}
-
 	if (isdecimal(nextc)) {
 		ep1 = newnode();
 		ep1->type = NUM;
 		ep1->v.num = getnum();
 		return(ep1);
 	}
-	syntax("unexpected character");
+	esyntax("unexpected character");
 	return NULL; /* pro forma return */
 }
 
 
 EPNODE *
 rconst(			/* reduce a constant expression */
-    register EPNODE	 *epar
+    EPNODE	 *epar
 )
 {
-    register EPNODE  *ep;
+    EPNODE  *ep;
 
     ep = newnode();
     ep->type = NUM;
     errno = 0;
     ep->v.num = evalue(epar);
-    if (errno == EDOM || errno == ERANGE)
-	syntax("bad constant expression");
-    epfree(epar);
+    if ((errno == EDOM) | (errno == ERANGE))
+	esyntax("bad constant expression");
+    epfree(epar,1);
  
     return(ep);
 }
@@ -801,10 +871,10 @@ rconst(			/* reduce a constant expression */
 
 int
 isconstvar(			/* is ep linked to a constant expression? */
-    register EPNODE	 *ep
+    EPNODE	 *ep
 )
 {
-    register EPNODE  *ep1;
+    EPNODE  *ep1;
 
     if (esupport&E_FUNCTION && ep->type == FUNC) {
 	if (!isconstfun(ep->v.kid))
@@ -827,11 +897,11 @@ isconstvar(			/* is ep linked to a constant expression
 
 int
 isconstfun(			/* is ep linked to a constant function? */
-    register EPNODE	 *ep
+    EPNODE	 *ep
 )
 {
-    register EPNODE  *dp;
-    register LIBR  *lp;
+    EPNODE  *dp;
+    ELIBR  *lp;
 
     if (ep->type != VAR)
 	return(0);