From 07582c158fcd4793f4c068dbaf291682e74c150b Mon Sep 17 00:00:00 2001 From: "Eric S. Raymond" Date: Mon, 7 Nov 1994 16:51:20 +0000 Subject: [PATCH] Added undocumented LN function. --- cupl.doc | 10 +- cupl.h | 16 ++- execute.c | 84 +++++++++++++++ monitor.c | 301 ++++++++++++++++++++++++++++++++++++++++++++++++++---- 4 files changed, 388 insertions(+), 23 deletions(-) diff --git a/cupl.doc b/cupl.doc index b2700aa..46a3c3b 100644 --- a/cupl.doc +++ b/cupl.doc @@ -8,6 +8,12 @@ references (by section-page number in the original manual) to notes which follow the appendix transcription. These notes are also excerpts from the manual. +There is one correction in the text. The original manual listed both LN +and LOG as built-in functions and wrote "LOG(a)" is "natural log of a". We +believe this is incorrect; we have changed the "LOG(a)" in the original to +LN(a) and inserted a new "LOG(a)" entry. This implementation's LOG function +is, accordingly, log_10() and not log_e(). + There are a few typographical changes to fit it into the ASCII character set. The differences: @@ -160,7 +166,9 @@ ELEMENTS OF THE LANGUAGE expression FLOOR(a) greatest integer not numerical valued exceeding a expression - LOG(a) natural og of a numerical valued + LN(a) natural log of a numerical valued + expression + LOG(a) log to base 10 of a numerical valued expression SQRT(a) positive square root numerical valued of a expression diff --git a/cupl.h b/cupl.h index 98d570f..5b06362 100644 --- a/cupl.h +++ b/cupl.h @@ -94,6 +94,11 @@ extern lvar *idlist; #define for_symbols(s) for (s = idlist; s; s = s->next) +/* subscripting operations */ +#define SUB(v, i, j) (v.elements + i * v.width + j) +#define SUBI(v, n) (n / v.width) +#define SUBJ(v, n) (n % v.width) + #define NOMEM "out of memory\n" /* miscellaneous */ @@ -111,6 +116,7 @@ extern void warn(char *msg, ...); extern void make_scalar(value *v, scalar n); extern value copy_value(value); +extern value allocate_value(int rank, int i, int j); extern void deallocate_value(value *); void cupl_reset_write(); @@ -130,13 +136,21 @@ extern value cupl_atan(value); extern value cupl_cos(value); extern value cupl_exp(value); extern value cupl_floor(value); +extern value cupl_ln(value); extern value cupl_log(value); extern value cupl_sqrt(value); extern value cupl_max(value, value); extern value cupl_min(value, value); extern value cupl_rand(value); -extern value cupl_sgm(value right); +extern value cupl_det(value); +extern value cupl_dot(value, value); +extern value cupl_inv(value); +extern value cupl_posmax(value); +extern value cupl_posmin(value); +extern value cupl_sgm(value); +extern value cupl_trc(value); +extern value cupl_trn(value); extern bool cupl_eq(value, value); extern bool cupl_lt(value, value); diff --git a/execute.c b/execute.c index f0a846c..d338512 100644 --- a/execute.c +++ b/execute.c @@ -234,6 +234,13 @@ static value cupl_eval(node *tree) RETURN_WRAP(tree, tree->car, tree->cdr, result) return(result); + case LN: + rightside = EVAL_WRAP(cupl_eval(tree->cdr)); + result = cupl_ln(rightside); + deallocate_value(&rightside); + RETURN_WRAP(tree, tree->car, tree->cdr, result) + return(result); + case SQRT: rightside = EVAL_WRAP(cupl_eval(tree->cdr)); result = cupl_sqrt(rightside); @@ -268,6 +275,41 @@ static value cupl_eval(node *tree) * Matrix functions */ + case DET: + rightside = EVAL_WRAP(cupl_eval(tree->cdr)); + result = cupl_det(rightside); + deallocate_value(&rightside); + RETURN_WRAP(tree, tree->car, tree->cdr, result) + return(result); + + case DOT: + leftside = EVAL_WRAP(cupl_eval(tree->car)); + rightside = EVAL_WRAP(cupl_eval(tree->cdr)); + result = cupl_dot(leftside, rightside); + RETURN_WRAP(tree, tree->car, tree->cdr, result) + return(result); + + case INV: + rightside = EVAL_WRAP(cupl_eval(tree->cdr)); + result = cupl_inv(rightside); + deallocate_value(&rightside); + RETURN_WRAP(tree, tree->car, tree->cdr, result) + return(result); + + case POSMAX: + rightside = EVAL_WRAP(cupl_eval(tree->cdr)); + result = cupl_posmax(rightside); + deallocate_value(&rightside); + RETURN_WRAP(tree, tree->car, tree->cdr, result) + return(result); + + case POSMIN: + rightside = EVAL_WRAP(cupl_eval(tree->cdr)); + result = cupl_posmin(rightside); + deallocate_value(&rightside); + RETURN_WRAP(tree, tree->car, tree->cdr, result) + return(result); + case SGM: rightside = EVAL_WRAP(cupl_eval(tree->cdr)); result = cupl_sgm(rightside); @@ -275,6 +317,20 @@ static value cupl_eval(node *tree) RETURN_WRAP(tree, tree->car, tree->cdr, result) return(result); + case TRC: + rightside = EVAL_WRAP(cupl_eval(tree->cdr)); + result = cupl_trc(rightside); + deallocate_value(&rightside); + RETURN_WRAP(tree, tree->car, tree->cdr, result) + return(result); + + case TRN: + rightside = EVAL_WRAP(cupl_eval(tree->cdr)); + result = cupl_trn(rightside); + deallocate_value(&rightside); + RETURN_WRAP(tree, tree->car, tree->cdr, result) + return(result); + /* * Relations. * @@ -345,6 +401,14 @@ static value cupl_eval(node *tree) RETURN_WRAP(tree, tree->car, tree->cdr, result) return(result); + /* + * Allocation. + */ + + case ALLOCATE: + /* FIXME: implement ALLOCATE */ + die("ALLOCATE is not implemented\n"); + /* * Control structures. */ @@ -359,6 +423,26 @@ static value cupl_eval(node *tree) result.rank = FAIL; break; + case FOR: + /* FIXME: implement FOR */ + die("FOR is not implemented\n"); + + case TIMES: + /* FIXME: implement TIMES */ + die("TIMES is not implemented\n"); + + case BLOCK: + /* FIXME: implement BLOCK */ + die("BLOCK skip is not implemented\n"); + + case GO: + /* FIXME: implement GO TO */ + die("GO TO is not implemented\n"); + + case IF: + /* FIXME: implement IF */ + die("GO TO is not implemented\n"); + case STOP: longjmp(jmpbuf, 1); /* no fall through */ diff --git a/monitor.c b/monitor.c index 879606b..8a81e55 100644 --- a/monitor.c +++ b/monitor.c @@ -9,6 +9,7 @@ SYNOPSIS void make_scalar(value *v, scalar i) void copy_value(value v); + value allocate_value(int rank, int i, int j) void deallocate_value(value *v) void cupl_reset_write() @@ -27,13 +28,27 @@ SYNOPSIS value cupl_cos(value) value cupl_exp(value) value cupl_floor(value) + value cupl_ln(value) value cupl_log(value) value cupl_sqrt(value) value cupl_max(value, value) value cupl_min(value, value) value cupl_rand(value) + bool cupl_eq(value, value) + bool cupl_lt(value, value) + bool cupl_gt(value, value) + bool cupl_le(value, value) + bool cupl_ge(value, value) + + value cupl_det(value); + value cupl_dot(value, vlue); + value cupl_inv(value); + value cupl_posmax(value); + value cupl_posmin(value); value cupl_sgm(value right) + value cupl_trc(value right) + value cupl_trn(value right) DESCRIPTION Runtime support. This is segregated from the execute() code in case anyone @@ -107,12 +122,25 @@ value copy_value(value v) return(newvalue); } +value allocate_value(int rank, int i, int j) +/* allocate a value of given shape */ +{ + value v; + + v.rank = rank; + v.width = j; v.depth = i; + v.elements = (scalar *)calloc(sizeof(scalar), i * j); + + return(v); +} + void deallocate_value(value *v) /* destroy a value copy, only if its reference count is 1 */ { (void) free(v->elements); } + /**************************************************************************** * * Output support @@ -144,6 +172,7 @@ static void needspace(int w) void cupl_scalar_write(char *name, scalar quant) /* write a numeric or skip a field in CUPL style */ { + /* FIXME: cupl_scalar_write screws up line wraps */ if (name) { needspace(2 * fieldwidth); @@ -179,7 +208,7 @@ value cupl_add(value left, value right) /* add two CUPL values */ { if (!CONGRUENT(left, right)) - die("addition failed, operands of different sizes\n"); + die("addition failed, operands of different sizes or ranks\n"); else { value result; @@ -196,7 +225,7 @@ value cupl_subtract(value left, value right) /* subtract two CUPL values */ { if (!CONGRUENT(left, right)) - die("subtract failed, operands of different sizes\n"); + die("subtract failed, operands of different sizes or ranks\n"); else { value result; @@ -212,6 +241,7 @@ value cupl_subtract(value left, value right) value cupl_multiply(value left, value right) /* multiply two CUPL values */ { + /* FIXME: implement cupl_multiply for non-scalar operands */ if (left.rank == 0 && right.rank == 0) { value result; @@ -227,6 +257,7 @@ value cupl_multiply(value left, value right) value cupl_divide(value left, value right) /* divide two CUPL values */ { + /* FIXME: implement cupl_divide for non-scalar operands */ if (left.rank == 0 && right.rank == 0) { value result; @@ -343,7 +374,7 @@ value cupl_floor(value right) } } -value cupl_log(value right) +value cupl_ln(value right) /* apply ln function */ { value result; @@ -358,6 +389,21 @@ value cupl_log(value right) } } +value cupl_log(value right) +/* apply log10 function */ +{ + value result; + + if (right.rank != 0) + die("LOG is only defined for scalar arguments\n"); + else + { + make_scalar(&result, 0); + result.elements[0] = log10(right.elements[0]); + return(result); + } +} + value cupl_sqrt(value right) /* apply square root */ { @@ -383,30 +429,28 @@ value cupl_max(value left, value right) /* apply max function */ { value result; + int n; - if (right.rank != 0 || left.rank != 0) - die("MAX is not yet supported for vectors and matrices\n"); - else - { - make_scalar(&result, 0); - result.elements[0] = max(right.elements[0], left.elements[0]); - return(result); - } + make_scalar(&result, right.elements[0]); + for (n = 0; n < right.width * right.depth; n++) + if (result.elements[0] < right.elements[n]) + result.elements[0] = right.elements[n]; + + return(result); } value cupl_min(value left, value right) /* apply min functionnot yet supported for vectors and matrices */ { value result; + int n; - if (right.rank != 0 || left.rank != 0) - die("MIN is \n"); - else - { - make_scalar(&result, 0); - result.elements[0] = min(right.elements[0], left.elements[0]); - return(result); - } + make_scalar(&result, right.elements[0]); + for (n = 0; n < right.width * right.depth; n++) + if (result.elements[0] > right.elements[n]) + result.elements[0] = right.elements[n]; + + return(result); } value cupl_rand(value right) @@ -425,13 +469,200 @@ value cupl_rand(value right) } } +/**************************************************************************** + * + * Relations + * + ****************************************************************************/ + +static bool normalize(scalar *ep1, scalar *ep2) +/* normalize a pair of scalars for conversion purposes */ +{ + int mult; + + if (*ep1 == 0 && *ep2 == 0) + return; + + /* + * Note: this is not precisely the roundoff method used by original + * CUPL, though it should (famous last words) be equivalent. + */ + if (fabs(*ep1) > fabs(*ep2)) + mult = pow(10, log10(*ep2) + 14); + else + mult = pow(10, log10(*ep1) + 14); + + *ep1 *= mult; + *ep2 *= mult; + + *ep1 = floor(*ep1); + *ep2 = floor(*ep2); +} + +bool cupl_eq(value v1, value v2) +/* test any two CUPL values for pairwise equality */ +{ + if (!CONGRUENT(v1, v2)) + die("comparison failed, operands of different sizes or ranks\n"); + else + { + int n; + + for (n = 0; n < v2.width * v2.depth; n++) + { + scalar e1 = v1.elements[n], e2 = v1.elements[n]; + + normalize(&e1, &e2); + + if (e1 != e2) + return(FALSE); + } + + return(TRUE); + } +} + +bool cupl_le(value v1, value v2) +/* are all elements of v1 less than or equal to their correspondents in v2? */ +{ + if (!CONGRUENT(v1, v2)) + die("LE failed, operands of different sizes or ranks\n"); + else + { + bool equal = TRUE; + int n; + + + for (n = 0; n < v2.width * v2.depth; n++) + { + scalar e1 = v1.elements[n], e2 = v1.elements[n]; + + normalize(&e1, &e2); + + if (e1 > e2) + return(FALSE); + } + + return(TRUE); + } +} + +bool cupl_ge(value v1, value v2) +/* are all els of v1 greater than or equal to their correspondents in v2? */ +{ + if (!CONGRUENT(v1, v2)) + die("GE failed, operands of different sizes or ranks\n"); + else + { + int n; + + for (n = 0; n < v2.width * v2.depth; n++) + { + scalar e1 = v1.elements[n], e2 = v1.elements[n]; + + normalize(&e1, &e2); + + if (e1 < e2) + return(FALSE); + } + + return(TRUE); + } +} + +bool cupl_lt(value v1, value v2) +/* strange CUPL definition #1, see appendix A */ +{ + return(cupl_le(v1, v2) && !cupl_eq(v1, v2)); +} + +bool cupl_gt(value v1, value v2) +/* strange CUPL definition #2, see appendix A */ +{ + return(cupl_ge(v1, v2) && !cupl_eq(v1, v2)); +} + /**************************************************************************** * * Matrix functions * ****************************************************************************/ -/* FIXME: most not yet implemented */ +value cupl_det(value right) +{ + /* FIXME: implement DET */ + die("the determinant function is not yet implemented"); +} + +value cupl_dot(value left, value right) +/* compute inner or dot product of two vectors */ +{ + if (!CONGRUENT(left, right) || right.rank != 1) + die("DOT failed, operands of different sizes or ranks\n"); + else + { + value result; + int n; + + make_scalar(&result, 0); + for (n = 0; n < left.width * left.depth; n++) + result.elements[0] += left.elements[n] * right.elements[n]; + return(result); + } +} + +value cupl_inv(value right) +/* compute the inverse of a matrix */ +{ + /* FIXME: implement INV */ + die("the inverse function is not yet implemented\n"); +} + +value cupl_posmax(value right) +/* row position of maximum element of argument */ +{ + if (right.rank == 0) + die("POSMAX argument is a scalar"); + else + { + value result; + int n; + scalar maxel = right.elements[0];; + + make_scalar(&result, 0); + for (n = 0; n < right.width * right.depth; n++) + if (maxel < right.elements[n]) + { + result.elements[0] = SUBI(right, n); + maxel = right.elements[0]; + } + + return(result); + } +} + +value cupl_posmin(value right) +/* row position of minimum element of argument */ +{ + if (right.rank == 0) + die("POSMIN argument is a scalar"); + else + { + value result; + int n; + scalar minel = right.elements[0];; + + make_scalar(&result, 0); + for (n = 0; n < right.width * right.depth; n++) + if (minel > right.elements[n]) + { + result.elements[0] = SUBI(right, n); + minel = right.elements[n]; + } + + return(result); + } +} value cupl_sgm(value right) /* compute the sum of the elements of a matrix */ @@ -445,6 +676,34 @@ value cupl_sgm(value right) return(result); } -/* monitor.c ends here */ +value cupl_trc(value right) +/* compute sum of elements on main diagonals */ +{ + if (right.rank != 2 || right.width != right.depth) + die("TRACE failed, operands is not a square matrix\n"); + else + { + value result; + int n; + + make_scalar(&result, 0); + for (n = 0; n < right.width; n++) + result.elements[0] += SUB(right, n, n)[0]; + return(result); + } +} + +value cupl_trn(value right) +/* compute the transpose of a matrix */ +{ + value result = allocate_value(right.rank, right.depth, right.width); + int i, j; + for (i = 0; i < right.depth; i++) + for (j = 0; j < right.width; j++) + SUB(result, j, i)[0] = SUB(right, i, j)[0]; + return(result); +} + +/* monitor.c ends here */ -- 2.43.0