alpar@1: /* glpmpl03.c */ alpar@1: alpar@1: /*********************************************************************** alpar@1: * This code is part of GLPK (GNU Linear Programming Kit). alpar@1: * alpar@1: * Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, alpar@1: * 2009, 2010 Andrew Makhorin, Department for Applied Informatics, alpar@1: * Moscow Aviation Institute, Moscow, Russia. All rights reserved. alpar@1: * E-mail: . alpar@1: * alpar@1: * GLPK is free software: you can redistribute it and/or modify it alpar@1: * under the terms of the GNU General Public License as published by alpar@1: * the Free Software Foundation, either version 3 of the License, or alpar@1: * (at your option) any later version. alpar@1: * alpar@1: * GLPK is distributed in the hope that it will be useful, but WITHOUT alpar@1: * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY alpar@1: * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public alpar@1: * License for more details. alpar@1: * alpar@1: * You should have received a copy of the GNU General Public License alpar@1: * along with GLPK. If not, see . alpar@1: ***********************************************************************/ alpar@1: alpar@1: #define _GLPSTD_ERRNO alpar@1: #define _GLPSTD_STDIO alpar@1: #include "glpenv.h" alpar@1: #include "glpmpl.h" alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * FLOATING-POINT NUMBERS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_add - floating-point addition. alpar@1: -- alpar@1: -- This routine computes the sum x + y. */ alpar@1: alpar@1: double fp_add(MPL *mpl, double x, double y) alpar@1: { if (x > 0.0 && y > 0.0 && x > + 0.999 * DBL_MAX - y || alpar@1: x < 0.0 && y < 0.0 && x < - 0.999 * DBL_MAX - y) alpar@1: error(mpl, "%.*g + %.*g; floating-point overflow", alpar@1: DBL_DIG, x, DBL_DIG, y); alpar@1: return x + y; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_sub - floating-point subtraction. alpar@1: -- alpar@1: -- This routine computes the difference x - y. */ alpar@1: alpar@1: double fp_sub(MPL *mpl, double x, double y) alpar@1: { if (x > 0.0 && y < 0.0 && x > + 0.999 * DBL_MAX + y || alpar@1: x < 0.0 && y > 0.0 && x < - 0.999 * DBL_MAX + y) alpar@1: error(mpl, "%.*g - %.*g; floating-point overflow", alpar@1: DBL_DIG, x, DBL_DIG, y); alpar@1: return x - y; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_less - floating-point non-negative subtraction. alpar@1: -- alpar@1: -- This routine computes the non-negative difference max(0, x - y). */ alpar@1: alpar@1: double fp_less(MPL *mpl, double x, double y) alpar@1: { if (x < y) return 0.0; alpar@1: if (x > 0.0 && y < 0.0 && x > + 0.999 * DBL_MAX + y) alpar@1: error(mpl, "%.*g less %.*g; floating-point overflow", alpar@1: DBL_DIG, x, DBL_DIG, y); alpar@1: return x - y; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_mul - floating-point multiplication. alpar@1: -- alpar@1: -- This routine computes the product x * y. */ alpar@1: alpar@1: double fp_mul(MPL *mpl, double x, double y) alpar@1: { if (fabs(y) > 1.0 && fabs(x) > (0.999 * DBL_MAX) / fabs(y)) alpar@1: error(mpl, "%.*g * %.*g; floating-point overflow", alpar@1: DBL_DIG, x, DBL_DIG, y); alpar@1: return x * y; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_div - floating-point division. alpar@1: -- alpar@1: -- This routine computes the quotient x / y. */ alpar@1: alpar@1: double fp_div(MPL *mpl, double x, double y) alpar@1: { if (fabs(y) < DBL_MIN) alpar@1: error(mpl, "%.*g / %.*g; floating-point zero divide", alpar@1: DBL_DIG, x, DBL_DIG, y); alpar@1: if (fabs(y) < 1.0 && fabs(x) > (0.999 * DBL_MAX) * fabs(y)) alpar@1: error(mpl, "%.*g / %.*g; floating-point overflow", alpar@1: DBL_DIG, x, DBL_DIG, y); alpar@1: return x / y; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_idiv - floating-point quotient of exact division. alpar@1: -- alpar@1: -- This routine computes the quotient of exact division x div y. */ alpar@1: alpar@1: double fp_idiv(MPL *mpl, double x, double y) alpar@1: { if (fabs(y) < DBL_MIN) alpar@1: error(mpl, "%.*g div %.*g; floating-point zero divide", alpar@1: DBL_DIG, x, DBL_DIG, y); alpar@1: if (fabs(y) < 1.0 && fabs(x) > (0.999 * DBL_MAX) * fabs(y)) alpar@1: error(mpl, "%.*g div %.*g; floating-point overflow", alpar@1: DBL_DIG, x, DBL_DIG, y); alpar@1: x /= y; alpar@1: return x > 0.0 ? floor(x) : x < 0.0 ? ceil(x) : 0.0; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_mod - floating-point remainder of exact division. alpar@1: -- alpar@1: -- This routine computes the remainder of exact division x mod y. alpar@1: -- alpar@1: -- NOTE: By definition x mod y = x - y * floor(x / y). */ alpar@1: alpar@1: double fp_mod(MPL *mpl, double x, double y) alpar@1: { double r; alpar@1: xassert(mpl == mpl); alpar@1: if (x == 0.0) alpar@1: r = 0.0; alpar@1: else if (y == 0.0) alpar@1: r = x; alpar@1: else alpar@1: { r = fmod(fabs(x), fabs(y)); alpar@1: if (r != 0.0) alpar@1: { if (x < 0.0) r = - r; alpar@1: if (x > 0.0 && y < 0.0 || x < 0.0 && y > 0.0) r += y; alpar@1: } alpar@1: } alpar@1: return r; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_power - floating-point exponentiation (raise to power). alpar@1: -- alpar@1: -- This routine computes the exponentiation x ** y. */ alpar@1: alpar@1: double fp_power(MPL *mpl, double x, double y) alpar@1: { double r; alpar@1: if (x == 0.0 && y <= 0.0 || x < 0.0 && y != floor(y)) alpar@1: error(mpl, "%.*g ** %.*g; result undefined", alpar@1: DBL_DIG, x, DBL_DIG, y); alpar@1: if (x == 0.0) goto eval; alpar@1: if (fabs(x) > 1.0 && y > +1.0 && alpar@1: +log(fabs(x)) > (0.999 * log(DBL_MAX)) / y || alpar@1: fabs(x) < 1.0 && y < -1.0 && alpar@1: +log(fabs(x)) < (0.999 * log(DBL_MAX)) / y) alpar@1: error(mpl, "%.*g ** %.*g; floating-point overflow", alpar@1: DBL_DIG, x, DBL_DIG, y); alpar@1: if (fabs(x) > 1.0 && y < -1.0 && alpar@1: -log(fabs(x)) < (0.999 * log(DBL_MAX)) / y || alpar@1: fabs(x) < 1.0 && y > +1.0 && alpar@1: -log(fabs(x)) > (0.999 * log(DBL_MAX)) / y) alpar@1: r = 0.0; alpar@1: else alpar@1: eval: r = pow(x, y); alpar@1: return r; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_exp - floating-point base-e exponential. alpar@1: -- alpar@1: -- This routine computes the base-e exponential e ** x. */ alpar@1: alpar@1: double fp_exp(MPL *mpl, double x) alpar@1: { if (x > 0.999 * log(DBL_MAX)) alpar@1: error(mpl, "exp(%.*g); floating-point overflow", DBL_DIG, x); alpar@1: return exp(x); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_log - floating-point natural logarithm. alpar@1: -- alpar@1: -- This routine computes the natural logarithm log x. */ alpar@1: alpar@1: double fp_log(MPL *mpl, double x) alpar@1: { if (x <= 0.0) alpar@1: error(mpl, "log(%.*g); non-positive argument", DBL_DIG, x); alpar@1: return log(x); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_log10 - floating-point common (decimal) logarithm. alpar@1: -- alpar@1: -- This routine computes the common (decimal) logarithm lg x. */ alpar@1: alpar@1: double fp_log10(MPL *mpl, double x) alpar@1: { if (x <= 0.0) alpar@1: error(mpl, "log10(%.*g); non-positive argument", DBL_DIG, x); alpar@1: return log10(x); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_sqrt - floating-point square root. alpar@1: -- alpar@1: -- This routine computes the square root x ** 0.5. */ alpar@1: alpar@1: double fp_sqrt(MPL *mpl, double x) alpar@1: { if (x < 0.0) alpar@1: error(mpl, "sqrt(%.*g); negative argument", DBL_DIG, x); alpar@1: return sqrt(x); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_sin - floating-point trigonometric sine. alpar@1: -- alpar@1: -- This routine computes the trigonometric sine sin(x). */ alpar@1: alpar@1: double fp_sin(MPL *mpl, double x) alpar@1: { if (!(-1e6 <= x && x <= +1e6)) alpar@1: error(mpl, "sin(%.*g); argument too large", DBL_DIG, x); alpar@1: return sin(x); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_cos - floating-point trigonometric cosine. alpar@1: -- alpar@1: -- This routine computes the trigonometric cosine cos(x). */ alpar@1: alpar@1: double fp_cos(MPL *mpl, double x) alpar@1: { if (!(-1e6 <= x && x <= +1e6)) alpar@1: error(mpl, "cos(%.*g); argument too large", DBL_DIG, x); alpar@1: return cos(x); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_atan - floating-point trigonometric arctangent. alpar@1: -- alpar@1: -- This routine computes the trigonometric arctangent atan(x). */ alpar@1: alpar@1: double fp_atan(MPL *mpl, double x) alpar@1: { xassert(mpl == mpl); alpar@1: return atan(x); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_atan2 - floating-point trigonometric arctangent. alpar@1: -- alpar@1: -- This routine computes the trigonometric arctangent atan(y / x). */ alpar@1: alpar@1: double fp_atan2(MPL *mpl, double y, double x) alpar@1: { xassert(mpl == mpl); alpar@1: return atan2(y, x); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_round - round floating-point value to n fractional digits. alpar@1: -- alpar@1: -- This routine rounds given floating-point value x to n fractional alpar@1: -- digits with the formula: alpar@1: -- alpar@1: -- round(x, n) = floor(x * 10^n + 0.5) / 10^n. alpar@1: -- alpar@1: -- The parameter n is assumed to be integer. */ alpar@1: alpar@1: double fp_round(MPL *mpl, double x, double n) alpar@1: { double ten_to_n; alpar@1: if (n != floor(n)) alpar@1: error(mpl, "round(%.*g, %.*g); non-integer second argument", alpar@1: DBL_DIG, x, DBL_DIG, n); alpar@1: if (n <= DBL_DIG + 2) alpar@1: { ten_to_n = pow(10.0, n); alpar@1: if (fabs(x) < (0.999 * DBL_MAX) / ten_to_n) alpar@1: { x = floor(x * ten_to_n + 0.5); alpar@1: if (x != 0.0) x /= ten_to_n; alpar@1: } alpar@1: } alpar@1: return x; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_trunc - truncate floating-point value to n fractional digits. alpar@1: -- alpar@1: -- This routine truncates given floating-point value x to n fractional alpar@1: -- digits with the formula: alpar@1: -- alpar@1: -- ( floor(x * 10^n) / 10^n, if x >= 0 alpar@1: -- trunc(x, n) = < alpar@1: -- ( ceil(x * 10^n) / 10^n, if x < 0 alpar@1: -- alpar@1: -- The parameter n is assumed to be integer. */ alpar@1: alpar@1: double fp_trunc(MPL *mpl, double x, double n) alpar@1: { double ten_to_n; alpar@1: if (n != floor(n)) alpar@1: error(mpl, "trunc(%.*g, %.*g); non-integer second argument", alpar@1: DBL_DIG, x, DBL_DIG, n); alpar@1: if (n <= DBL_DIG + 2) alpar@1: { ten_to_n = pow(10.0, n); alpar@1: if (fabs(x) < (0.999 * DBL_MAX) / ten_to_n) alpar@1: { x = (x >= 0.0 ? floor(x * ten_to_n) : ceil(x * ten_to_n)); alpar@1: if (x != 0.0) x /= ten_to_n; alpar@1: } alpar@1: } alpar@1: return x; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * PSEUDO-RANDOM NUMBER GENERATORS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_irand224 - pseudo-random integer in the range [0, 2^24). alpar@1: -- alpar@1: -- This routine returns a next pseudo-random integer (converted to alpar@1: -- floating-point) which is uniformly distributed between 0 and 2^24-1, alpar@1: -- inclusive. */ alpar@1: alpar@1: #define two_to_the_24 0x1000000 alpar@1: alpar@1: double fp_irand224(MPL *mpl) alpar@1: { return alpar@1: (double)rng_unif_rand(mpl->rand, two_to_the_24); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_uniform01 - pseudo-random number in the range [0, 1). alpar@1: -- alpar@1: -- This routine returns a next pseudo-random number which is uniformly alpar@1: -- distributed in the range [0, 1). */ alpar@1: alpar@1: #define two_to_the_31 ((unsigned int)0x80000000) alpar@1: alpar@1: double fp_uniform01(MPL *mpl) alpar@1: { return alpar@1: (double)rng_next_rand(mpl->rand) / (double)two_to_the_31; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_uniform - pseudo-random number in the range [a, b). alpar@1: -- alpar@1: -- This routine returns a next pseudo-random number which is uniformly alpar@1: -- distributed in the range [a, b). */ alpar@1: alpar@1: double fp_uniform(MPL *mpl, double a, double b) alpar@1: { double x; alpar@1: if (a >= b) alpar@1: error(mpl, "Uniform(%.*g, %.*g); invalid range", alpar@1: DBL_DIG, a, DBL_DIG, b); alpar@1: x = fp_uniform01(mpl); alpar@1: #if 0 alpar@1: x = a * (1.0 - x) + b * x; alpar@1: #else alpar@1: x = fp_add(mpl, a * (1.0 - x), b * x); alpar@1: #endif alpar@1: return x; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_normal01 - Gaussian random variate with mu = 0 and sigma = 1. alpar@1: -- alpar@1: -- This routine returns a Gaussian random variate with zero mean and alpar@1: -- unit standard deviation. The polar (Box-Mueller) method is used. alpar@1: -- alpar@1: -- This code is a modified version of the routine gsl_ran_gaussian from alpar@1: -- the GNU Scientific Library Version 1.0. */ alpar@1: alpar@1: double fp_normal01(MPL *mpl) alpar@1: { double x, y, r2; alpar@1: do alpar@1: { /* choose x, y in uniform square (-1,-1) to (+1,+1) */ alpar@1: x = -1.0 + 2.0 * fp_uniform01(mpl); alpar@1: y = -1.0 + 2.0 * fp_uniform01(mpl); alpar@1: /* see if it is in the unit circle */ alpar@1: r2 = x * x + y * y; alpar@1: } while (r2 > 1.0 || r2 == 0.0); alpar@1: /* Box-Muller transform */ alpar@1: return y * sqrt(-2.0 * log (r2) / r2); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fp_normal - Gaussian random variate with specified mu and sigma. alpar@1: -- alpar@1: -- This routine returns a Gaussian random variate with mean mu and alpar@1: -- standard deviation sigma. */ alpar@1: alpar@1: double fp_normal(MPL *mpl, double mu, double sigma) alpar@1: { double x; alpar@1: #if 0 alpar@1: x = mu + sigma * fp_normal01(mpl); alpar@1: #else alpar@1: x = fp_add(mpl, mu, fp_mul(mpl, sigma, fp_normal01(mpl))); alpar@1: #endif alpar@1: return x; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * SEGMENTED CHARACTER STRINGS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- create_string - create character string. alpar@1: -- alpar@1: -- This routine creates a segmented character string, which is exactly alpar@1: -- equivalent to specified character string. */ alpar@1: alpar@1: STRING *create_string alpar@1: ( MPL *mpl, alpar@1: char buf[MAX_LENGTH+1] /* not changed */ alpar@1: ) alpar@1: #if 0 alpar@1: { STRING *head, *tail; alpar@1: int i, j; alpar@1: xassert(buf != NULL); alpar@1: xassert(strlen(buf) <= MAX_LENGTH); alpar@1: head = tail = dmp_get_atom(mpl->strings, sizeof(STRING)); alpar@1: for (i = j = 0; ; i++) alpar@1: { if ((tail->seg[j++] = buf[i]) == '\0') break; alpar@1: if (j == STRSEG_SIZE) alpar@1: tail = (tail->next = dmp_get_atom(mpl->strings, sizeof(STRING))), j = 0; alpar@1: } alpar@1: tail->next = NULL; alpar@1: return head; alpar@1: } alpar@1: #else alpar@1: { STRING *str; alpar@1: xassert(strlen(buf) <= MAX_LENGTH); alpar@1: str = dmp_get_atom(mpl->strings, strlen(buf)+1); alpar@1: strcpy(str, buf); alpar@1: return str; alpar@1: } alpar@1: #endif alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- copy_string - make copy of character string. alpar@1: -- alpar@1: -- This routine returns an exact copy of segmented character string. */ alpar@1: alpar@1: STRING *copy_string alpar@1: ( MPL *mpl, alpar@1: STRING *str /* not changed */ alpar@1: ) alpar@1: #if 0 alpar@1: { STRING *head, *tail; alpar@1: xassert(str != NULL); alpar@1: head = tail = dmp_get_atom(mpl->strings, sizeof(STRING)); alpar@1: for (; str != NULL; str = str->next) alpar@1: { memcpy(tail->seg, str->seg, STRSEG_SIZE); alpar@1: if (str->next != NULL) alpar@1: tail = (tail->next = dmp_get_atom(mpl->strings, sizeof(STRING))); alpar@1: } alpar@1: tail->next = NULL; alpar@1: return head; alpar@1: } alpar@1: #else alpar@1: { xassert(mpl == mpl); alpar@1: return create_string(mpl, str); alpar@1: } alpar@1: #endif alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- compare_strings - compare one character string with another. alpar@1: -- alpar@1: -- This routine compares one segmented character strings with another alpar@1: -- and returns the result of comparison as follows: alpar@1: -- alpar@1: -- = 0 - both strings are identical; alpar@1: -- < 0 - the first string precedes the second one; alpar@1: -- > 0 - the first string follows the second one. */ alpar@1: alpar@1: int compare_strings alpar@1: ( MPL *mpl, alpar@1: STRING *str1, /* not changed */ alpar@1: STRING *str2 /* not changed */ alpar@1: ) alpar@1: #if 0 alpar@1: { int j, c1, c2; alpar@1: xassert(mpl == mpl); alpar@1: for (;; str1 = str1->next, str2 = str2->next) alpar@1: { xassert(str1 != NULL); alpar@1: xassert(str2 != NULL); alpar@1: for (j = 0; j < STRSEG_SIZE; j++) alpar@1: { c1 = (unsigned char)str1->seg[j]; alpar@1: c2 = (unsigned char)str2->seg[j]; alpar@1: if (c1 < c2) return -1; alpar@1: if (c1 > c2) return +1; alpar@1: if (c1 == '\0') goto done; alpar@1: } alpar@1: } alpar@1: done: return 0; alpar@1: } alpar@1: #else alpar@1: { xassert(mpl == mpl); alpar@1: return strcmp(str1, str2); alpar@1: } alpar@1: #endif alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- fetch_string - extract content of character string. alpar@1: -- alpar@1: -- This routine returns a character string, which is exactly equivalent alpar@1: -- to specified segmented character string. */ alpar@1: alpar@1: char *fetch_string alpar@1: ( MPL *mpl, alpar@1: STRING *str, /* not changed */ alpar@1: char buf[MAX_LENGTH+1] /* modified */ alpar@1: ) alpar@1: #if 0 alpar@1: { int i, j; alpar@1: xassert(mpl == mpl); alpar@1: xassert(buf != NULL); alpar@1: for (i = 0; ; str = str->next) alpar@1: { xassert(str != NULL); alpar@1: for (j = 0; j < STRSEG_SIZE; j++) alpar@1: if ((buf[i++] = str->seg[j]) == '\0') goto done; alpar@1: } alpar@1: done: xassert(strlen(buf) <= MAX_LENGTH); alpar@1: return buf; alpar@1: } alpar@1: #else alpar@1: { xassert(mpl == mpl); alpar@1: return strcpy(buf, str); alpar@1: } alpar@1: #endif alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- delete_string - delete character string. alpar@1: -- alpar@1: -- This routine deletes specified segmented character string. */ alpar@1: alpar@1: void delete_string alpar@1: ( MPL *mpl, alpar@1: STRING *str /* destroyed */ alpar@1: ) alpar@1: #if 0 alpar@1: { STRING *temp; alpar@1: xassert(str != NULL); alpar@1: while (str != NULL) alpar@1: { temp = str; alpar@1: str = str->next; alpar@1: dmp_free_atom(mpl->strings, temp, sizeof(STRING)); alpar@1: } alpar@1: return; alpar@1: } alpar@1: #else alpar@1: { dmp_free_atom(mpl->strings, str, strlen(str)+1); alpar@1: return; alpar@1: } alpar@1: #endif alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * SYMBOLS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- create_symbol_num - create symbol of numeric type. alpar@1: -- alpar@1: -- This routine creates a symbol, which has a numeric value specified alpar@1: -- as floating-point number. */ alpar@1: alpar@1: SYMBOL *create_symbol_num(MPL *mpl, double num) alpar@1: { SYMBOL *sym; alpar@1: sym = dmp_get_atom(mpl->symbols, sizeof(SYMBOL)); alpar@1: sym->num = num; alpar@1: sym->str = NULL; alpar@1: return sym; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- create_symbol_str - create symbol of abstract type. alpar@1: -- alpar@1: -- This routine creates a symbol, which has an abstract value specified alpar@1: -- as segmented character string. */ alpar@1: alpar@1: SYMBOL *create_symbol_str alpar@1: ( MPL *mpl, alpar@1: STRING *str /* destroyed */ alpar@1: ) alpar@1: { SYMBOL *sym; alpar@1: xassert(str != NULL); alpar@1: sym = dmp_get_atom(mpl->symbols, sizeof(SYMBOL)); alpar@1: sym->num = 0.0; alpar@1: sym->str = str; alpar@1: return sym; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- copy_symbol - make copy of symbol. alpar@1: -- alpar@1: -- This routine returns an exact copy of symbol. */ alpar@1: alpar@1: SYMBOL *copy_symbol alpar@1: ( MPL *mpl, alpar@1: SYMBOL *sym /* not changed */ alpar@1: ) alpar@1: { SYMBOL *copy; alpar@1: xassert(sym != NULL); alpar@1: copy = dmp_get_atom(mpl->symbols, sizeof(SYMBOL)); alpar@1: if (sym->str == NULL) alpar@1: { copy->num = sym->num; alpar@1: copy->str = NULL; alpar@1: } alpar@1: else alpar@1: { copy->num = 0.0; alpar@1: copy->str = copy_string(mpl, sym->str); alpar@1: } alpar@1: return copy; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- compare_symbols - compare one symbol with another. alpar@1: -- alpar@1: -- This routine compares one symbol with another and returns the result alpar@1: -- of comparison as follows: alpar@1: -- alpar@1: -- = 0 - both symbols are identical; alpar@1: -- < 0 - the first symbol precedes the second one; alpar@1: -- > 0 - the first symbol follows the second one. alpar@1: -- alpar@1: -- Note that the linear order, in which symbols follow each other, is alpar@1: -- implementation-dependent. It may be not an alphabetical order. */ alpar@1: alpar@1: int compare_symbols alpar@1: ( MPL *mpl, alpar@1: SYMBOL *sym1, /* not changed */ alpar@1: SYMBOL *sym2 /* not changed */ alpar@1: ) alpar@1: { xassert(sym1 != NULL); alpar@1: xassert(sym2 != NULL); alpar@1: /* let all numeric quantities precede all symbolic quantities */ alpar@1: if (sym1->str == NULL && sym2->str == NULL) alpar@1: { if (sym1->num < sym2->num) return -1; alpar@1: if (sym1->num > sym2->num) return +1; alpar@1: return 0; alpar@1: } alpar@1: if (sym1->str == NULL) return -1; alpar@1: if (sym2->str == NULL) return +1; alpar@1: return compare_strings(mpl, sym1->str, sym2->str); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- delete_symbol - delete symbol. alpar@1: -- alpar@1: -- This routine deletes specified symbol. */ alpar@1: alpar@1: void delete_symbol alpar@1: ( MPL *mpl, alpar@1: SYMBOL *sym /* destroyed */ alpar@1: ) alpar@1: { xassert(sym != NULL); alpar@1: if (sym->str != NULL) delete_string(mpl, sym->str); alpar@1: dmp_free_atom(mpl->symbols, sym, sizeof(SYMBOL)); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- format_symbol - format symbol for displaying or printing. alpar@1: -- alpar@1: -- This routine converts specified symbol to a charater string, which alpar@1: -- is suitable for displaying or printing. alpar@1: -- alpar@1: -- The resultant string is never longer than 255 characters. If it gets alpar@1: -- longer, it is truncated from the right and appended by dots. */ alpar@1: alpar@1: char *format_symbol alpar@1: ( MPL *mpl, alpar@1: SYMBOL *sym /* not changed */ alpar@1: ) alpar@1: { char *buf = mpl->sym_buf; alpar@1: xassert(sym != NULL); alpar@1: if (sym->str == NULL) alpar@1: sprintf(buf, "%.*g", DBL_DIG, sym->num); alpar@1: else alpar@1: { char str[MAX_LENGTH+1]; alpar@1: int quoted, j, len; alpar@1: fetch_string(mpl, sym->str, str); alpar@1: if (!(isalpha((unsigned char)str[0]) || str[0] == '_')) alpar@1: quoted = 1; alpar@1: else alpar@1: { quoted = 0; alpar@1: for (j = 1; str[j] != '\0'; j++) alpar@1: { if (!(isalnum((unsigned char)str[j]) || alpar@1: strchr("+-._", (unsigned char)str[j]) != NULL)) alpar@1: { quoted = 1; alpar@1: break; alpar@1: } alpar@1: } alpar@1: } alpar@1: # define safe_append(c) \ alpar@1: (void)(len < 255 ? (buf[len++] = (char)(c)) : 0) alpar@1: buf[0] = '\0', len = 0; alpar@1: if (quoted) safe_append('\''); alpar@1: for (j = 0; str[j] != '\0'; j++) alpar@1: { if (quoted && str[j] == '\'') safe_append('\''); alpar@1: safe_append(str[j]); alpar@1: } alpar@1: if (quoted) safe_append('\''); alpar@1: # undef safe_append alpar@1: buf[len] = '\0'; alpar@1: if (len == 255) strcpy(buf+252, "..."); alpar@1: } alpar@1: xassert(strlen(buf) <= 255); alpar@1: return buf; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- concat_symbols - concatenate one symbol with another. alpar@1: -- alpar@1: -- This routine concatenates values of two given symbols and assigns alpar@1: -- the resultant character string to a new symbol, which is returned on alpar@1: -- exit. Both original symbols are destroyed. */ alpar@1: alpar@1: SYMBOL *concat_symbols alpar@1: ( MPL *mpl, alpar@1: SYMBOL *sym1, /* destroyed */ alpar@1: SYMBOL *sym2 /* destroyed */ alpar@1: ) alpar@1: { char str1[MAX_LENGTH+1], str2[MAX_LENGTH+1]; alpar@1: xassert(MAX_LENGTH >= DBL_DIG + DBL_DIG); alpar@1: if (sym1->str == NULL) alpar@1: sprintf(str1, "%.*g", DBL_DIG, sym1->num); alpar@1: else alpar@1: fetch_string(mpl, sym1->str, str1); alpar@1: if (sym2->str == NULL) alpar@1: sprintf(str2, "%.*g", DBL_DIG, sym2->num); alpar@1: else alpar@1: fetch_string(mpl, sym2->str, str2); alpar@1: if (strlen(str1) + strlen(str2) > MAX_LENGTH) alpar@1: { char buf[255+1]; alpar@1: strcpy(buf, format_symbol(mpl, sym1)); alpar@1: xassert(strlen(buf) < sizeof(buf)); alpar@1: error(mpl, "%s & %s; resultant symbol exceeds %d characters", alpar@1: buf, format_symbol(mpl, sym2), MAX_LENGTH); alpar@1: } alpar@1: delete_symbol(mpl, sym1); alpar@1: delete_symbol(mpl, sym2); alpar@1: return create_symbol_str(mpl, create_string(mpl, strcat(str1, alpar@1: str2))); alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * N-TUPLES * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- create_tuple - create n-tuple. alpar@1: -- alpar@1: -- This routine creates a n-tuple, which initially has no components, alpar@1: -- i.e. which is 0-tuple. */ alpar@1: alpar@1: TUPLE *create_tuple(MPL *mpl) alpar@1: { TUPLE *tuple; alpar@1: xassert(mpl == mpl); alpar@1: tuple = NULL; alpar@1: return tuple; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- expand_tuple - append symbol to n-tuple. alpar@1: -- alpar@1: -- This routine expands n-tuple appending to it a given symbol, which alpar@1: -- becomes its new last component. */ alpar@1: alpar@1: TUPLE *expand_tuple alpar@1: ( MPL *mpl, alpar@1: TUPLE *tuple, /* destroyed */ alpar@1: SYMBOL *sym /* destroyed */ alpar@1: ) alpar@1: { TUPLE *tail, *temp; alpar@1: xassert(sym != NULL); alpar@1: /* create a new component */ alpar@1: tail = dmp_get_atom(mpl->tuples, sizeof(TUPLE)); alpar@1: tail->sym = sym; alpar@1: tail->next = NULL; alpar@1: /* and append it to the component list */ alpar@1: if (tuple == NULL) alpar@1: tuple = tail; alpar@1: else alpar@1: { for (temp = tuple; temp->next != NULL; temp = temp->next); alpar@1: temp->next = tail; alpar@1: } alpar@1: return tuple; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- tuple_dimen - determine dimension of n-tuple. alpar@1: -- alpar@1: -- This routine returns dimension of n-tuple, i.e. number of components alpar@1: -- in the n-tuple. */ alpar@1: alpar@1: int tuple_dimen alpar@1: ( MPL *mpl, alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { TUPLE *temp; alpar@1: int dim = 0; alpar@1: xassert(mpl == mpl); alpar@1: for (temp = tuple; temp != NULL; temp = temp->next) dim++; alpar@1: return dim; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- copy_tuple - make copy of n-tuple. alpar@1: -- alpar@1: -- This routine returns an exact copy of n-tuple. */ alpar@1: alpar@1: TUPLE *copy_tuple alpar@1: ( MPL *mpl, alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { TUPLE *head, *tail; alpar@1: if (tuple == NULL) alpar@1: head = NULL; alpar@1: else alpar@1: { head = tail = dmp_get_atom(mpl->tuples, sizeof(TUPLE)); alpar@1: for (; tuple != NULL; tuple = tuple->next) alpar@1: { xassert(tuple->sym != NULL); alpar@1: tail->sym = copy_symbol(mpl, tuple->sym); alpar@1: if (tuple->next != NULL) alpar@1: tail = (tail->next = dmp_get_atom(mpl->tuples, sizeof(TUPLE))); alpar@1: } alpar@1: tail->next = NULL; alpar@1: } alpar@1: return head; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- compare_tuples - compare one n-tuple with another. alpar@1: -- alpar@1: -- This routine compares two given n-tuples, which must have the same alpar@1: -- dimension (not checked for the sake of efficiency), and returns one alpar@1: -- of the following codes: alpar@1: -- alpar@1: -- = 0 - both n-tuples are identical; alpar@1: -- < 0 - the first n-tuple precedes the second one; alpar@1: -- > 0 - the first n-tuple follows the second one. alpar@1: -- alpar@1: -- Note that the linear order, in which n-tuples follow each other, is alpar@1: -- implementation-dependent. It may be not an alphabetical order. */ alpar@1: alpar@1: int compare_tuples alpar@1: ( MPL *mpl, alpar@1: TUPLE *tuple1, /* not changed */ alpar@1: TUPLE *tuple2 /* not changed */ alpar@1: ) alpar@1: { TUPLE *item1, *item2; alpar@1: int ret; alpar@1: xassert(mpl == mpl); alpar@1: for (item1 = tuple1, item2 = tuple2; item1 != NULL; alpar@1: item1 = item1->next, item2 = item2->next) alpar@1: { xassert(item2 != NULL); alpar@1: xassert(item1->sym != NULL); alpar@1: xassert(item2->sym != NULL); alpar@1: ret = compare_symbols(mpl, item1->sym, item2->sym); alpar@1: if (ret != 0) return ret; alpar@1: } alpar@1: xassert(item2 == NULL); alpar@1: return 0; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- build_subtuple - build subtuple of given n-tuple. alpar@1: -- alpar@1: -- This routine builds subtuple, which consists of first dim components alpar@1: -- of given n-tuple. */ alpar@1: alpar@1: TUPLE *build_subtuple alpar@1: ( MPL *mpl, alpar@1: TUPLE *tuple, /* not changed */ alpar@1: int dim alpar@1: ) alpar@1: { TUPLE *head, *temp; alpar@1: int j; alpar@1: head = create_tuple(mpl); alpar@1: for (j = 1, temp = tuple; j <= dim; j++, temp = temp->next) alpar@1: { xassert(temp != NULL); alpar@1: head = expand_tuple(mpl, head, copy_symbol(mpl, temp->sym)); alpar@1: } alpar@1: return head; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- delete_tuple - delete n-tuple. alpar@1: -- alpar@1: -- This routine deletes specified n-tuple. */ alpar@1: alpar@1: void delete_tuple alpar@1: ( MPL *mpl, alpar@1: TUPLE *tuple /* destroyed */ alpar@1: ) alpar@1: { TUPLE *temp; alpar@1: while (tuple != NULL) alpar@1: { temp = tuple; alpar@1: tuple = temp->next; alpar@1: xassert(temp->sym != NULL); alpar@1: delete_symbol(mpl, temp->sym); alpar@1: dmp_free_atom(mpl->tuples, temp, sizeof(TUPLE)); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- format_tuple - format n-tuple for displaying or printing. alpar@1: -- alpar@1: -- This routine converts specified n-tuple to a character string, which alpar@1: -- is suitable for displaying or printing. alpar@1: -- alpar@1: -- The resultant string is never longer than 255 characters. If it gets alpar@1: -- longer, it is truncated from the right and appended by dots. */ alpar@1: alpar@1: char *format_tuple alpar@1: ( MPL *mpl, alpar@1: int c, alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { TUPLE *temp; alpar@1: int dim, j, len; alpar@1: char *buf = mpl->tup_buf, str[255+1], *save; alpar@1: # define safe_append(c) \ alpar@1: (void)(len < 255 ? (buf[len++] = (char)(c)) : 0) alpar@1: buf[0] = '\0', len = 0; alpar@1: dim = tuple_dimen(mpl, tuple); alpar@1: if (c == '[' && dim > 0) safe_append('['); alpar@1: if (c == '(' && dim > 1) safe_append('('); alpar@1: for (temp = tuple; temp != NULL; temp = temp->next) alpar@1: { if (temp != tuple) safe_append(','); alpar@1: xassert(temp->sym != NULL); alpar@1: save = mpl->sym_buf; alpar@1: mpl->sym_buf = str; alpar@1: format_symbol(mpl, temp->sym); alpar@1: mpl->sym_buf = save; alpar@1: xassert(strlen(str) < sizeof(str)); alpar@1: for (j = 0; str[j] != '\0'; j++) safe_append(str[j]); alpar@1: } alpar@1: if (c == '[' && dim > 0) safe_append(']'); alpar@1: if (c == '(' && dim > 1) safe_append(')'); alpar@1: # undef safe_append alpar@1: buf[len] = '\0'; alpar@1: if (len == 255) strcpy(buf+252, "..."); alpar@1: xassert(strlen(buf) <= 255); alpar@1: return buf; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * ELEMENTAL SETS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- create_elemset - create elemental set. alpar@1: -- alpar@1: -- This routine creates an elemental set, whose members are n-tuples of alpar@1: -- specified dimension. Being created the set is initially empty. */ alpar@1: alpar@1: ELEMSET *create_elemset(MPL *mpl, int dim) alpar@1: { ELEMSET *set; alpar@1: xassert(dim > 0); alpar@1: set = create_array(mpl, A_NONE, dim); alpar@1: return set; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- find_tuple - check if elemental set contains given n-tuple. alpar@1: -- alpar@1: -- This routine finds given n-tuple in specified elemental set in order alpar@1: -- to check if the set contains that n-tuple. If the n-tuple is found, alpar@1: -- the routine returns pointer to corresponding array member. Otherwise alpar@1: -- null pointer is returned. */ alpar@1: alpar@1: MEMBER *find_tuple alpar@1: ( MPL *mpl, alpar@1: ELEMSET *set, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { xassert(set != NULL); alpar@1: xassert(set->type == A_NONE); alpar@1: xassert(set->dim == tuple_dimen(mpl, tuple)); alpar@1: return find_member(mpl, set, tuple); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- add_tuple - add new n-tuple to elemental set. alpar@1: -- alpar@1: -- This routine adds given n-tuple to specified elemental set. alpar@1: -- alpar@1: -- For the sake of efficiency this routine doesn't check whether the alpar@1: -- set already contains the same n-tuple or not. Therefore the calling alpar@1: -- program should use the routine find_tuple (if necessary) in order to alpar@1: -- make sure that the given n-tuple is not contained in the set, since alpar@1: -- duplicate n-tuples within the same set are not allowed. */ alpar@1: alpar@1: MEMBER *add_tuple alpar@1: ( MPL *mpl, alpar@1: ELEMSET *set, /* modified */ alpar@1: TUPLE *tuple /* destroyed */ alpar@1: ) alpar@1: { MEMBER *memb; alpar@1: xassert(set != NULL); alpar@1: xassert(set->type == A_NONE); alpar@1: xassert(set->dim == tuple_dimen(mpl, tuple)); alpar@1: memb = add_member(mpl, set, tuple); alpar@1: memb->value.none = NULL; alpar@1: return memb; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- check_then_add - check and add new n-tuple to elemental set. alpar@1: -- alpar@1: -- This routine is equivalent to the routine add_tuple except that it alpar@1: -- does check for duplicate n-tuples. */ alpar@1: alpar@1: MEMBER *check_then_add alpar@1: ( MPL *mpl, alpar@1: ELEMSET *set, /* modified */ alpar@1: TUPLE *tuple /* destroyed */ alpar@1: ) alpar@1: { if (find_tuple(mpl, set, tuple) != NULL) alpar@1: error(mpl, "duplicate tuple %s detected", format_tuple(mpl, alpar@1: '(', tuple)); alpar@1: return add_tuple(mpl, set, tuple); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- copy_elemset - make copy of elemental set. alpar@1: -- alpar@1: -- This routine makes an exact copy of elemental set. */ alpar@1: alpar@1: ELEMSET *copy_elemset alpar@1: ( MPL *mpl, alpar@1: ELEMSET *set /* not changed */ alpar@1: ) alpar@1: { ELEMSET *copy; alpar@1: MEMBER *memb; alpar@1: xassert(set != NULL); alpar@1: xassert(set->type == A_NONE); alpar@1: xassert(set->dim > 0); alpar@1: copy = create_elemset(mpl, set->dim); alpar@1: for (memb = set->head; memb != NULL; memb = memb->next) alpar@1: add_tuple(mpl, copy, copy_tuple(mpl, memb->tuple)); alpar@1: return copy; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- delete_elemset - delete elemental set. alpar@1: -- alpar@1: -- This routine deletes specified elemental set. */ alpar@1: alpar@1: void delete_elemset alpar@1: ( MPL *mpl, alpar@1: ELEMSET *set /* destroyed */ alpar@1: ) alpar@1: { xassert(set != NULL); alpar@1: xassert(set->type == A_NONE); alpar@1: delete_array(mpl, set); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- arelset_size - compute size of "arithmetic" elemental set. alpar@1: -- alpar@1: -- This routine computes the size of "arithmetic" elemental set, which alpar@1: -- is specified in the form of arithmetic progression: alpar@1: -- alpar@1: -- { t0 .. tf by dt }. alpar@1: -- alpar@1: -- The size is computed using the formula: alpar@1: -- alpar@1: -- n = max(0, floor((tf - t0) / dt) + 1). */ alpar@1: alpar@1: int arelset_size(MPL *mpl, double t0, double tf, double dt) alpar@1: { double temp; alpar@1: if (dt == 0.0) alpar@1: error(mpl, "%.*g .. %.*g by %.*g; zero stride not allowed", alpar@1: DBL_DIG, t0, DBL_DIG, tf, DBL_DIG, dt); alpar@1: if (tf > 0.0 && t0 < 0.0 && tf > + 0.999 * DBL_MAX + t0) alpar@1: temp = +DBL_MAX; alpar@1: else if (tf < 0.0 && t0 > 0.0 && tf < - 0.999 * DBL_MAX + t0) alpar@1: temp = -DBL_MAX; alpar@1: else alpar@1: temp = tf - t0; alpar@1: if (fabs(dt) < 1.0 && fabs(temp) > (0.999 * DBL_MAX) * fabs(dt)) alpar@1: { if (temp > 0.0 && dt > 0.0 || temp < 0.0 && dt < 0.0) alpar@1: temp = +DBL_MAX; alpar@1: else alpar@1: temp = 0.0; alpar@1: } alpar@1: else alpar@1: { temp = floor(temp / dt) + 1.0; alpar@1: if (temp < 0.0) temp = 0.0; alpar@1: } alpar@1: xassert(temp >= 0.0); alpar@1: if (temp > (double)(INT_MAX - 1)) alpar@1: error(mpl, "%.*g .. %.*g by %.*g; set too large", alpar@1: DBL_DIG, t0, DBL_DIG, tf, DBL_DIG, dt); alpar@1: return (int)(temp + 0.5); alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- arelset_member - compute member of "arithmetic" elemental set. alpar@1: -- alpar@1: -- This routine returns a numeric value of symbol, which is equivalent alpar@1: -- to j-th member of given "arithmetic" elemental set specified in the alpar@1: -- form of arithmetic progression: alpar@1: -- alpar@1: -- { t0 .. tf by dt }. alpar@1: -- alpar@1: -- The symbol value is computed with the formula: alpar@1: -- alpar@1: -- j-th member = t0 + (j - 1) * dt, alpar@1: -- alpar@1: -- The number j must satisfy to the restriction 1 <= j <= n, where n is alpar@1: -- the set size computed by the routine arelset_size. */ alpar@1: alpar@1: double arelset_member(MPL *mpl, double t0, double tf, double dt, int j) alpar@1: { xassert(1 <= j && j <= arelset_size(mpl, t0, tf, dt)); alpar@1: return t0 + (double)(j - 1) * dt; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- create_arelset - create "arithmetic" elemental set. alpar@1: -- alpar@1: -- This routine creates "arithmetic" elemental set, which is specified alpar@1: -- in the form of arithmetic progression: alpar@1: -- alpar@1: -- { t0 .. tf by dt }. alpar@1: -- alpar@1: -- Components of this set are 1-tuples. */ alpar@1: alpar@1: ELEMSET *create_arelset(MPL *mpl, double t0, double tf, double dt) alpar@1: { ELEMSET *set; alpar@1: int j, n; alpar@1: set = create_elemset(mpl, 1); alpar@1: n = arelset_size(mpl, t0, tf, dt); alpar@1: for (j = 1; j <= n; j++) alpar@1: { add_tuple alpar@1: ( mpl, alpar@1: set, alpar@1: expand_tuple alpar@1: ( mpl, alpar@1: create_tuple(mpl), alpar@1: create_symbol_num alpar@1: ( mpl, alpar@1: arelset_member(mpl, t0, tf, dt, j) alpar@1: ) alpar@1: ) alpar@1: ); alpar@1: } alpar@1: return set; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- set_union - union of two elemental sets. alpar@1: -- alpar@1: -- This routine computes the union: alpar@1: -- alpar@1: -- X U Y = { j | (j in X) or (j in Y) }, alpar@1: -- alpar@1: -- where X and Y are given elemental sets (destroyed on exit). */ alpar@1: alpar@1: ELEMSET *set_union alpar@1: ( MPL *mpl, alpar@1: ELEMSET *X, /* destroyed */ alpar@1: ELEMSET *Y /* destroyed */ alpar@1: ) alpar@1: { MEMBER *memb; alpar@1: xassert(X != NULL); alpar@1: xassert(X->type == A_NONE); alpar@1: xassert(X->dim > 0); alpar@1: xassert(Y != NULL); alpar@1: xassert(Y->type == A_NONE); alpar@1: xassert(Y->dim > 0); alpar@1: xassert(X->dim == Y->dim); alpar@1: for (memb = Y->head; memb != NULL; memb = memb->next) alpar@1: { if (find_tuple(mpl, X, memb->tuple) == NULL) alpar@1: add_tuple(mpl, X, copy_tuple(mpl, memb->tuple)); alpar@1: } alpar@1: delete_elemset(mpl, Y); alpar@1: return X; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- set_diff - difference between two elemental sets. alpar@1: -- alpar@1: -- This routine computes the difference: alpar@1: -- alpar@1: -- X \ Y = { j | (j in X) and (j not in Y) }, alpar@1: -- alpar@1: -- where X and Y are given elemental sets (destroyed on exit). */ alpar@1: alpar@1: ELEMSET *set_diff alpar@1: ( MPL *mpl, alpar@1: ELEMSET *X, /* destroyed */ alpar@1: ELEMSET *Y /* destroyed */ alpar@1: ) alpar@1: { ELEMSET *Z; alpar@1: MEMBER *memb; alpar@1: xassert(X != NULL); alpar@1: xassert(X->type == A_NONE); alpar@1: xassert(X->dim > 0); alpar@1: xassert(Y != NULL); alpar@1: xassert(Y->type == A_NONE); alpar@1: xassert(Y->dim > 0); alpar@1: xassert(X->dim == Y->dim); alpar@1: Z = create_elemset(mpl, X->dim); alpar@1: for (memb = X->head; memb != NULL; memb = memb->next) alpar@1: { if (find_tuple(mpl, Y, memb->tuple) == NULL) alpar@1: add_tuple(mpl, Z, copy_tuple(mpl, memb->tuple)); alpar@1: } alpar@1: delete_elemset(mpl, X); alpar@1: delete_elemset(mpl, Y); alpar@1: return Z; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- set_symdiff - symmetric difference between two elemental sets. alpar@1: -- alpar@1: -- This routine computes the symmetric difference: alpar@1: -- alpar@1: -- X (+) Y = (X \ Y) U (Y \ X), alpar@1: -- alpar@1: -- where X and Y are given elemental sets (destroyed on exit). */ alpar@1: alpar@1: ELEMSET *set_symdiff alpar@1: ( MPL *mpl, alpar@1: ELEMSET *X, /* destroyed */ alpar@1: ELEMSET *Y /* destroyed */ alpar@1: ) alpar@1: { ELEMSET *Z; alpar@1: MEMBER *memb; alpar@1: xassert(X != NULL); alpar@1: xassert(X->type == A_NONE); alpar@1: xassert(X->dim > 0); alpar@1: xassert(Y != NULL); alpar@1: xassert(Y->type == A_NONE); alpar@1: xassert(Y->dim > 0); alpar@1: xassert(X->dim == Y->dim); alpar@1: /* Z := X \ Y */ alpar@1: Z = create_elemset(mpl, X->dim); alpar@1: for (memb = X->head; memb != NULL; memb = memb->next) alpar@1: { if (find_tuple(mpl, Y, memb->tuple) == NULL) alpar@1: add_tuple(mpl, Z, copy_tuple(mpl, memb->tuple)); alpar@1: } alpar@1: /* Z := Z U (Y \ X) */ alpar@1: for (memb = Y->head; memb != NULL; memb = memb->next) alpar@1: { if (find_tuple(mpl, X, memb->tuple) == NULL) alpar@1: add_tuple(mpl, Z, copy_tuple(mpl, memb->tuple)); alpar@1: } alpar@1: delete_elemset(mpl, X); alpar@1: delete_elemset(mpl, Y); alpar@1: return Z; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- set_inter - intersection of two elemental sets. alpar@1: -- alpar@1: -- This routine computes the intersection: alpar@1: -- alpar@1: -- X ^ Y = { j | (j in X) and (j in Y) }, alpar@1: -- alpar@1: -- where X and Y are given elemental sets (destroyed on exit). */ alpar@1: alpar@1: ELEMSET *set_inter alpar@1: ( MPL *mpl, alpar@1: ELEMSET *X, /* destroyed */ alpar@1: ELEMSET *Y /* destroyed */ alpar@1: ) alpar@1: { ELEMSET *Z; alpar@1: MEMBER *memb; alpar@1: xassert(X != NULL); alpar@1: xassert(X->type == A_NONE); alpar@1: xassert(X->dim > 0); alpar@1: xassert(Y != NULL); alpar@1: xassert(Y->type == A_NONE); alpar@1: xassert(Y->dim > 0); alpar@1: xassert(X->dim == Y->dim); alpar@1: Z = create_elemset(mpl, X->dim); alpar@1: for (memb = X->head; memb != NULL; memb = memb->next) alpar@1: { if (find_tuple(mpl, Y, memb->tuple) != NULL) alpar@1: add_tuple(mpl, Z, copy_tuple(mpl, memb->tuple)); alpar@1: } alpar@1: delete_elemset(mpl, X); alpar@1: delete_elemset(mpl, Y); alpar@1: return Z; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- set_cross - cross (Cartesian) product of two elemental sets. alpar@1: -- alpar@1: -- This routine computes the cross (Cartesian) product: alpar@1: -- alpar@1: -- X x Y = { (i,j) | (i in X) and (j in Y) }, alpar@1: -- alpar@1: -- where X and Y are given elemental sets (destroyed on exit). */ alpar@1: alpar@1: ELEMSET *set_cross alpar@1: ( MPL *mpl, alpar@1: ELEMSET *X, /* destroyed */ alpar@1: ELEMSET *Y /* destroyed */ alpar@1: ) alpar@1: { ELEMSET *Z; alpar@1: MEMBER *memx, *memy; alpar@1: TUPLE *tuple, *temp; alpar@1: xassert(X != NULL); alpar@1: xassert(X->type == A_NONE); alpar@1: xassert(X->dim > 0); alpar@1: xassert(Y != NULL); alpar@1: xassert(Y->type == A_NONE); alpar@1: xassert(Y->dim > 0); alpar@1: Z = create_elemset(mpl, X->dim + Y->dim); alpar@1: for (memx = X->head; memx != NULL; memx = memx->next) alpar@1: { for (memy = Y->head; memy != NULL; memy = memy->next) alpar@1: { tuple = copy_tuple(mpl, memx->tuple); alpar@1: for (temp = memy->tuple; temp != NULL; temp = temp->next) alpar@1: tuple = expand_tuple(mpl, tuple, copy_symbol(mpl, alpar@1: temp->sym)); alpar@1: add_tuple(mpl, Z, tuple); alpar@1: } alpar@1: } alpar@1: delete_elemset(mpl, X); alpar@1: delete_elemset(mpl, Y); alpar@1: return Z; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * ELEMENTAL VARIABLES * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /* (there are no specific routines for elemental variables) */ alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * LINEAR FORMS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- constant_term - create constant term. alpar@1: -- alpar@1: -- This routine creates the linear form, which is a constant term. */ alpar@1: alpar@1: FORMULA *constant_term(MPL *mpl, double coef) alpar@1: { FORMULA *form; alpar@1: if (coef == 0.0) alpar@1: form = NULL; alpar@1: else alpar@1: { form = dmp_get_atom(mpl->formulae, sizeof(FORMULA)); alpar@1: form->coef = coef; alpar@1: form->var = NULL; alpar@1: form->next = NULL; alpar@1: } alpar@1: return form; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- single_variable - create single variable. alpar@1: -- alpar@1: -- This routine creates the linear form, which is a single elemental alpar@1: -- variable. */ alpar@1: alpar@1: FORMULA *single_variable alpar@1: ( MPL *mpl, alpar@1: ELEMVAR *var /* referenced */ alpar@1: ) alpar@1: { FORMULA *form; alpar@1: xassert(var != NULL); alpar@1: form = dmp_get_atom(mpl->formulae, sizeof(FORMULA)); alpar@1: form->coef = 1.0; alpar@1: form->var = var; alpar@1: form->next = NULL; alpar@1: return form; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- copy_formula - make copy of linear form. alpar@1: -- alpar@1: -- This routine returns an exact copy of linear form. */ alpar@1: alpar@1: FORMULA *copy_formula alpar@1: ( MPL *mpl, alpar@1: FORMULA *form /* not changed */ alpar@1: ) alpar@1: { FORMULA *head, *tail; alpar@1: if (form == NULL) alpar@1: head = NULL; alpar@1: else alpar@1: { head = tail = dmp_get_atom(mpl->formulae, sizeof(FORMULA)); alpar@1: for (; form != NULL; form = form->next) alpar@1: { tail->coef = form->coef; alpar@1: tail->var = form->var; alpar@1: if (form->next != NULL) alpar@1: tail = (tail->next = dmp_get_atom(mpl->formulae, sizeof(FORMULA))); alpar@1: } alpar@1: tail->next = NULL; alpar@1: } alpar@1: return head; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- delete_formula - delete linear form. alpar@1: -- alpar@1: -- This routine deletes specified linear form. */ alpar@1: alpar@1: void delete_formula alpar@1: ( MPL *mpl, alpar@1: FORMULA *form /* destroyed */ alpar@1: ) alpar@1: { FORMULA *temp; alpar@1: while (form != NULL) alpar@1: { temp = form; alpar@1: form = form->next; alpar@1: dmp_free_atom(mpl->formulae, temp, sizeof(FORMULA)); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- linear_comb - linear combination of two linear forms. alpar@1: -- alpar@1: -- This routine computes the linear combination: alpar@1: -- alpar@1: -- a * fx + b * fy, alpar@1: -- alpar@1: -- where a and b are numeric coefficients, fx and fy are linear forms alpar@1: -- (destroyed on exit). */ alpar@1: alpar@1: FORMULA *linear_comb alpar@1: ( MPL *mpl, alpar@1: double a, FORMULA *fx, /* destroyed */ alpar@1: double b, FORMULA *fy /* destroyed */ alpar@1: ) alpar@1: { FORMULA *form = NULL, *term, *temp; alpar@1: double c0 = 0.0; alpar@1: for (term = fx; term != NULL; term = term->next) alpar@1: { if (term->var == NULL) alpar@1: c0 = fp_add(mpl, c0, fp_mul(mpl, a, term->coef)); alpar@1: else alpar@1: term->var->temp = alpar@1: fp_add(mpl, term->var->temp, fp_mul(mpl, a, term->coef)); alpar@1: } alpar@1: for (term = fy; term != NULL; term = term->next) alpar@1: { if (term->var == NULL) alpar@1: c0 = fp_add(mpl, c0, fp_mul(mpl, b, term->coef)); alpar@1: else alpar@1: term->var->temp = alpar@1: fp_add(mpl, term->var->temp, fp_mul(mpl, b, term->coef)); alpar@1: } alpar@1: for (term = fx; term != NULL; term = term->next) alpar@1: { if (term->var != NULL && term->var->temp != 0.0) alpar@1: { temp = dmp_get_atom(mpl->formulae, sizeof(FORMULA)); alpar@1: temp->coef = term->var->temp, temp->var = term->var; alpar@1: temp->next = form, form = temp; alpar@1: term->var->temp = 0.0; alpar@1: } alpar@1: } alpar@1: for (term = fy; term != NULL; term = term->next) alpar@1: { if (term->var != NULL && term->var->temp != 0.0) alpar@1: { temp = dmp_get_atom(mpl->formulae, sizeof(FORMULA)); alpar@1: temp->coef = term->var->temp, temp->var = term->var; alpar@1: temp->next = form, form = temp; alpar@1: term->var->temp = 0.0; alpar@1: } alpar@1: } alpar@1: if (c0 != 0.0) alpar@1: { temp = dmp_get_atom(mpl->formulae, sizeof(FORMULA)); alpar@1: temp->coef = c0, temp->var = NULL; alpar@1: temp->next = form, form = temp; alpar@1: } alpar@1: delete_formula(mpl, fx); alpar@1: delete_formula(mpl, fy); alpar@1: return form; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- remove_constant - remove constant term from linear form. alpar@1: -- alpar@1: -- This routine removes constant term from linear form and stores its alpar@1: -- value to given location. */ alpar@1: alpar@1: FORMULA *remove_constant alpar@1: ( MPL *mpl, alpar@1: FORMULA *form, /* destroyed */ alpar@1: double *coef /* modified */ alpar@1: ) alpar@1: { FORMULA *head = NULL, *temp; alpar@1: *coef = 0.0; alpar@1: while (form != NULL) alpar@1: { temp = form; alpar@1: form = form->next; alpar@1: if (temp->var == NULL) alpar@1: { /* constant term */ alpar@1: *coef = fp_add(mpl, *coef, temp->coef); alpar@1: dmp_free_atom(mpl->formulae, temp, sizeof(FORMULA)); alpar@1: } alpar@1: else alpar@1: { /* linear term */ alpar@1: temp->next = head; alpar@1: head = temp; alpar@1: } alpar@1: } alpar@1: return head; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- reduce_terms - reduce identical terms in linear form. alpar@1: -- alpar@1: -- This routine reduces identical terms in specified linear form. */ alpar@1: alpar@1: FORMULA *reduce_terms alpar@1: ( MPL *mpl, alpar@1: FORMULA *form /* destroyed */ alpar@1: ) alpar@1: { FORMULA *term, *next_term; alpar@1: double c0 = 0.0; alpar@1: for (term = form; term != NULL; term = term->next) alpar@1: { if (term->var == NULL) alpar@1: c0 = fp_add(mpl, c0, term->coef); alpar@1: else alpar@1: term->var->temp = fp_add(mpl, term->var->temp, term->coef); alpar@1: } alpar@1: next_term = form, form = NULL; alpar@1: for (term = next_term; term != NULL; term = next_term) alpar@1: { next_term = term->next; alpar@1: if (term->var == NULL && c0 != 0.0) alpar@1: { term->coef = c0, c0 = 0.0; alpar@1: term->next = form, form = term; alpar@1: } alpar@1: else if (term->var != NULL && term->var->temp != 0.0) alpar@1: { term->coef = term->var->temp, term->var->temp = 0.0; alpar@1: term->next = form, form = term; alpar@1: } alpar@1: else alpar@1: dmp_free_atom(mpl->formulae, term, sizeof(FORMULA)); alpar@1: } alpar@1: return form; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * ELEMENTAL CONSTRAINTS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /* (there are no specific routines for elemental constraints) */ alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * GENERIC VALUES * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- delete_value - delete generic value. alpar@1: -- alpar@1: -- This routine deletes specified generic value. alpar@1: -- alpar@1: -- NOTE: The generic value to be deleted must be valid. */ alpar@1: alpar@1: void delete_value alpar@1: ( MPL *mpl, alpar@1: int type, alpar@1: VALUE *value /* content destroyed */ alpar@1: ) alpar@1: { xassert(value != NULL); alpar@1: switch (type) alpar@1: { case A_NONE: alpar@1: value->none = NULL; alpar@1: break; alpar@1: case A_NUMERIC: alpar@1: value->num = 0.0; alpar@1: break; alpar@1: case A_SYMBOLIC: alpar@1: delete_symbol(mpl, value->sym), value->sym = NULL; alpar@1: break; alpar@1: case A_LOGICAL: alpar@1: value->bit = 0; alpar@1: break; alpar@1: case A_TUPLE: alpar@1: delete_tuple(mpl, value->tuple), value->tuple = NULL; alpar@1: break; alpar@1: case A_ELEMSET: alpar@1: delete_elemset(mpl, value->set), value->set = NULL; alpar@1: break; alpar@1: case A_ELEMVAR: alpar@1: value->var = NULL; alpar@1: break; alpar@1: case A_FORMULA: alpar@1: delete_formula(mpl, value->form), value->form = NULL; alpar@1: break; alpar@1: case A_ELEMCON: alpar@1: value->con = NULL; alpar@1: break; alpar@1: default: alpar@1: xassert(type != type); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * SYMBOLICALLY INDEXED ARRAYS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- create_array - create array. alpar@1: -- alpar@1: -- This routine creates an array of specified type and dimension. Being alpar@1: -- created the array is initially empty. alpar@1: -- alpar@1: -- The type indicator determines generic values, which can be assigned alpar@1: -- to the array members: alpar@1: -- alpar@1: -- A_NONE - none (members have no assigned values) alpar@1: -- A_NUMERIC - floating-point numbers alpar@1: -- A_SYMBOLIC - symbols alpar@1: -- A_ELEMSET - elemental sets alpar@1: -- A_ELEMVAR - elemental variables alpar@1: -- A_ELEMCON - elemental constraints alpar@1: -- alpar@1: -- The dimension may be 0, in which case the array consists of the only alpar@1: -- member (such arrays represent 0-dimensional objects). */ alpar@1: alpar@1: ARRAY *create_array(MPL *mpl, int type, int dim) alpar@1: { ARRAY *array; alpar@1: xassert(type == A_NONE || type == A_NUMERIC || alpar@1: type == A_SYMBOLIC || type == A_ELEMSET || alpar@1: type == A_ELEMVAR || type == A_ELEMCON); alpar@1: xassert(dim >= 0); alpar@1: array = dmp_get_atom(mpl->arrays, sizeof(ARRAY)); alpar@1: array->type = type; alpar@1: array->dim = dim; alpar@1: array->size = 0; alpar@1: array->head = NULL; alpar@1: array->tail = NULL; alpar@1: array->tree = NULL; alpar@1: array->prev = NULL; alpar@1: array->next = mpl->a_list; alpar@1: /* include the array in the global array list */ alpar@1: if (array->next != NULL) array->next->prev = array; alpar@1: mpl->a_list = array; alpar@1: return array; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- find_member - find array member with given n-tuple. alpar@1: -- alpar@1: -- This routine finds an array member, which has given n-tuple. If the alpar@1: -- array is short, the linear search is used. Otherwise the routine alpar@1: -- autimatically creates the search tree (i.e. the array index) to find alpar@1: -- members for logarithmic time. */ alpar@1: alpar@1: static int compare_member_tuples(void *info, const void *key1, alpar@1: const void *key2) alpar@1: { /* this is an auxiliary routine used to compare keys, which are alpar@1: n-tuples assigned to array members */ alpar@1: return compare_tuples((MPL *)info, (TUPLE *)key1, (TUPLE *)key2); alpar@1: } alpar@1: alpar@1: MEMBER *find_member alpar@1: ( MPL *mpl, alpar@1: ARRAY *array, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { MEMBER *memb; alpar@1: xassert(array != NULL); alpar@1: /* the n-tuple must have the same dimension as the array */ alpar@1: xassert(tuple_dimen(mpl, tuple) == array->dim); alpar@1: /* if the array is large enough, create the search tree and index alpar@1: all existing members of the array */ alpar@1: if (array->size > 30 && array->tree == NULL) alpar@1: { array->tree = avl_create_tree(compare_member_tuples, mpl); alpar@1: for (memb = array->head; memb != NULL; memb = memb->next) alpar@1: avl_set_node_link(avl_insert_node(array->tree, memb->tuple), alpar@1: (void *)memb); alpar@1: } alpar@1: /* find a member, which has the given tuple */ alpar@1: if (array->tree == NULL) alpar@1: { /* the search tree doesn't exist; use the linear search */ alpar@1: for (memb = array->head; memb != NULL; memb = memb->next) alpar@1: if (compare_tuples(mpl, memb->tuple, tuple) == 0) break; alpar@1: } alpar@1: else alpar@1: { /* the search tree exists; use the binary search */ alpar@1: AVLNODE *node; alpar@1: node = avl_find_node(array->tree, tuple); alpar@1: memb = (MEMBER *)(node == NULL ? NULL : avl_get_node_link(node)); alpar@1: } alpar@1: return memb; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- add_member - add new member to array. alpar@1: -- alpar@1: -- This routine creates a new member with given n-tuple and adds it to alpar@1: -- specified array. alpar@1: -- alpar@1: -- For the sake of efficiency this routine doesn't check whether the alpar@1: -- array already contains a member with the given n-tuple or not. Thus, alpar@1: -- if necessary, the calling program should use the routine find_member alpar@1: -- in order to be sure that the array contains no member with the same alpar@1: -- n-tuple, because members with duplicate n-tuples are not allowed. alpar@1: -- alpar@1: -- This routine assigns no generic value to the new member, because the alpar@1: -- calling program must do that. */ alpar@1: alpar@1: MEMBER *add_member alpar@1: ( MPL *mpl, alpar@1: ARRAY *array, /* modified */ alpar@1: TUPLE *tuple /* destroyed */ alpar@1: ) alpar@1: { MEMBER *memb; alpar@1: xassert(array != NULL); alpar@1: /* the n-tuple must have the same dimension as the array */ alpar@1: xassert(tuple_dimen(mpl, tuple) == array->dim); alpar@1: /* create new member */ alpar@1: memb = dmp_get_atom(mpl->members, sizeof(MEMBER)); alpar@1: memb->tuple = tuple; alpar@1: memb->next = NULL; alpar@1: memset(&memb->value, '?', sizeof(VALUE)); alpar@1: /* and append it to the member list */ alpar@1: array->size++; alpar@1: if (array->head == NULL) alpar@1: array->head = memb; alpar@1: else alpar@1: array->tail->next = memb; alpar@1: array->tail = memb; alpar@1: /* if the search tree exists, index the new member */ alpar@1: if (array->tree != NULL) alpar@1: avl_set_node_link(avl_insert_node(array->tree, memb->tuple), alpar@1: (void *)memb); alpar@1: return memb; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- delete_array - delete array. alpar@1: -- alpar@1: -- This routine deletes specified array. alpar@1: -- alpar@1: -- Generic values assigned to the array members are not deleted by this alpar@1: -- routine. The calling program itself must delete all assigned generic alpar@1: -- values before deleting the array. */ alpar@1: alpar@1: void delete_array alpar@1: ( MPL *mpl, alpar@1: ARRAY *array /* destroyed */ alpar@1: ) alpar@1: { MEMBER *memb; alpar@1: xassert(array != NULL); alpar@1: /* delete all existing array members */ alpar@1: while (array->head != NULL) alpar@1: { memb = array->head; alpar@1: array->head = memb->next; alpar@1: delete_tuple(mpl, memb->tuple); alpar@1: dmp_free_atom(mpl->members, memb, sizeof(MEMBER)); alpar@1: } alpar@1: /* if the search tree exists, also delete it */ alpar@1: if (array->tree != NULL) avl_delete_tree(array->tree); alpar@1: /* remove the array from the global array list */ alpar@1: if (array->prev == NULL) alpar@1: mpl->a_list = array->next; alpar@1: else alpar@1: array->prev->next = array->next; alpar@1: if (array->next == NULL) alpar@1: ; alpar@1: else alpar@1: array->next->prev = array->prev; alpar@1: /* delete the array descriptor */ alpar@1: dmp_free_atom(mpl->arrays, array, sizeof(ARRAY)); alpar@1: return; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * DOMAINS AND DUMMY INDICES * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- assign_dummy_index - assign new value to dummy index. alpar@1: -- alpar@1: -- This routine assigns new value to specified dummy index and, that is alpar@1: -- important, invalidates all temporary resultant values, which depends alpar@1: -- on that dummy index. */ alpar@1: alpar@1: void assign_dummy_index alpar@1: ( MPL *mpl, alpar@1: DOMAIN_SLOT *slot, /* modified */ alpar@1: SYMBOL *value /* not changed */ alpar@1: ) alpar@1: { CODE *leaf, *code; alpar@1: xassert(slot != NULL); alpar@1: xassert(value != NULL); alpar@1: /* delete the current value assigned to the dummy index */ alpar@1: if (slot->value != NULL) alpar@1: { /* if the current value and the new one are identical, actual alpar@1: assignment is not needed */ alpar@1: if (compare_symbols(mpl, slot->value, value) == 0) goto done; alpar@1: /* delete a symbol, which is the current value */ alpar@1: delete_symbol(mpl, slot->value), slot->value = NULL; alpar@1: } alpar@1: /* now walk through all the pseudo-codes with op = O_INDEX, which alpar@1: refer to the dummy index to be changed (these pseudo-codes are alpar@1: leaves in the forest of *all* expressions in the database) */ alpar@1: for (leaf = slot->list; leaf != NULL; leaf = leaf->arg.index. alpar@1: next) alpar@1: { xassert(leaf->op == O_INDEX); alpar@1: /* invalidate all resultant values, which depend on the dummy alpar@1: index, walking from the current leaf toward the root of the alpar@1: corresponding expression tree */ alpar@1: for (code = leaf; code != NULL; code = code->up) alpar@1: { if (code->valid) alpar@1: { /* invalidate and delete resultant value */ alpar@1: code->valid = 0; alpar@1: delete_value(mpl, code->type, &code->value); alpar@1: } alpar@1: } alpar@1: } alpar@1: /* assign new value to the dummy index */ alpar@1: slot->value = copy_symbol(mpl, value); alpar@1: done: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- update_dummy_indices - update current values of dummy indices. alpar@1: -- alpar@1: -- This routine assigns components of "backup" n-tuple to dummy indices alpar@1: -- of specified domain block. If no "backup" n-tuple is defined for the alpar@1: -- domain block, values of the dummy indices remain untouched. */ alpar@1: alpar@1: void update_dummy_indices alpar@1: ( MPL *mpl, alpar@1: DOMAIN_BLOCK *block /* not changed */ alpar@1: ) alpar@1: { DOMAIN_SLOT *slot; alpar@1: TUPLE *temp; alpar@1: if (block->backup != NULL) alpar@1: { for (slot = block->list, temp = block->backup; slot != NULL; alpar@1: slot = slot->next, temp = temp->next) alpar@1: { xassert(temp != NULL); alpar@1: xassert(temp->sym != NULL); alpar@1: assign_dummy_index(mpl, slot, temp->sym); alpar@1: } alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- enter_domain_block - enter domain block. alpar@1: -- alpar@1: -- Let specified domain block have the form: alpar@1: -- alpar@1: -- { ..., (j1, j2, ..., jn) in J, ... } alpar@1: -- alpar@1: -- where j1, j2, ..., jn are dummy indices, J is a basic set. alpar@1: -- alpar@1: -- This routine does the following: alpar@1: -- alpar@1: -- 1. Checks if the given n-tuple is a member of the basic set J. Note alpar@1: -- that J being *out of the scope* of the domain block cannot depend alpar@1: -- on the dummy indices in the same and inner domain blocks, so it alpar@1: -- can be computed before the dummy indices are assigned new values. alpar@1: -- If this check fails, the routine returns with non-zero code. alpar@1: -- alpar@1: -- 2. Saves current values of the dummy indices j1, j2, ..., jn. alpar@1: -- alpar@1: -- 3. Assigns new values, which are components of the given n-tuple, to alpar@1: -- the dummy indices j1, j2, ..., jn. If dimension of the n-tuple is alpar@1: -- larger than n, its extra components n+1, n+2, ... are not used. alpar@1: -- alpar@1: -- 4. Calls the formal routine func which either enters the next domain alpar@1: -- block or evaluates some code within the domain scope. alpar@1: -- alpar@1: -- 5. Restores former values of the dummy indices j1, j2, ..., jn. alpar@1: -- alpar@1: -- Since current values assigned to the dummy indices on entry to this alpar@1: -- routine are restored on exit, the formal routine func is allowed to alpar@1: -- call this routine recursively. */ alpar@1: alpar@1: int enter_domain_block alpar@1: ( MPL *mpl, alpar@1: DOMAIN_BLOCK *block, /* not changed */ alpar@1: TUPLE *tuple, /* not changed */ alpar@1: void *info, void (*func)(MPL *mpl, void *info) alpar@1: ) alpar@1: { TUPLE *backup; alpar@1: int ret = 0; alpar@1: /* check if the given n-tuple is a member of the basic set */ alpar@1: xassert(block->code != NULL); alpar@1: if (!is_member(mpl, block->code, tuple)) alpar@1: { ret = 1; alpar@1: goto done; alpar@1: } alpar@1: /* save reference to "backup" n-tuple, which was used to assign alpar@1: current values of the dummy indices (it is sufficient to save alpar@1: reference, not value, because that n-tuple is defined in some alpar@1: outer level of recursion and therefore cannot be changed on alpar@1: this and deeper recursive calls) */ alpar@1: backup = block->backup; alpar@1: /* set up new "backup" n-tuple, which defines new values of the alpar@1: dummy indices */ alpar@1: block->backup = tuple; alpar@1: /* assign new values to the dummy indices */ alpar@1: update_dummy_indices(mpl, block); alpar@1: /* call the formal routine that does the rest part of the job */ alpar@1: func(mpl, info); alpar@1: /* restore reference to the former "backup" n-tuple */ alpar@1: block->backup = backup; alpar@1: /* restore former values of the dummy indices; note that if the alpar@1: domain block just escaped has no other active instances which alpar@1: may exist due to recursion (it is indicated by a null pointer alpar@1: to the former n-tuple), former values of the dummy indices are alpar@1: undefined; therefore in this case the routine keeps currently alpar@1: assigned values of the dummy indices that involves keeping all alpar@1: dependent temporary results and thereby, if this domain block alpar@1: is not used recursively, allows improving efficiency */ alpar@1: update_dummy_indices(mpl, block); alpar@1: done: return ret; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_within_domain - perform evaluation within domain scope. alpar@1: -- alpar@1: -- This routine assigns new values (symbols) to all dummy indices of alpar@1: -- specified domain and calls the formal routine func, which is used to alpar@1: -- evaluate some code in the domain scope. Each free dummy index in the alpar@1: -- domain is assigned a value specified in the corresponding component alpar@1: -- of given n-tuple. Non-free dummy indices are assigned values, which alpar@1: -- are computed by this routine. alpar@1: -- alpar@1: -- Number of components in the given n-tuple must be the same as number alpar@1: -- of free indices in the domain. alpar@1: -- alpar@1: -- If the given n-tuple is not a member of the domain set, the routine alpar@1: -- func is not called, and non-zero code is returned. alpar@1: -- alpar@1: -- For the sake of convenience it is allowed to specify domain as NULL alpar@1: -- (then n-tuple also must be 0-tuple, i.e. empty), in which case this alpar@1: -- routine just calls the routine func and returns zero. alpar@1: -- alpar@1: -- This routine allows recursive calls from the routine func providing alpar@1: -- correct values of dummy indices for each instance. alpar@1: -- alpar@1: -- NOTE: The n-tuple passed to this routine must not be changed by any alpar@1: -- other routines called from the formal routine func until this alpar@1: -- routine has returned. */ alpar@1: alpar@1: struct eval_domain_info alpar@1: { /* working info used by the routine eval_within_domain */ alpar@1: DOMAIN *domain; alpar@1: /* domain, which has to be entered */ alpar@1: DOMAIN_BLOCK *block; alpar@1: /* domain block, which is currently processed */ alpar@1: TUPLE *tuple; alpar@1: /* tail of original n-tuple, whose components have to be assigned alpar@1: to free dummy indices in the current domain block */ alpar@1: void *info; alpar@1: /* transit pointer passed to the formal routine func */ alpar@1: void (*func)(MPL *mpl, void *info); alpar@1: /* routine, which has to be executed in the domain scope */ alpar@1: int failure; alpar@1: /* this flag indicates that given n-tuple is not a member of the alpar@1: domain set */ alpar@1: }; alpar@1: alpar@1: static void eval_domain_func(MPL *mpl, void *_my_info) alpar@1: { /* this routine recursively enters into the domain scope and then alpar@1: calls the routine func */ alpar@1: struct eval_domain_info *my_info = _my_info; alpar@1: if (my_info->block != NULL) alpar@1: { /* the current domain block to be entered exists */ alpar@1: DOMAIN_BLOCK *block; alpar@1: DOMAIN_SLOT *slot; alpar@1: TUPLE *tuple = NULL, *temp = NULL; alpar@1: /* save pointer to the current domain block */ alpar@1: block = my_info->block; alpar@1: /* and get ready to enter the next block (if it exists) */ alpar@1: my_info->block = block->next; alpar@1: /* construct temporary n-tuple, whose components correspond to alpar@1: dummy indices (slots) of the current domain; components of alpar@1: the temporary n-tuple that correspond to free dummy indices alpar@1: are assigned references (not values!) to symbols specified alpar@1: in the corresponding components of the given n-tuple, while alpar@1: other components that correspond to non-free dummy indices alpar@1: are assigned symbolic values computed here */ alpar@1: for (slot = block->list; slot != NULL; slot = slot->next) alpar@1: { /* create component that corresponds to the current slot */ alpar@1: if (tuple == NULL) alpar@1: tuple = temp = dmp_get_atom(mpl->tuples, sizeof(TUPLE)); alpar@1: else alpar@1: temp = (temp->next = dmp_get_atom(mpl->tuples, sizeof(TUPLE))); alpar@1: if (slot->code == NULL) alpar@1: { /* dummy index is free; take reference to symbol, which alpar@1: is specified in the corresponding component of given alpar@1: n-tuple */ alpar@1: xassert(my_info->tuple != NULL); alpar@1: temp->sym = my_info->tuple->sym; alpar@1: xassert(temp->sym != NULL); alpar@1: my_info->tuple = my_info->tuple->next; alpar@1: } alpar@1: else alpar@1: { /* dummy index is non-free; compute symbolic value to be alpar@1: temporarily assigned to the dummy index */ alpar@1: temp->sym = eval_symbolic(mpl, slot->code); alpar@1: } alpar@1: } alpar@1: temp->next = NULL; alpar@1: /* enter the current domain block */ alpar@1: if (enter_domain_block(mpl, block, tuple, my_info, alpar@1: eval_domain_func)) my_info->failure = 1; alpar@1: /* delete temporary n-tuple as well as symbols that correspond alpar@1: to non-free dummy indices (they were computed here) */ alpar@1: for (slot = block->list; slot != NULL; slot = slot->next) alpar@1: { xassert(tuple != NULL); alpar@1: temp = tuple; alpar@1: tuple = tuple->next; alpar@1: if (slot->code != NULL) alpar@1: { /* dummy index is non-free; delete symbolic value */ alpar@1: delete_symbol(mpl, temp->sym); alpar@1: } alpar@1: /* delete component that corresponds to the current slot */ alpar@1: dmp_free_atom(mpl->tuples, temp, sizeof(TUPLE)); alpar@1: } alpar@1: } alpar@1: else alpar@1: { /* there are no more domain blocks, i.e. we have reached the alpar@1: domain scope */ alpar@1: xassert(my_info->tuple == NULL); alpar@1: /* check optional predicate specified for the domain */ alpar@1: if (my_info->domain->code != NULL && !eval_logical(mpl, alpar@1: my_info->domain->code)) alpar@1: { /* the predicate is false */ alpar@1: my_info->failure = 2; alpar@1: } alpar@1: else alpar@1: { /* the predicate is true; do the job */ alpar@1: my_info->func(mpl, my_info->info); alpar@1: } alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: int eval_within_domain alpar@1: ( MPL *mpl, alpar@1: DOMAIN *domain, /* not changed */ alpar@1: TUPLE *tuple, /* not changed */ alpar@1: void *info, void (*func)(MPL *mpl, void *info) alpar@1: ) alpar@1: { /* this routine performs evaluation within domain scope */ alpar@1: struct eval_domain_info _my_info, *my_info = &_my_info; alpar@1: if (domain == NULL) alpar@1: { xassert(tuple == NULL); alpar@1: func(mpl, info); alpar@1: my_info->failure = 0; alpar@1: } alpar@1: else alpar@1: { xassert(tuple != NULL); alpar@1: my_info->domain = domain; alpar@1: my_info->block = domain->list; alpar@1: my_info->tuple = tuple; alpar@1: my_info->info = info; alpar@1: my_info->func = func; alpar@1: my_info->failure = 0; alpar@1: /* enter the very first domain block */ alpar@1: eval_domain_func(mpl, my_info); alpar@1: } alpar@1: return my_info->failure; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- loop_within_domain - perform iterations within domain scope. alpar@1: -- alpar@1: -- This routine iteratively assigns new values (symbols) to the dummy alpar@1: -- indices of specified domain by enumerating all n-tuples, which are alpar@1: -- members of the domain set, and for every n-tuple it calls the formal alpar@1: -- routine func to evaluate some code within the domain scope. alpar@1: -- alpar@1: -- If the routine func returns non-zero, enumeration within the domain alpar@1: -- is prematurely terminated. alpar@1: -- alpar@1: -- For the sake of convenience it is allowed to specify domain as NULL, alpar@1: -- in which case this routine just calls the routine func only once and alpar@1: -- returns zero. alpar@1: -- alpar@1: -- This routine allows recursive calls from the routine func providing alpar@1: -- correct values of dummy indices for each instance. */ alpar@1: alpar@1: struct loop_domain_info alpar@1: { /* working info used by the routine loop_within_domain */ alpar@1: DOMAIN *domain; alpar@1: /* domain, which has to be entered */ alpar@1: DOMAIN_BLOCK *block; alpar@1: /* domain block, which is currently processed */ alpar@1: int looping; alpar@1: /* clearing this flag leads to terminating enumeration */ alpar@1: void *info; alpar@1: /* transit pointer passed to the formal routine func */ alpar@1: int (*func)(MPL *mpl, void *info); alpar@1: /* routine, which needs to be executed in the domain scope */ alpar@1: }; alpar@1: alpar@1: static void loop_domain_func(MPL *mpl, void *_my_info) alpar@1: { /* this routine enumerates all n-tuples in the basic set of the alpar@1: current domain block, enters recursively into the domain scope alpar@1: for every n-tuple, and then calls the routine func */ alpar@1: struct loop_domain_info *my_info = _my_info; alpar@1: if (my_info->block != NULL) alpar@1: { /* the current domain block to be entered exists */ alpar@1: DOMAIN_BLOCK *block; alpar@1: DOMAIN_SLOT *slot; alpar@1: TUPLE *bound; alpar@1: /* save pointer to the current domain block */ alpar@1: block = my_info->block; alpar@1: /* and get ready to enter the next block (if it exists) */ alpar@1: my_info->block = block->next; alpar@1: /* compute symbolic values, at which non-free dummy indices of alpar@1: the current domain block are bound; since that values don't alpar@1: depend on free dummy indices of the current block, they can alpar@1: be computed once out of the enumeration loop */ alpar@1: bound = create_tuple(mpl); alpar@1: for (slot = block->list; slot != NULL; slot = slot->next) alpar@1: { if (slot->code != NULL) alpar@1: bound = expand_tuple(mpl, bound, eval_symbolic(mpl, alpar@1: slot->code)); alpar@1: } alpar@1: /* start enumeration */ alpar@1: xassert(block->code != NULL); alpar@1: if (block->code->op == O_DOTS) alpar@1: { /* the basic set is "arithmetic", in which case it doesn't alpar@1: need to be computed explicitly */ alpar@1: TUPLE *tuple; alpar@1: int n, j; alpar@1: double t0, tf, dt; alpar@1: /* compute "parameters" of the basic set */ alpar@1: t0 = eval_numeric(mpl, block->code->arg.arg.x); alpar@1: tf = eval_numeric(mpl, block->code->arg.arg.y); alpar@1: if (block->code->arg.arg.z == NULL) alpar@1: dt = 1.0; alpar@1: else alpar@1: dt = eval_numeric(mpl, block->code->arg.arg.z); alpar@1: /* determine cardinality of the basic set */ alpar@1: n = arelset_size(mpl, t0, tf, dt); alpar@1: /* create dummy 1-tuple for members of the basic set */ alpar@1: tuple = expand_tuple(mpl, create_tuple(mpl), alpar@1: create_symbol_num(mpl, 0.0)); alpar@1: /* in case of "arithmetic" set there is exactly one dummy alpar@1: index, which cannot be non-free */ alpar@1: xassert(bound == NULL); alpar@1: /* walk through 1-tuples of the basic set */ alpar@1: for (j = 1; j <= n && my_info->looping; j++) alpar@1: { /* construct dummy 1-tuple for the current member */ alpar@1: tuple->sym->num = arelset_member(mpl, t0, tf, dt, j); alpar@1: /* enter the current domain block */ alpar@1: enter_domain_block(mpl, block, tuple, my_info, alpar@1: loop_domain_func); alpar@1: } alpar@1: /* delete dummy 1-tuple */ alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: else alpar@1: { /* the basic set is of general kind, in which case it needs alpar@1: to be explicitly computed */ alpar@1: ELEMSET *set; alpar@1: MEMBER *memb; alpar@1: TUPLE *temp1, *temp2; alpar@1: /* compute the basic set */ alpar@1: set = eval_elemset(mpl, block->code); alpar@1: /* walk through all n-tuples of the basic set */ alpar@1: for (memb = set->head; memb != NULL && my_info->looping; alpar@1: memb = memb->next) alpar@1: { /* all components of the current n-tuple that correspond alpar@1: to non-free dummy indices must be feasible; otherwise alpar@1: the n-tuple is not in the basic set */ alpar@1: temp1 = memb->tuple; alpar@1: temp2 = bound; alpar@1: for (slot = block->list; slot != NULL; slot = slot->next) alpar@1: { xassert(temp1 != NULL); alpar@1: if (slot->code != NULL) alpar@1: { /* non-free dummy index */ alpar@1: xassert(temp2 != NULL); alpar@1: if (compare_symbols(mpl, temp1->sym, temp2->sym) alpar@1: != 0) alpar@1: { /* the n-tuple is not in the basic set */ alpar@1: goto skip; alpar@1: } alpar@1: temp2 = temp2->next; alpar@1: } alpar@1: temp1 = temp1->next; alpar@1: } alpar@1: xassert(temp1 == NULL); alpar@1: xassert(temp2 == NULL); alpar@1: /* enter the current domain block */ alpar@1: enter_domain_block(mpl, block, memb->tuple, my_info, alpar@1: loop_domain_func); alpar@1: skip: ; alpar@1: } alpar@1: /* delete the basic set */ alpar@1: delete_elemset(mpl, set); alpar@1: } alpar@1: /* delete symbolic values binding non-free dummy indices */ alpar@1: delete_tuple(mpl, bound); alpar@1: /* restore pointer to the current domain block */ alpar@1: my_info->block = block; alpar@1: } alpar@1: else alpar@1: { /* there are no more domain blocks, i.e. we have reached the alpar@1: domain scope */ alpar@1: /* check optional predicate specified for the domain */ alpar@1: if (my_info->domain->code != NULL && !eval_logical(mpl, alpar@1: my_info->domain->code)) alpar@1: { /* the predicate is false */ alpar@1: /* nop */; alpar@1: } alpar@1: else alpar@1: { /* the predicate is true; do the job */ alpar@1: my_info->looping = !my_info->func(mpl, my_info->info); alpar@1: } alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: void loop_within_domain alpar@1: ( MPL *mpl, alpar@1: DOMAIN *domain, /* not changed */ alpar@1: void *info, int (*func)(MPL *mpl, void *info) alpar@1: ) alpar@1: { /* this routine performs iterations within domain scope */ alpar@1: struct loop_domain_info _my_info, *my_info = &_my_info; alpar@1: if (domain == NULL) alpar@1: func(mpl, info); alpar@1: else alpar@1: { my_info->domain = domain; alpar@1: my_info->block = domain->list; alpar@1: my_info->looping = 1; alpar@1: my_info->info = info; alpar@1: my_info->func = func; alpar@1: /* enter the very first domain block */ alpar@1: loop_domain_func(mpl, my_info); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- out_of_domain - raise domain exception. alpar@1: -- alpar@1: -- This routine is called when a reference is made to a member of some alpar@1: -- model object, but its n-tuple is out of the object domain. */ alpar@1: alpar@1: void out_of_domain alpar@1: ( MPL *mpl, alpar@1: char *name, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { xassert(name != NULL); alpar@1: xassert(tuple != NULL); alpar@1: error(mpl, "%s%s out of domain", name, format_tuple(mpl, '[', alpar@1: tuple)); alpar@1: /* no return */ alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- get_domain_tuple - obtain current n-tuple from domain. alpar@1: -- alpar@1: -- This routine constructs n-tuple, whose components are current values alpar@1: -- assigned to *free* dummy indices of specified domain. alpar@1: -- alpar@1: -- For the sake of convenience it is allowed to specify domain as NULL, alpar@1: -- in which case this routine returns 0-tuple. alpar@1: -- alpar@1: -- NOTE: This routine must not be called out of domain scope. */ alpar@1: alpar@1: TUPLE *get_domain_tuple alpar@1: ( MPL *mpl, alpar@1: DOMAIN *domain /* not changed */ alpar@1: ) alpar@1: { DOMAIN_BLOCK *block; alpar@1: DOMAIN_SLOT *slot; alpar@1: TUPLE *tuple; alpar@1: tuple = create_tuple(mpl); alpar@1: if (domain != NULL) alpar@1: { for (block = domain->list; block != NULL; block = block->next) alpar@1: { for (slot = block->list; slot != NULL; slot = slot->next) alpar@1: { if (slot->code == NULL) alpar@1: { xassert(slot->value != NULL); alpar@1: tuple = expand_tuple(mpl, tuple, copy_symbol(mpl, alpar@1: slot->value)); alpar@1: } alpar@1: } alpar@1: } alpar@1: } alpar@1: return tuple; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean_domain - clean domain. alpar@1: -- alpar@1: -- This routine cleans specified domain that assumes deleting all stuff alpar@1: -- dynamically allocated during the generation phase. */ alpar@1: alpar@1: void clean_domain(MPL *mpl, DOMAIN *domain) alpar@1: { DOMAIN_BLOCK *block; alpar@1: DOMAIN_SLOT *slot; alpar@1: /* if no domain is specified, do nothing */ alpar@1: if (domain == NULL) goto done; alpar@1: /* clean all domain blocks */ alpar@1: for (block = domain->list; block != NULL; block = block->next) alpar@1: { /* clean all domain slots */ alpar@1: for (slot = block->list; slot != NULL; slot = slot->next) alpar@1: { /* clean pseudo-code for computing bound value */ alpar@1: clean_code(mpl, slot->code); alpar@1: /* delete symbolic value assigned to dummy index */ alpar@1: if (slot->value != NULL) alpar@1: delete_symbol(mpl, slot->value), slot->value = NULL; alpar@1: } alpar@1: /* clean pseudo-code for computing basic set */ alpar@1: clean_code(mpl, block->code); alpar@1: } alpar@1: /* clean pseudo-code for computing domain predicate */ alpar@1: clean_code(mpl, domain->code); alpar@1: done: return; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * MODEL SETS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- check_elem_set - check elemental set assigned to set member. alpar@1: -- alpar@1: -- This routine checks if given elemental set being assigned to member alpar@1: -- of specified model set satisfies to all restrictions. alpar@1: -- alpar@1: -- NOTE: This routine must not be called out of domain scope. */ alpar@1: alpar@1: void check_elem_set alpar@1: ( MPL *mpl, alpar@1: SET *set, /* not changed */ alpar@1: TUPLE *tuple, /* not changed */ alpar@1: ELEMSET *refer /* not changed */ alpar@1: ) alpar@1: { WITHIN *within; alpar@1: MEMBER *memb; alpar@1: int eqno; alpar@1: /* elemental set must be within all specified supersets */ alpar@1: for (within = set->within, eqno = 1; within != NULL; within = alpar@1: within->next, eqno++) alpar@1: { xassert(within->code != NULL); alpar@1: for (memb = refer->head; memb != NULL; memb = memb->next) alpar@1: { if (!is_member(mpl, within->code, memb->tuple)) alpar@1: { char buf[255+1]; alpar@1: strcpy(buf, format_tuple(mpl, '(', memb->tuple)); alpar@1: xassert(strlen(buf) < sizeof(buf)); alpar@1: error(mpl, "%s%s contains %s which not within specified " alpar@1: "set; see (%d)", set->name, format_tuple(mpl, '[', alpar@1: tuple), buf, eqno); alpar@1: } alpar@1: } alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- take_member_set - obtain elemental set assigned to set member. alpar@1: -- alpar@1: -- This routine obtains a reference to elemental set assigned to given alpar@1: -- member of specified model set and returns it on exit. alpar@1: -- alpar@1: -- NOTE: This routine must not be called out of domain scope. */ alpar@1: alpar@1: ELEMSET *take_member_set /* returns reference, not value */ alpar@1: ( MPL *mpl, alpar@1: SET *set, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { MEMBER *memb; alpar@1: ELEMSET *refer; alpar@1: /* find member in the set array */ alpar@1: memb = find_member(mpl, set->array, tuple); alpar@1: if (memb != NULL) alpar@1: { /* member exists, so just take the reference */ alpar@1: refer = memb->value.set; alpar@1: } alpar@1: else if (set->assign != NULL) alpar@1: { /* compute value using assignment expression */ alpar@1: refer = eval_elemset(mpl, set->assign); alpar@1: add: /* check that the elemental set satisfies to all restrictions, alpar@1: assign it to new member, and add the member to the array */ alpar@1: check_elem_set(mpl, set, tuple, refer); alpar@1: memb = add_member(mpl, set->array, copy_tuple(mpl, tuple)); alpar@1: memb->value.set = refer; alpar@1: } alpar@1: else if (set->option != NULL) alpar@1: { /* compute default elemental set */ alpar@1: refer = eval_elemset(mpl, set->option); alpar@1: goto add; alpar@1: } alpar@1: else alpar@1: { /* no value (elemental set) is provided */ alpar@1: error(mpl, "no value for %s%s", set->name, format_tuple(mpl, alpar@1: '[', tuple)); alpar@1: } alpar@1: return refer; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_member_set - evaluate elemental set assigned to set member. alpar@1: -- alpar@1: -- This routine evaluates a reference to elemental set assigned to given alpar@1: -- member of specified model set and returns it on exit. */ alpar@1: alpar@1: struct eval_set_info alpar@1: { /* working info used by the routine eval_member_set */ alpar@1: SET *set; alpar@1: /* model set */ alpar@1: TUPLE *tuple; alpar@1: /* n-tuple, which defines set member */ alpar@1: MEMBER *memb; alpar@1: /* normally this pointer is NULL; the routine uses this pointer alpar@1: to check data provided in the data section, in which case it alpar@1: points to a member currently checked; this check is performed alpar@1: automatically only once when a reference to any member occurs alpar@1: for the first time */ alpar@1: ELEMSET *refer; alpar@1: /* evaluated reference to elemental set */ alpar@1: }; alpar@1: alpar@1: static void eval_set_func(MPL *mpl, void *_info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: struct eval_set_info *info = _info; alpar@1: if (info->memb != NULL) alpar@1: { /* checking call; check elemental set being assigned */ alpar@1: check_elem_set(mpl, info->set, info->memb->tuple, alpar@1: info->memb->value.set); alpar@1: } alpar@1: else alpar@1: { /* normal call; evaluate member, which has given n-tuple */ alpar@1: info->refer = take_member_set(mpl, info->set, info->tuple); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: #if 1 /* 12/XII-2008 */ alpar@1: static void saturate_set(MPL *mpl, SET *set) alpar@1: { GADGET *gadget = set->gadget; alpar@1: ELEMSET *data; alpar@1: MEMBER *elem, *memb; alpar@1: TUPLE *tuple, *work[20]; alpar@1: int i; alpar@1: xprintf("Generating %s...\n", set->name); alpar@1: eval_whole_set(mpl, gadget->set); alpar@1: /* gadget set must have exactly one member */ alpar@1: xassert(gadget->set->array != NULL); alpar@1: xassert(gadget->set->array->head != NULL); alpar@1: xassert(gadget->set->array->head == gadget->set->array->tail); alpar@1: data = gadget->set->array->head->value.set; alpar@1: xassert(data->type == A_NONE); alpar@1: xassert(data->dim == gadget->set->dimen); alpar@1: /* walk thru all elements of the plain set */ alpar@1: for (elem = data->head; elem != NULL; elem = elem->next) alpar@1: { /* create a copy of n-tuple */ alpar@1: tuple = copy_tuple(mpl, elem->tuple); alpar@1: /* rearrange component of the n-tuple */ alpar@1: for (i = 0; i < gadget->set->dimen; i++) alpar@1: work[i] = NULL; alpar@1: for (i = 0; tuple != NULL; tuple = tuple->next) alpar@1: work[gadget->ind[i++]-1] = tuple; alpar@1: xassert(i == gadget->set->dimen); alpar@1: for (i = 0; i < gadget->set->dimen; i++) alpar@1: { xassert(work[i] != NULL); alpar@1: work[i]->next = work[i+1]; alpar@1: } alpar@1: /* construct subscript list from first set->dim components */ alpar@1: if (set->dim == 0) alpar@1: tuple = NULL; alpar@1: else alpar@1: tuple = work[0], work[set->dim-1]->next = NULL; alpar@1: /* find corresponding member of the set to be initialized */ alpar@1: memb = find_member(mpl, set->array, tuple); alpar@1: if (memb == NULL) alpar@1: { /* not found; add new member to the set and assign it empty alpar@1: elemental set */ alpar@1: memb = add_member(mpl, set->array, tuple); alpar@1: memb->value.set = create_elemset(mpl, set->dimen); alpar@1: } alpar@1: else alpar@1: { /* found; free subscript list */ alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: /* construct new n-tuple from rest set->dimen components */ alpar@1: tuple = work[set->dim]; alpar@1: xassert(set->dim + set->dimen == gadget->set->dimen); alpar@1: work[gadget->set->dimen-1]->next = NULL; alpar@1: /* and add it to the elemental set assigned to the member alpar@1: (no check for duplicates is needed) */ alpar@1: add_tuple(mpl, memb->value.set, tuple); alpar@1: } alpar@1: /* the set has been saturated with data */ alpar@1: set->data = 1; alpar@1: return; alpar@1: } alpar@1: #endif alpar@1: alpar@1: ELEMSET *eval_member_set /* returns reference, not value */ alpar@1: ( MPL *mpl, alpar@1: SET *set, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { /* this routine evaluates set member */ alpar@1: struct eval_set_info _info, *info = &_info; alpar@1: xassert(set->dim == tuple_dimen(mpl, tuple)); alpar@1: info->set = set; alpar@1: info->tuple = tuple; alpar@1: #if 1 /* 12/XII-2008 */ alpar@1: if (set->gadget != NULL && set->data == 0) alpar@1: { /* initialize the set with data from a plain set */ alpar@1: saturate_set(mpl, set); alpar@1: } alpar@1: #endif alpar@1: if (set->data == 1) alpar@1: { /* check data, which are provided in the data section, but not alpar@1: checked yet */ alpar@1: /* save pointer to the last array member; note that during the alpar@1: check new members may be added beyond the last member due to alpar@1: references to the same parameter from default expression as alpar@1: well as from expressions that define restricting supersets; alpar@1: however, values assigned to the new members will be checked alpar@1: by other routine, so we don't need to check them here */ alpar@1: MEMBER *tail = set->array->tail; alpar@1: /* change the data status to prevent infinite recursive loop alpar@1: due to references to the same set during the check */ alpar@1: set->data = 2; alpar@1: /* check elemental sets assigned to array members in the data alpar@1: section until the marked member has been reached */ alpar@1: for (info->memb = set->array->head; info->memb != NULL; alpar@1: info->memb = info->memb->next) alpar@1: { if (eval_within_domain(mpl, set->domain, info->memb->tuple, alpar@1: info, eval_set_func)) alpar@1: out_of_domain(mpl, set->name, info->memb->tuple); alpar@1: if (info->memb == tail) break; alpar@1: } alpar@1: /* the check has been finished */ alpar@1: } alpar@1: /* evaluate member, which has given n-tuple */ alpar@1: info->memb = NULL; alpar@1: if (eval_within_domain(mpl, info->set->domain, info->tuple, info, alpar@1: eval_set_func)) alpar@1: out_of_domain(mpl, set->name, info->tuple); alpar@1: /* bring evaluated reference to the calling program */ alpar@1: return info->refer; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_whole_set - evaluate model set over entire domain. alpar@1: -- alpar@1: -- This routine evaluates all members of specified model set over entire alpar@1: -- domain. */ alpar@1: alpar@1: static int whole_set_func(MPL *mpl, void *info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: SET *set = (SET *)info; alpar@1: TUPLE *tuple = get_domain_tuple(mpl, set->domain); alpar@1: eval_member_set(mpl, set, tuple); alpar@1: delete_tuple(mpl, tuple); alpar@1: return 0; alpar@1: } alpar@1: alpar@1: void eval_whole_set(MPL *mpl, SET *set) alpar@1: { loop_within_domain(mpl, set->domain, set, whole_set_func); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean set - clean model set. alpar@1: -- alpar@1: -- This routine cleans specified model set that assumes deleting all alpar@1: -- stuff dynamically allocated during the generation phase. */ alpar@1: alpar@1: void clean_set(MPL *mpl, SET *set) alpar@1: { WITHIN *within; alpar@1: MEMBER *memb; alpar@1: /* clean subscript domain */ alpar@1: clean_domain(mpl, set->domain); alpar@1: /* clean pseudo-code for computing supersets */ alpar@1: for (within = set->within; within != NULL; within = within->next) alpar@1: clean_code(mpl, within->code); alpar@1: /* clean pseudo-code for computing assigned value */ alpar@1: clean_code(mpl, set->assign); alpar@1: /* clean pseudo-code for computing default value */ alpar@1: clean_code(mpl, set->option); alpar@1: /* reset data status flag */ alpar@1: set->data = 0; alpar@1: /* delete content array */ alpar@1: for (memb = set->array->head; memb != NULL; memb = memb->next) alpar@1: delete_value(mpl, set->array->type, &memb->value); alpar@1: delete_array(mpl, set->array), set->array = NULL; alpar@1: return; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * MODEL PARAMETERS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- check_value_num - check numeric value assigned to parameter member. alpar@1: -- alpar@1: -- This routine checks if numeric value being assigned to some member alpar@1: -- of specified numeric model parameter satisfies to all restrictions. alpar@1: -- alpar@1: -- NOTE: This routine must not be called out of domain scope. */ alpar@1: alpar@1: void check_value_num alpar@1: ( MPL *mpl, alpar@1: PARAMETER *par, /* not changed */ alpar@1: TUPLE *tuple, /* not changed */ alpar@1: double value alpar@1: ) alpar@1: { CONDITION *cond; alpar@1: WITHIN *in; alpar@1: int eqno; alpar@1: /* the value must satisfy to the parameter type */ alpar@1: switch (par->type) alpar@1: { case A_NUMERIC: alpar@1: break; alpar@1: case A_INTEGER: alpar@1: if (value != floor(value)) alpar@1: error(mpl, "%s%s = %.*g not integer", par->name, alpar@1: format_tuple(mpl, '[', tuple), DBL_DIG, value); alpar@1: break; alpar@1: case A_BINARY: alpar@1: if (!(value == 0.0 || value == 1.0)) alpar@1: error(mpl, "%s%s = %.*g not binary", par->name, alpar@1: format_tuple(mpl, '[', tuple), DBL_DIG, value); alpar@1: break; alpar@1: default: alpar@1: xassert(par != par); alpar@1: } alpar@1: /* the value must satisfy to all specified conditions */ alpar@1: for (cond = par->cond, eqno = 1; cond != NULL; cond = cond->next, alpar@1: eqno++) alpar@1: { double bound; alpar@1: char *rho; alpar@1: xassert(cond->code != NULL); alpar@1: bound = eval_numeric(mpl, cond->code); alpar@1: switch (cond->rho) alpar@1: { case O_LT: alpar@1: if (!(value < bound)) alpar@1: { rho = "<"; alpar@1: err: error(mpl, "%s%s = %.*g not %s %.*g; see (%d)", alpar@1: par->name, format_tuple(mpl, '[', tuple), DBL_DIG, alpar@1: value, rho, DBL_DIG, bound, eqno); alpar@1: } alpar@1: break; alpar@1: case O_LE: alpar@1: if (!(value <= bound)) { rho = "<="; goto err; } alpar@1: break; alpar@1: case O_EQ: alpar@1: if (!(value == bound)) { rho = "="; goto err; } alpar@1: break; alpar@1: case O_GE: alpar@1: if (!(value >= bound)) { rho = ">="; goto err; } alpar@1: break; alpar@1: case O_GT: alpar@1: if (!(value > bound)) { rho = ">"; goto err; } alpar@1: break; alpar@1: case O_NE: alpar@1: if (!(value != bound)) { rho = "<>"; goto err; } alpar@1: break; alpar@1: default: alpar@1: xassert(cond != cond); alpar@1: } alpar@1: } alpar@1: /* the value must be in all specified supersets */ alpar@1: for (in = par->in, eqno = 1; in != NULL; in = in->next, eqno++) alpar@1: { TUPLE *dummy; alpar@1: xassert(in->code != NULL); alpar@1: xassert(in->code->dim == 1); alpar@1: dummy = expand_tuple(mpl, create_tuple(mpl), alpar@1: create_symbol_num(mpl, value)); alpar@1: if (!is_member(mpl, in->code, dummy)) alpar@1: error(mpl, "%s%s = %.*g not in specified set; see (%d)", alpar@1: par->name, format_tuple(mpl, '[', tuple), DBL_DIG, alpar@1: value, eqno); alpar@1: delete_tuple(mpl, dummy); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- take_member_num - obtain num. value assigned to parameter member. alpar@1: -- alpar@1: -- This routine obtains a numeric value assigned to member of specified alpar@1: -- numeric model parameter and returns it on exit. alpar@1: -- alpar@1: -- NOTE: This routine must not be called out of domain scope. */ alpar@1: alpar@1: double take_member_num alpar@1: ( MPL *mpl, alpar@1: PARAMETER *par, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { MEMBER *memb; alpar@1: double value; alpar@1: /* find member in the parameter array */ alpar@1: memb = find_member(mpl, par->array, tuple); alpar@1: if (memb != NULL) alpar@1: { /* member exists, so just take its value */ alpar@1: value = memb->value.num; alpar@1: } alpar@1: else if (par->assign != NULL) alpar@1: { /* compute value using assignment expression */ alpar@1: value = eval_numeric(mpl, par->assign); alpar@1: add: /* check that the value satisfies to all restrictions, assign alpar@1: it to new member, and add the member to the array */ alpar@1: check_value_num(mpl, par, tuple, value); alpar@1: memb = add_member(mpl, par->array, copy_tuple(mpl, tuple)); alpar@1: memb->value.num = value; alpar@1: } alpar@1: else if (par->option != NULL) alpar@1: { /* compute default value */ alpar@1: value = eval_numeric(mpl, par->option); alpar@1: goto add; alpar@1: } alpar@1: else if (par->defval != NULL) alpar@1: { /* take default value provided in the data section */ alpar@1: if (par->defval->str != NULL) alpar@1: error(mpl, "cannot convert %s to floating-point number", alpar@1: format_symbol(mpl, par->defval)); alpar@1: value = par->defval->num; alpar@1: goto add; alpar@1: } alpar@1: else alpar@1: { /* no value is provided */ alpar@1: error(mpl, "no value for %s%s", par->name, format_tuple(mpl, alpar@1: '[', tuple)); alpar@1: } alpar@1: return value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_member_num - evaluate num. value assigned to parameter member. alpar@1: -- alpar@1: -- This routine evaluates a numeric value assigned to given member of alpar@1: -- specified numeric model parameter and returns it on exit. */ alpar@1: alpar@1: struct eval_num_info alpar@1: { /* working info used by the routine eval_member_num */ alpar@1: PARAMETER *par; alpar@1: /* model parameter */ alpar@1: TUPLE *tuple; alpar@1: /* n-tuple, which defines parameter member */ alpar@1: MEMBER *memb; alpar@1: /* normally this pointer is NULL; the routine uses this pointer alpar@1: to check data provided in the data section, in which case it alpar@1: points to a member currently checked; this check is performed alpar@1: automatically only once when a reference to any member occurs alpar@1: for the first time */ alpar@1: double value; alpar@1: /* evaluated numeric value */ alpar@1: }; alpar@1: alpar@1: static void eval_num_func(MPL *mpl, void *_info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: struct eval_num_info *info = _info; alpar@1: if (info->memb != NULL) alpar@1: { /* checking call; check numeric value being assigned */ alpar@1: check_value_num(mpl, info->par, info->memb->tuple, alpar@1: info->memb->value.num); alpar@1: } alpar@1: else alpar@1: { /* normal call; evaluate member, which has given n-tuple */ alpar@1: info->value = take_member_num(mpl, info->par, info->tuple); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: double eval_member_num alpar@1: ( MPL *mpl, alpar@1: PARAMETER *par, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { /* this routine evaluates numeric parameter member */ alpar@1: struct eval_num_info _info, *info = &_info; alpar@1: xassert(par->type == A_NUMERIC || par->type == A_INTEGER || alpar@1: par->type == A_BINARY); alpar@1: xassert(par->dim == tuple_dimen(mpl, tuple)); alpar@1: info->par = par; alpar@1: info->tuple = tuple; alpar@1: if (par->data == 1) alpar@1: { /* check data, which are provided in the data section, but not alpar@1: checked yet */ alpar@1: /* save pointer to the last array member; note that during the alpar@1: check new members may be added beyond the last member due to alpar@1: references to the same parameter from default expression as alpar@1: well as from expressions that define restricting conditions; alpar@1: however, values assigned to the new members will be checked alpar@1: by other routine, so we don't need to check them here */ alpar@1: MEMBER *tail = par->array->tail; alpar@1: /* change the data status to prevent infinite recursive loop alpar@1: due to references to the same parameter during the check */ alpar@1: par->data = 2; alpar@1: /* check values assigned to array members in the data section alpar@1: until the marked member has been reached */ alpar@1: for (info->memb = par->array->head; info->memb != NULL; alpar@1: info->memb = info->memb->next) alpar@1: { if (eval_within_domain(mpl, par->domain, info->memb->tuple, alpar@1: info, eval_num_func)) alpar@1: out_of_domain(mpl, par->name, info->memb->tuple); alpar@1: if (info->memb == tail) break; alpar@1: } alpar@1: /* the check has been finished */ alpar@1: } alpar@1: /* evaluate member, which has given n-tuple */ alpar@1: info->memb = NULL; alpar@1: if (eval_within_domain(mpl, info->par->domain, info->tuple, info, alpar@1: eval_num_func)) alpar@1: out_of_domain(mpl, par->name, info->tuple); alpar@1: /* bring evaluated value to the calling program */ alpar@1: return info->value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- check_value_sym - check symbolic value assigned to parameter member. alpar@1: -- alpar@1: -- This routine checks if symbolic value being assigned to some member alpar@1: -- of specified symbolic model parameter satisfies to all restrictions. alpar@1: -- alpar@1: -- NOTE: This routine must not be called out of domain scope. */ alpar@1: alpar@1: void check_value_sym alpar@1: ( MPL *mpl, alpar@1: PARAMETER *par, /* not changed */ alpar@1: TUPLE *tuple, /* not changed */ alpar@1: SYMBOL *value /* not changed */ alpar@1: ) alpar@1: { CONDITION *cond; alpar@1: WITHIN *in; alpar@1: int eqno; alpar@1: /* the value must satisfy to all specified conditions */ alpar@1: for (cond = par->cond, eqno = 1; cond != NULL; cond = cond->next, alpar@1: eqno++) alpar@1: { SYMBOL *bound; alpar@1: char buf[255+1]; alpar@1: xassert(cond->code != NULL); alpar@1: bound = eval_symbolic(mpl, cond->code); alpar@1: switch (cond->rho) alpar@1: { alpar@1: #if 1 /* 13/VIII-2008 */ alpar@1: case O_LT: alpar@1: if (!(compare_symbols(mpl, value, bound) < 0)) alpar@1: { strcpy(buf, format_symbol(mpl, bound)); alpar@1: xassert(strlen(buf) < sizeof(buf)); alpar@1: error(mpl, "%s%s = %s not < %s", alpar@1: par->name, format_tuple(mpl, '[', tuple), alpar@1: format_symbol(mpl, value), buf, eqno); alpar@1: } alpar@1: break; alpar@1: case O_LE: alpar@1: if (!(compare_symbols(mpl, value, bound) <= 0)) alpar@1: { strcpy(buf, format_symbol(mpl, bound)); alpar@1: xassert(strlen(buf) < sizeof(buf)); alpar@1: error(mpl, "%s%s = %s not <= %s", alpar@1: par->name, format_tuple(mpl, '[', tuple), alpar@1: format_symbol(mpl, value), buf, eqno); alpar@1: } alpar@1: break; alpar@1: #endif alpar@1: case O_EQ: alpar@1: if (!(compare_symbols(mpl, value, bound) == 0)) alpar@1: { strcpy(buf, format_symbol(mpl, bound)); alpar@1: xassert(strlen(buf) < sizeof(buf)); alpar@1: error(mpl, "%s%s = %s not = %s", alpar@1: par->name, format_tuple(mpl, '[', tuple), alpar@1: format_symbol(mpl, value), buf, eqno); alpar@1: } alpar@1: break; alpar@1: #if 1 /* 13/VIII-2008 */ alpar@1: case O_GE: alpar@1: if (!(compare_symbols(mpl, value, bound) >= 0)) alpar@1: { strcpy(buf, format_symbol(mpl, bound)); alpar@1: xassert(strlen(buf) < sizeof(buf)); alpar@1: error(mpl, "%s%s = %s not >= %s", alpar@1: par->name, format_tuple(mpl, '[', tuple), alpar@1: format_symbol(mpl, value), buf, eqno); alpar@1: } alpar@1: break; alpar@1: case O_GT: alpar@1: if (!(compare_symbols(mpl, value, bound) > 0)) alpar@1: { strcpy(buf, format_symbol(mpl, bound)); alpar@1: xassert(strlen(buf) < sizeof(buf)); alpar@1: error(mpl, "%s%s = %s not > %s", alpar@1: par->name, format_tuple(mpl, '[', tuple), alpar@1: format_symbol(mpl, value), buf, eqno); alpar@1: } alpar@1: break; alpar@1: #endif alpar@1: case O_NE: alpar@1: if (!(compare_symbols(mpl, value, bound) != 0)) alpar@1: { strcpy(buf, format_symbol(mpl, bound)); alpar@1: xassert(strlen(buf) < sizeof(buf)); alpar@1: error(mpl, "%s%s = %s not <> %s", alpar@1: par->name, format_tuple(mpl, '[', tuple), alpar@1: format_symbol(mpl, value), buf, eqno); alpar@1: } alpar@1: break; alpar@1: default: alpar@1: xassert(cond != cond); alpar@1: } alpar@1: delete_symbol(mpl, bound); alpar@1: } alpar@1: /* the value must be in all specified supersets */ alpar@1: for (in = par->in, eqno = 1; in != NULL; in = in->next, eqno++) alpar@1: { TUPLE *dummy; alpar@1: xassert(in->code != NULL); alpar@1: xassert(in->code->dim == 1); alpar@1: dummy = expand_tuple(mpl, create_tuple(mpl), copy_symbol(mpl, alpar@1: value)); alpar@1: if (!is_member(mpl, in->code, dummy)) alpar@1: error(mpl, "%s%s = %s not in specified set; see (%d)", alpar@1: par->name, format_tuple(mpl, '[', tuple), alpar@1: format_symbol(mpl, value), eqno); alpar@1: delete_tuple(mpl, dummy); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- take_member_sym - obtain symb. value assigned to parameter member. alpar@1: -- alpar@1: -- This routine obtains a symbolic value assigned to member of specified alpar@1: -- symbolic model parameter and returns it on exit. alpar@1: -- alpar@1: -- NOTE: This routine must not be called out of domain scope. */ alpar@1: alpar@1: SYMBOL *take_member_sym /* returns value, not reference */ alpar@1: ( MPL *mpl, alpar@1: PARAMETER *par, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { MEMBER *memb; alpar@1: SYMBOL *value; alpar@1: /* find member in the parameter array */ alpar@1: memb = find_member(mpl, par->array, tuple); alpar@1: if (memb != NULL) alpar@1: { /* member exists, so just take its value */ alpar@1: value = copy_symbol(mpl, memb->value.sym); alpar@1: } alpar@1: else if (par->assign != NULL) alpar@1: { /* compute value using assignment expression */ alpar@1: value = eval_symbolic(mpl, par->assign); alpar@1: add: /* check that the value satisfies to all restrictions, assign alpar@1: it to new member, and add the member to the array */ alpar@1: check_value_sym(mpl, par, tuple, value); alpar@1: memb = add_member(mpl, par->array, copy_tuple(mpl, tuple)); alpar@1: memb->value.sym = copy_symbol(mpl, value); alpar@1: } alpar@1: else if (par->option != NULL) alpar@1: { /* compute default value */ alpar@1: value = eval_symbolic(mpl, par->option); alpar@1: goto add; alpar@1: } alpar@1: else if (par->defval != NULL) alpar@1: { /* take default value provided in the data section */ alpar@1: value = copy_symbol(mpl, par->defval); alpar@1: goto add; alpar@1: } alpar@1: else alpar@1: { /* no value is provided */ alpar@1: error(mpl, "no value for %s%s", par->name, format_tuple(mpl, alpar@1: '[', tuple)); alpar@1: } alpar@1: return value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_member_sym - evaluate symb. value assigned to parameter member. alpar@1: -- alpar@1: -- This routine evaluates a symbolic value assigned to given member of alpar@1: -- specified symbolic model parameter and returns it on exit. */ alpar@1: alpar@1: struct eval_sym_info alpar@1: { /* working info used by the routine eval_member_sym */ alpar@1: PARAMETER *par; alpar@1: /* model parameter */ alpar@1: TUPLE *tuple; alpar@1: /* n-tuple, which defines parameter member */ alpar@1: MEMBER *memb; alpar@1: /* normally this pointer is NULL; the routine uses this pointer alpar@1: to check data provided in the data section, in which case it alpar@1: points to a member currently checked; this check is performed alpar@1: automatically only once when a reference to any member occurs alpar@1: for the first time */ alpar@1: SYMBOL *value; alpar@1: /* evaluated symbolic value */ alpar@1: }; alpar@1: alpar@1: static void eval_sym_func(MPL *mpl, void *_info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: struct eval_sym_info *info = _info; alpar@1: if (info->memb != NULL) alpar@1: { /* checking call; check symbolic value being assigned */ alpar@1: check_value_sym(mpl, info->par, info->memb->tuple, alpar@1: info->memb->value.sym); alpar@1: } alpar@1: else alpar@1: { /* normal call; evaluate member, which has given n-tuple */ alpar@1: info->value = take_member_sym(mpl, info->par, info->tuple); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: SYMBOL *eval_member_sym /* returns value, not reference */ alpar@1: ( MPL *mpl, alpar@1: PARAMETER *par, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { /* this routine evaluates symbolic parameter member */ alpar@1: struct eval_sym_info _info, *info = &_info; alpar@1: xassert(par->type == A_SYMBOLIC); alpar@1: xassert(par->dim == tuple_dimen(mpl, tuple)); alpar@1: info->par = par; alpar@1: info->tuple = tuple; alpar@1: if (par->data == 1) alpar@1: { /* check data, which are provided in the data section, but not alpar@1: checked yet */ alpar@1: /* save pointer to the last array member; note that during the alpar@1: check new members may be added beyond the last member due to alpar@1: references to the same parameter from default expression as alpar@1: well as from expressions that define restricting conditions; alpar@1: however, values assigned to the new members will be checked alpar@1: by other routine, so we don't need to check them here */ alpar@1: MEMBER *tail = par->array->tail; alpar@1: /* change the data status to prevent infinite recursive loop alpar@1: due to references to the same parameter during the check */ alpar@1: par->data = 2; alpar@1: /* check values assigned to array members in the data section alpar@1: until the marked member has been reached */ alpar@1: for (info->memb = par->array->head; info->memb != NULL; alpar@1: info->memb = info->memb->next) alpar@1: { if (eval_within_domain(mpl, par->domain, info->memb->tuple, alpar@1: info, eval_sym_func)) alpar@1: out_of_domain(mpl, par->name, info->memb->tuple); alpar@1: if (info->memb == tail) break; alpar@1: } alpar@1: /* the check has been finished */ alpar@1: } alpar@1: /* evaluate member, which has given n-tuple */ alpar@1: info->memb = NULL; alpar@1: if (eval_within_domain(mpl, info->par->domain, info->tuple, info, alpar@1: eval_sym_func)) alpar@1: out_of_domain(mpl, par->name, info->tuple); alpar@1: /* bring evaluated value to the calling program */ alpar@1: return info->value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_whole_par - evaluate model parameter over entire domain. alpar@1: -- alpar@1: -- This routine evaluates all members of specified model parameter over alpar@1: -- entire domain. */ alpar@1: alpar@1: static int whole_par_func(MPL *mpl, void *info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: PARAMETER *par = (PARAMETER *)info; alpar@1: TUPLE *tuple = get_domain_tuple(mpl, par->domain); alpar@1: switch (par->type) alpar@1: { case A_NUMERIC: alpar@1: case A_INTEGER: alpar@1: case A_BINARY: alpar@1: eval_member_num(mpl, par, tuple); alpar@1: break; alpar@1: case A_SYMBOLIC: alpar@1: delete_symbol(mpl, eval_member_sym(mpl, par, tuple)); alpar@1: break; alpar@1: default: alpar@1: xassert(par != par); alpar@1: } alpar@1: delete_tuple(mpl, tuple); alpar@1: return 0; alpar@1: } alpar@1: alpar@1: void eval_whole_par(MPL *mpl, PARAMETER *par) alpar@1: { loop_within_domain(mpl, par->domain, par, whole_par_func); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean_parameter - clean model parameter. alpar@1: -- alpar@1: -- This routine cleans specified model parameter that assumes deleting alpar@1: -- all stuff dynamically allocated during the generation phase. */ alpar@1: alpar@1: void clean_parameter(MPL *mpl, PARAMETER *par) alpar@1: { CONDITION *cond; alpar@1: WITHIN *in; alpar@1: MEMBER *memb; alpar@1: /* clean subscript domain */ alpar@1: clean_domain(mpl, par->domain); alpar@1: /* clean pseudo-code for computing restricting conditions */ alpar@1: for (cond = par->cond; cond != NULL; cond = cond->next) alpar@1: clean_code(mpl, cond->code); alpar@1: /* clean pseudo-code for computing restricting supersets */ alpar@1: for (in = par->in; in != NULL; in = in->next) alpar@1: clean_code(mpl, in->code); alpar@1: /* clean pseudo-code for computing assigned value */ alpar@1: clean_code(mpl, par->assign); alpar@1: /* clean pseudo-code for computing default value */ alpar@1: clean_code(mpl, par->option); alpar@1: /* reset data status flag */ alpar@1: par->data = 0; alpar@1: /* delete default symbolic value */ alpar@1: if (par->defval != NULL) alpar@1: delete_symbol(mpl, par->defval), par->defval = NULL; alpar@1: /* delete content array */ alpar@1: for (memb = par->array->head; memb != NULL; memb = memb->next) alpar@1: delete_value(mpl, par->array->type, &memb->value); alpar@1: delete_array(mpl, par->array), par->array = NULL; alpar@1: return; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * MODEL VARIABLES * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- take_member_var - obtain reference to elemental variable. alpar@1: -- alpar@1: -- This routine obtains a reference to elemental variable assigned to alpar@1: -- given member of specified model variable and returns it on exit. If alpar@1: -- necessary, new elemental variable is created. alpar@1: -- alpar@1: -- NOTE: This routine must not be called out of domain scope. */ alpar@1: alpar@1: ELEMVAR *take_member_var /* returns reference */ alpar@1: ( MPL *mpl, alpar@1: VARIABLE *var, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { MEMBER *memb; alpar@1: ELEMVAR *refer; alpar@1: /* find member in the variable array */ alpar@1: memb = find_member(mpl, var->array, tuple); alpar@1: if (memb != NULL) alpar@1: { /* member exists, so just take the reference */ alpar@1: refer = memb->value.var; alpar@1: } alpar@1: else alpar@1: { /* member is referenced for the first time and therefore does alpar@1: not exist; create new elemental variable, assign it to new alpar@1: member, and add the member to the variable array */ alpar@1: memb = add_member(mpl, var->array, copy_tuple(mpl, tuple)); alpar@1: refer = (memb->value.var = alpar@1: dmp_get_atom(mpl->elemvars, sizeof(ELEMVAR))); alpar@1: refer->j = 0; alpar@1: refer->var = var; alpar@1: refer->memb = memb; alpar@1: /* compute lower bound */ alpar@1: if (var->lbnd == NULL) alpar@1: refer->lbnd = 0.0; alpar@1: else alpar@1: refer->lbnd = eval_numeric(mpl, var->lbnd); alpar@1: /* compute upper bound */ alpar@1: if (var->ubnd == NULL) alpar@1: refer->ubnd = 0.0; alpar@1: else if (var->ubnd == var->lbnd) alpar@1: refer->ubnd = refer->lbnd; alpar@1: else alpar@1: refer->ubnd = eval_numeric(mpl, var->ubnd); alpar@1: /* nullify working quantity */ alpar@1: refer->temp = 0.0; alpar@1: #if 1 /* 15/V-2010 */ alpar@1: /* solution has not been obtained by the solver yet */ alpar@1: refer->stat = 0; alpar@1: refer->prim = refer->dual = 0.0; alpar@1: #endif alpar@1: } alpar@1: return refer; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_member_var - evaluate reference to elemental variable. alpar@1: -- alpar@1: -- This routine evaluates a reference to elemental variable assigned to alpar@1: -- member of specified model variable and returns it on exit. */ alpar@1: alpar@1: struct eval_var_info alpar@1: { /* working info used by the routine eval_member_var */ alpar@1: VARIABLE *var; alpar@1: /* model variable */ alpar@1: TUPLE *tuple; alpar@1: /* n-tuple, which defines variable member */ alpar@1: ELEMVAR *refer; alpar@1: /* evaluated reference to elemental variable */ alpar@1: }; alpar@1: alpar@1: static void eval_var_func(MPL *mpl, void *_info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: struct eval_var_info *info = _info; alpar@1: info->refer = take_member_var(mpl, info->var, info->tuple); alpar@1: return; alpar@1: } alpar@1: alpar@1: ELEMVAR *eval_member_var /* returns reference */ alpar@1: ( MPL *mpl, alpar@1: VARIABLE *var, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { /* this routine evaluates variable member */ alpar@1: struct eval_var_info _info, *info = &_info; alpar@1: xassert(var->dim == tuple_dimen(mpl, tuple)); alpar@1: info->var = var; alpar@1: info->tuple = tuple; alpar@1: /* evaluate member, which has given n-tuple */ alpar@1: if (eval_within_domain(mpl, info->var->domain, info->tuple, info, alpar@1: eval_var_func)) alpar@1: out_of_domain(mpl, var->name, info->tuple); alpar@1: /* bring evaluated reference to the calling program */ alpar@1: return info->refer; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_whole_var - evaluate model variable over entire domain. alpar@1: -- alpar@1: -- This routine evaluates all members of specified model variable over alpar@1: -- entire domain. */ alpar@1: alpar@1: static int whole_var_func(MPL *mpl, void *info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: VARIABLE *var = (VARIABLE *)info; alpar@1: TUPLE *tuple = get_domain_tuple(mpl, var->domain); alpar@1: eval_member_var(mpl, var, tuple); alpar@1: delete_tuple(mpl, tuple); alpar@1: return 0; alpar@1: } alpar@1: alpar@1: void eval_whole_var(MPL *mpl, VARIABLE *var) alpar@1: { loop_within_domain(mpl, var->domain, var, whole_var_func); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean_variable - clean model variable. alpar@1: -- alpar@1: -- This routine cleans specified model variable that assumes deleting alpar@1: -- all stuff dynamically allocated during the generation phase. */ alpar@1: alpar@1: void clean_variable(MPL *mpl, VARIABLE *var) alpar@1: { MEMBER *memb; alpar@1: /* clean subscript domain */ alpar@1: clean_domain(mpl, var->domain); alpar@1: /* clean code for computing lower bound */ alpar@1: clean_code(mpl, var->lbnd); alpar@1: /* clean code for computing upper bound */ alpar@1: if (var->ubnd != var->lbnd) clean_code(mpl, var->ubnd); alpar@1: /* delete content array */ alpar@1: for (memb = var->array->head; memb != NULL; memb = memb->next) alpar@1: dmp_free_atom(mpl->elemvars, memb->value.var, sizeof(ELEMVAR)); alpar@1: delete_array(mpl, var->array), var->array = NULL; alpar@1: return; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * MODEL CONSTRAINTS AND OBJECTIVES * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- take_member_con - obtain reference to elemental constraint. alpar@1: -- alpar@1: -- This routine obtains a reference to elemental constraint assigned alpar@1: -- to given member of specified model constraint and returns it on exit. alpar@1: -- If necessary, new elemental constraint is created. alpar@1: -- alpar@1: -- NOTE: This routine must not be called out of domain scope. */ alpar@1: alpar@1: ELEMCON *take_member_con /* returns reference */ alpar@1: ( MPL *mpl, alpar@1: CONSTRAINT *con, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { MEMBER *memb; alpar@1: ELEMCON *refer; alpar@1: /* find member in the constraint array */ alpar@1: memb = find_member(mpl, con->array, tuple); alpar@1: if (memb != NULL) alpar@1: { /* member exists, so just take the reference */ alpar@1: refer = memb->value.con; alpar@1: } alpar@1: else alpar@1: { /* member is referenced for the first time and therefore does alpar@1: not exist; create new elemental constraint, assign it to new alpar@1: member, and add the member to the constraint array */ alpar@1: memb = add_member(mpl, con->array, copy_tuple(mpl, tuple)); alpar@1: refer = (memb->value.con = alpar@1: dmp_get_atom(mpl->elemcons, sizeof(ELEMCON))); alpar@1: refer->i = 0; alpar@1: refer->con = con; alpar@1: refer->memb = memb; alpar@1: /* compute linear form */ alpar@1: xassert(con->code != NULL); alpar@1: refer->form = eval_formula(mpl, con->code); alpar@1: /* compute lower and upper bounds */ alpar@1: if (con->lbnd == NULL && con->ubnd == NULL) alpar@1: { /* objective has no bounds */ alpar@1: double temp; alpar@1: xassert(con->type == A_MINIMIZE || con->type == A_MAXIMIZE); alpar@1: /* carry the constant term to the right-hand side */ alpar@1: refer->form = remove_constant(mpl, refer->form, &temp); alpar@1: refer->lbnd = refer->ubnd = - temp; alpar@1: } alpar@1: else if (con->lbnd != NULL && con->ubnd == NULL) alpar@1: { /* constraint a * x + b >= c * y + d is transformed to the alpar@1: standard form a * x - c * y >= d - b */ alpar@1: double temp; alpar@1: xassert(con->type == A_CONSTRAINT); alpar@1: refer->form = linear_comb(mpl, alpar@1: +1.0, refer->form, alpar@1: -1.0, eval_formula(mpl, con->lbnd)); alpar@1: refer->form = remove_constant(mpl, refer->form, &temp); alpar@1: refer->lbnd = - temp; alpar@1: refer->ubnd = 0.0; alpar@1: } alpar@1: else if (con->lbnd == NULL && con->ubnd != NULL) alpar@1: { /* constraint a * x + b <= c * y + d is transformed to the alpar@1: standard form a * x - c * y <= d - b */ alpar@1: double temp; alpar@1: xassert(con->type == A_CONSTRAINT); alpar@1: refer->form = linear_comb(mpl, alpar@1: +1.0, refer->form, alpar@1: -1.0, eval_formula(mpl, con->ubnd)); alpar@1: refer->form = remove_constant(mpl, refer->form, &temp); alpar@1: refer->lbnd = 0.0; alpar@1: refer->ubnd = - temp; alpar@1: } alpar@1: else if (con->lbnd == con->ubnd) alpar@1: { /* constraint a * x + b = c * y + d is transformed to the alpar@1: standard form a * x - c * y = d - b */ alpar@1: double temp; alpar@1: xassert(con->type == A_CONSTRAINT); alpar@1: refer->form = linear_comb(mpl, alpar@1: +1.0, refer->form, alpar@1: -1.0, eval_formula(mpl, con->lbnd)); alpar@1: refer->form = remove_constant(mpl, refer->form, &temp); alpar@1: refer->lbnd = refer->ubnd = - temp; alpar@1: } alpar@1: else alpar@1: { /* ranged constraint c <= a * x + b <= d is transformed to alpar@1: the standard form c - b <= a * x <= d - b */ alpar@1: double temp, temp1, temp2; alpar@1: xassert(con->type == A_CONSTRAINT); alpar@1: refer->form = remove_constant(mpl, refer->form, &temp); alpar@1: xassert(remove_constant(mpl, eval_formula(mpl, con->lbnd), alpar@1: &temp1) == NULL); alpar@1: xassert(remove_constant(mpl, eval_formula(mpl, con->ubnd), alpar@1: &temp2) == NULL); alpar@1: refer->lbnd = fp_sub(mpl, temp1, temp); alpar@1: refer->ubnd = fp_sub(mpl, temp2, temp); alpar@1: } alpar@1: #if 1 /* 15/V-2010 */ alpar@1: /* solution has not been obtained by the solver yet */ alpar@1: refer->stat = 0; alpar@1: refer->prim = refer->dual = 0.0; alpar@1: #endif alpar@1: } alpar@1: return refer; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_member_con - evaluate reference to elemental constraint. alpar@1: -- alpar@1: -- This routine evaluates a reference to elemental constraint assigned alpar@1: -- to member of specified model constraint and returns it on exit. */ alpar@1: alpar@1: struct eval_con_info alpar@1: { /* working info used by the routine eval_member_con */ alpar@1: CONSTRAINT *con; alpar@1: /* model constraint */ alpar@1: TUPLE *tuple; alpar@1: /* n-tuple, which defines constraint member */ alpar@1: ELEMCON *refer; alpar@1: /* evaluated reference to elemental constraint */ alpar@1: }; alpar@1: alpar@1: static void eval_con_func(MPL *mpl, void *_info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: struct eval_con_info *info = _info; alpar@1: info->refer = take_member_con(mpl, info->con, info->tuple); alpar@1: return; alpar@1: } alpar@1: alpar@1: ELEMCON *eval_member_con /* returns reference */ alpar@1: ( MPL *mpl, alpar@1: CONSTRAINT *con, /* not changed */ alpar@1: TUPLE *tuple /* not changed */ alpar@1: ) alpar@1: { /* this routine evaluates constraint member */ alpar@1: struct eval_con_info _info, *info = &_info; alpar@1: xassert(con->dim == tuple_dimen(mpl, tuple)); alpar@1: info->con = con; alpar@1: info->tuple = tuple; alpar@1: /* evaluate member, which has given n-tuple */ alpar@1: if (eval_within_domain(mpl, info->con->domain, info->tuple, info, alpar@1: eval_con_func)) alpar@1: out_of_domain(mpl, con->name, info->tuple); alpar@1: /* bring evaluated reference to the calling program */ alpar@1: return info->refer; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_whole_con - evaluate model constraint over entire domain. alpar@1: -- alpar@1: -- This routine evaluates all members of specified model constraint over alpar@1: -- entire domain. */ alpar@1: alpar@1: static int whole_con_func(MPL *mpl, void *info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: CONSTRAINT *con = (CONSTRAINT *)info; alpar@1: TUPLE *tuple = get_domain_tuple(mpl, con->domain); alpar@1: eval_member_con(mpl, con, tuple); alpar@1: delete_tuple(mpl, tuple); alpar@1: return 0; alpar@1: } alpar@1: alpar@1: void eval_whole_con(MPL *mpl, CONSTRAINT *con) alpar@1: { loop_within_domain(mpl, con->domain, con, whole_con_func); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean_constraint - clean model constraint. alpar@1: -- alpar@1: -- This routine cleans specified model constraint that assumes deleting alpar@1: -- all stuff dynamically allocated during the generation phase. */ alpar@1: alpar@1: void clean_constraint(MPL *mpl, CONSTRAINT *con) alpar@1: { MEMBER *memb; alpar@1: /* clean subscript domain */ alpar@1: clean_domain(mpl, con->domain); alpar@1: /* clean code for computing main linear form */ alpar@1: clean_code(mpl, con->code); alpar@1: /* clean code for computing lower bound */ alpar@1: clean_code(mpl, con->lbnd); alpar@1: /* clean code for computing upper bound */ alpar@1: if (con->ubnd != con->lbnd) clean_code(mpl, con->ubnd); alpar@1: /* delete content array */ alpar@1: for (memb = con->array->head; memb != NULL; memb = memb->next) alpar@1: { delete_formula(mpl, memb->value.con->form); alpar@1: dmp_free_atom(mpl->elemcons, memb->value.con, sizeof(ELEMCON)); alpar@1: } alpar@1: delete_array(mpl, con->array), con->array = NULL; alpar@1: return; alpar@1: } alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * PSEUDO-CODE * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_numeric - evaluate pseudo-code to determine numeric value. alpar@1: -- alpar@1: -- This routine evaluates specified pseudo-code to determine resultant alpar@1: -- numeric value, which is returned on exit. */ alpar@1: alpar@1: struct iter_num_info alpar@1: { /* working info used by the routine iter_num_func */ alpar@1: CODE *code; alpar@1: /* pseudo-code for iterated operation to be performed */ alpar@1: double value; alpar@1: /* resultant value */ alpar@1: }; alpar@1: alpar@1: static int iter_num_func(MPL *mpl, void *_info) alpar@1: { /* this is auxiliary routine used to perform iterated operation alpar@1: on numeric "integrand" within domain scope */ alpar@1: struct iter_num_info *info = _info; alpar@1: double temp; alpar@1: temp = eval_numeric(mpl, info->code->arg.loop.x); alpar@1: switch (info->code->op) alpar@1: { case O_SUM: alpar@1: /* summation over domain */ alpar@1: info->value = fp_add(mpl, info->value, temp); alpar@1: break; alpar@1: case O_PROD: alpar@1: /* multiplication over domain */ alpar@1: info->value = fp_mul(mpl, info->value, temp); alpar@1: break; alpar@1: case O_MINIMUM: alpar@1: /* minimum over domain */ alpar@1: if (info->value > temp) info->value = temp; alpar@1: break; alpar@1: case O_MAXIMUM: alpar@1: /* maximum over domain */ alpar@1: if (info->value < temp) info->value = temp; alpar@1: break; alpar@1: default: alpar@1: xassert(info != info); alpar@1: } alpar@1: return 0; alpar@1: } alpar@1: alpar@1: double eval_numeric(MPL *mpl, CODE *code) alpar@1: { double value; alpar@1: xassert(code != NULL); alpar@1: xassert(code->type == A_NUMERIC); alpar@1: xassert(code->dim == 0); alpar@1: /* if the operation has a side effect, invalidate and delete the alpar@1: resultant value */ alpar@1: if (code->vflag && code->valid) alpar@1: { code->valid = 0; alpar@1: delete_value(mpl, code->type, &code->value); alpar@1: } alpar@1: /* if resultant value is valid, no evaluation is needed */ alpar@1: if (code->valid) alpar@1: { value = code->value.num; alpar@1: goto done; alpar@1: } alpar@1: /* evaluate pseudo-code recursively */ alpar@1: switch (code->op) alpar@1: { case O_NUMBER: alpar@1: /* take floating-point number */ alpar@1: value = code->arg.num; alpar@1: break; alpar@1: case O_MEMNUM: alpar@1: /* take member of numeric parameter */ alpar@1: { TUPLE *tuple; alpar@1: ARG_LIST *e; alpar@1: tuple = create_tuple(mpl); alpar@1: for (e = code->arg.par.list; e != NULL; e = e->next) alpar@1: tuple = expand_tuple(mpl, tuple, eval_symbolic(mpl, alpar@1: e->x)); alpar@1: value = eval_member_num(mpl, code->arg.par.par, tuple); alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: break; alpar@1: case O_MEMVAR: alpar@1: /* take computed value of elemental variable */ alpar@1: { TUPLE *tuple; alpar@1: ARG_LIST *e; alpar@1: #if 1 /* 15/V-2010 */ alpar@1: ELEMVAR *var; alpar@1: #endif alpar@1: tuple = create_tuple(mpl); alpar@1: for (e = code->arg.var.list; e != NULL; e = e->next) alpar@1: tuple = expand_tuple(mpl, tuple, eval_symbolic(mpl, alpar@1: e->x)); alpar@1: #if 0 /* 15/V-2010 */ alpar@1: value = eval_member_var(mpl, code->arg.var.var, tuple) alpar@1: ->value; alpar@1: #else alpar@1: var = eval_member_var(mpl, code->arg.var.var, tuple); alpar@1: switch (code->arg.var.suff) alpar@1: { case DOT_LB: alpar@1: if (var->var->lbnd == NULL) alpar@1: value = -DBL_MAX; alpar@1: else alpar@1: value = var->lbnd; alpar@1: break; alpar@1: case DOT_UB: alpar@1: if (var->var->ubnd == NULL) alpar@1: value = +DBL_MAX; alpar@1: else alpar@1: value = var->ubnd; alpar@1: break; alpar@1: case DOT_STATUS: alpar@1: value = var->stat; alpar@1: break; alpar@1: case DOT_VAL: alpar@1: value = var->prim; alpar@1: break; alpar@1: case DOT_DUAL: alpar@1: value = var->dual; alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: #endif alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: break; alpar@1: #if 1 /* 15/V-2010 */ alpar@1: case O_MEMCON: alpar@1: /* take computed value of elemental constraint */ alpar@1: { TUPLE *tuple; alpar@1: ARG_LIST *e; alpar@1: ELEMCON *con; alpar@1: tuple = create_tuple(mpl); alpar@1: for (e = code->arg.con.list; e != NULL; e = e->next) alpar@1: tuple = expand_tuple(mpl, tuple, eval_symbolic(mpl, alpar@1: e->x)); alpar@1: con = eval_member_con(mpl, code->arg.con.con, tuple); alpar@1: switch (code->arg.con.suff) alpar@1: { case DOT_LB: alpar@1: if (con->con->lbnd == NULL) alpar@1: value = -DBL_MAX; alpar@1: else alpar@1: value = con->lbnd; alpar@1: break; alpar@1: case DOT_UB: alpar@1: if (con->con->ubnd == NULL) alpar@1: value = +DBL_MAX; alpar@1: else alpar@1: value = con->ubnd; alpar@1: break; alpar@1: case DOT_STATUS: alpar@1: value = con->stat; alpar@1: break; alpar@1: case DOT_VAL: alpar@1: value = con->prim; alpar@1: break; alpar@1: case DOT_DUAL: alpar@1: value = con->dual; alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: break; alpar@1: #endif alpar@1: case O_IRAND224: alpar@1: /* pseudo-random in [0, 2^24-1] */ alpar@1: value = fp_irand224(mpl); alpar@1: break; alpar@1: case O_UNIFORM01: alpar@1: /* pseudo-random in [0, 1) */ alpar@1: value = fp_uniform01(mpl); alpar@1: break; alpar@1: case O_NORMAL01: alpar@1: /* gaussian random, mu = 0, sigma = 1 */ alpar@1: value = fp_normal01(mpl); alpar@1: break; alpar@1: case O_GMTIME: alpar@1: /* current calendar time */ alpar@1: value = fn_gmtime(mpl); alpar@1: break; alpar@1: case O_CVTNUM: alpar@1: /* conversion to numeric */ alpar@1: { SYMBOL *sym; alpar@1: sym = eval_symbolic(mpl, code->arg.arg.x); alpar@1: #if 0 /* 23/XI-2008 */ alpar@1: if (sym->str != NULL) alpar@1: error(mpl, "cannot convert %s to floating-point numbe" alpar@1: "r", format_symbol(mpl, sym)); alpar@1: value = sym->num; alpar@1: #else alpar@1: if (sym->str == NULL) alpar@1: value = sym->num; alpar@1: else alpar@1: { if (str2num(sym->str, &value)) alpar@1: error(mpl, "cannot convert %s to floating-point nu" alpar@1: "mber", format_symbol(mpl, sym)); alpar@1: } alpar@1: #endif alpar@1: delete_symbol(mpl, sym); alpar@1: } alpar@1: break; alpar@1: case O_PLUS: alpar@1: /* unary plus */ alpar@1: value = + eval_numeric(mpl, code->arg.arg.x); alpar@1: break; alpar@1: case O_MINUS: alpar@1: /* unary minus */ alpar@1: value = - eval_numeric(mpl, code->arg.arg.x); alpar@1: break; alpar@1: case O_ABS: alpar@1: /* absolute value */ alpar@1: value = fabs(eval_numeric(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_CEIL: alpar@1: /* round upward ("ceiling of x") */ alpar@1: value = ceil(eval_numeric(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_FLOOR: alpar@1: /* round downward ("floor of x") */ alpar@1: value = floor(eval_numeric(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_EXP: alpar@1: /* base-e exponential */ alpar@1: value = fp_exp(mpl, eval_numeric(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_LOG: alpar@1: /* natural logarithm */ alpar@1: value = fp_log(mpl, eval_numeric(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_LOG10: alpar@1: /* common (decimal) logarithm */ alpar@1: value = fp_log10(mpl, eval_numeric(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_SQRT: alpar@1: /* square root */ alpar@1: value = fp_sqrt(mpl, eval_numeric(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_SIN: alpar@1: /* trigonometric sine */ alpar@1: value = fp_sin(mpl, eval_numeric(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_COS: alpar@1: /* trigonometric cosine */ alpar@1: value = fp_cos(mpl, eval_numeric(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_ATAN: alpar@1: /* trigonometric arctangent (one argument) */ alpar@1: value = fp_atan(mpl, eval_numeric(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_ATAN2: alpar@1: /* trigonometric arctangent (two arguments) */ alpar@1: value = fp_atan2(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_ROUND: alpar@1: /* round to nearest integer */ alpar@1: value = fp_round(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), 0.0); alpar@1: break; alpar@1: case O_ROUND2: alpar@1: /* round to n fractional digits */ alpar@1: value = fp_round(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_TRUNC: alpar@1: /* truncate to nearest integer */ alpar@1: value = fp_trunc(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), 0.0); alpar@1: break; alpar@1: case O_TRUNC2: alpar@1: /* truncate to n fractional digits */ alpar@1: value = fp_trunc(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_ADD: alpar@1: /* addition */ alpar@1: value = fp_add(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_SUB: alpar@1: /* subtraction */ alpar@1: value = fp_sub(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_LESS: alpar@1: /* non-negative subtraction */ alpar@1: value = fp_less(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_MUL: alpar@1: /* multiplication */ alpar@1: value = fp_mul(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_DIV: alpar@1: /* division */ alpar@1: value = fp_div(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_IDIV: alpar@1: /* quotient of exact division */ alpar@1: value = fp_idiv(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_MOD: alpar@1: /* remainder of exact division */ alpar@1: value = fp_mod(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_POWER: alpar@1: /* exponentiation (raise to power) */ alpar@1: value = fp_power(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_UNIFORM: alpar@1: /* pseudo-random in [a, b) */ alpar@1: value = fp_uniform(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_NORMAL: alpar@1: /* gaussian random, given mu and sigma */ alpar@1: value = fp_normal(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_CARD: alpar@1: { ELEMSET *set; alpar@1: set = eval_elemset(mpl, code->arg.arg.x); alpar@1: value = set->size; alpar@1: delete_array(mpl, set); alpar@1: } alpar@1: break; alpar@1: case O_LENGTH: alpar@1: { SYMBOL *sym; alpar@1: char str[MAX_LENGTH+1]; alpar@1: sym = eval_symbolic(mpl, code->arg.arg.x); alpar@1: if (sym->str == NULL) alpar@1: sprintf(str, "%.*g", DBL_DIG, sym->num); alpar@1: else alpar@1: fetch_string(mpl, sym->str, str); alpar@1: delete_symbol(mpl, sym); alpar@1: value = strlen(str); alpar@1: } alpar@1: break; alpar@1: case O_STR2TIME: alpar@1: { SYMBOL *sym; alpar@1: char str[MAX_LENGTH+1], fmt[MAX_LENGTH+1]; alpar@1: sym = eval_symbolic(mpl, code->arg.arg.x); alpar@1: if (sym->str == NULL) alpar@1: sprintf(str, "%.*g", DBL_DIG, sym->num); alpar@1: else alpar@1: fetch_string(mpl, sym->str, str); alpar@1: delete_symbol(mpl, sym); alpar@1: sym = eval_symbolic(mpl, code->arg.arg.y); alpar@1: if (sym->str == NULL) alpar@1: sprintf(fmt, "%.*g", DBL_DIG, sym->num); alpar@1: else alpar@1: fetch_string(mpl, sym->str, fmt); alpar@1: delete_symbol(mpl, sym); alpar@1: value = fn_str2time(mpl, str, fmt); alpar@1: } alpar@1: break; alpar@1: case O_FORK: alpar@1: /* if-then-else */ alpar@1: if (eval_logical(mpl, code->arg.arg.x)) alpar@1: value = eval_numeric(mpl, code->arg.arg.y); alpar@1: else if (code->arg.arg.z == NULL) alpar@1: value = 0.0; alpar@1: else alpar@1: value = eval_numeric(mpl, code->arg.arg.z); alpar@1: break; alpar@1: case O_MIN: alpar@1: /* minimal value (n-ary) */ alpar@1: { ARG_LIST *e; alpar@1: double temp; alpar@1: value = +DBL_MAX; alpar@1: for (e = code->arg.list; e != NULL; e = e->next) alpar@1: { temp = eval_numeric(mpl, e->x); alpar@1: if (value > temp) value = temp; alpar@1: } alpar@1: } alpar@1: break; alpar@1: case O_MAX: alpar@1: /* maximal value (n-ary) */ alpar@1: { ARG_LIST *e; alpar@1: double temp; alpar@1: value = -DBL_MAX; alpar@1: for (e = code->arg.list; e != NULL; e = e->next) alpar@1: { temp = eval_numeric(mpl, e->x); alpar@1: if (value < temp) value = temp; alpar@1: } alpar@1: } alpar@1: break; alpar@1: case O_SUM: alpar@1: /* summation over domain */ alpar@1: { struct iter_num_info _info, *info = &_info; alpar@1: info->code = code; alpar@1: info->value = 0.0; alpar@1: loop_within_domain(mpl, code->arg.loop.domain, info, alpar@1: iter_num_func); alpar@1: value = info->value; alpar@1: } alpar@1: break; alpar@1: case O_PROD: alpar@1: /* multiplication over domain */ alpar@1: { struct iter_num_info _info, *info = &_info; alpar@1: info->code = code; alpar@1: info->value = 1.0; alpar@1: loop_within_domain(mpl, code->arg.loop.domain, info, alpar@1: iter_num_func); alpar@1: value = info->value; alpar@1: } alpar@1: break; alpar@1: case O_MINIMUM: alpar@1: /* minimum over domain */ alpar@1: { struct iter_num_info _info, *info = &_info; alpar@1: info->code = code; alpar@1: info->value = +DBL_MAX; alpar@1: loop_within_domain(mpl, code->arg.loop.domain, info, alpar@1: iter_num_func); alpar@1: if (info->value == +DBL_MAX) alpar@1: error(mpl, "min{} over empty set; result undefined"); alpar@1: value = info->value; alpar@1: } alpar@1: break; alpar@1: case O_MAXIMUM: alpar@1: /* maximum over domain */ alpar@1: { struct iter_num_info _info, *info = &_info; alpar@1: info->code = code; alpar@1: info->value = -DBL_MAX; alpar@1: loop_within_domain(mpl, code->arg.loop.domain, info, alpar@1: iter_num_func); alpar@1: if (info->value == -DBL_MAX) alpar@1: error(mpl, "max{} over empty set; result undefined"); alpar@1: value = info->value; alpar@1: } alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: /* save resultant value */ alpar@1: xassert(!code->valid); alpar@1: code->valid = 1; alpar@1: code->value.num = value; alpar@1: done: return value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_symbolic - evaluate pseudo-code to determine symbolic value. alpar@1: -- alpar@1: -- This routine evaluates specified pseudo-code to determine resultant alpar@1: -- symbolic value, which is returned on exit. */ alpar@1: alpar@1: SYMBOL *eval_symbolic(MPL *mpl, CODE *code) alpar@1: { SYMBOL *value; alpar@1: xassert(code != NULL); alpar@1: xassert(code->type == A_SYMBOLIC); alpar@1: xassert(code->dim == 0); alpar@1: /* if the operation has a side effect, invalidate and delete the alpar@1: resultant value */ alpar@1: if (code->vflag && code->valid) alpar@1: { code->valid = 0; alpar@1: delete_value(mpl, code->type, &code->value); alpar@1: } alpar@1: /* if resultant value is valid, no evaluation is needed */ alpar@1: if (code->valid) alpar@1: { value = copy_symbol(mpl, code->value.sym); alpar@1: goto done; alpar@1: } alpar@1: /* evaluate pseudo-code recursively */ alpar@1: switch (code->op) alpar@1: { case O_STRING: alpar@1: /* take character string */ alpar@1: value = create_symbol_str(mpl, create_string(mpl, alpar@1: code->arg.str)); alpar@1: break; alpar@1: case O_INDEX: alpar@1: /* take dummy index */ alpar@1: xassert(code->arg.index.slot->value != NULL); alpar@1: value = copy_symbol(mpl, code->arg.index.slot->value); alpar@1: break; alpar@1: case O_MEMSYM: alpar@1: /* take member of symbolic parameter */ alpar@1: { TUPLE *tuple; alpar@1: ARG_LIST *e; alpar@1: tuple = create_tuple(mpl); alpar@1: for (e = code->arg.par.list; e != NULL; e = e->next) alpar@1: tuple = expand_tuple(mpl, tuple, eval_symbolic(mpl, alpar@1: e->x)); alpar@1: value = eval_member_sym(mpl, code->arg.par.par, tuple); alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: break; alpar@1: case O_CVTSYM: alpar@1: /* conversion to symbolic */ alpar@1: value = create_symbol_num(mpl, eval_numeric(mpl, alpar@1: code->arg.arg.x)); alpar@1: break; alpar@1: case O_CONCAT: alpar@1: /* concatenation */ alpar@1: value = concat_symbols(mpl, alpar@1: eval_symbolic(mpl, code->arg.arg.x), alpar@1: eval_symbolic(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_FORK: alpar@1: /* if-then-else */ alpar@1: if (eval_logical(mpl, code->arg.arg.x)) alpar@1: value = eval_symbolic(mpl, code->arg.arg.y); alpar@1: else if (code->arg.arg.z == NULL) alpar@1: value = create_symbol_num(mpl, 0.0); alpar@1: else alpar@1: value = eval_symbolic(mpl, code->arg.arg.z); alpar@1: break; alpar@1: case O_SUBSTR: alpar@1: case O_SUBSTR3: alpar@1: { double pos, len; alpar@1: char str[MAX_LENGTH+1]; alpar@1: value = eval_symbolic(mpl, code->arg.arg.x); alpar@1: if (value->str == NULL) alpar@1: sprintf(str, "%.*g", DBL_DIG, value->num); alpar@1: else alpar@1: fetch_string(mpl, value->str, str); alpar@1: delete_symbol(mpl, value); alpar@1: if (code->op == O_SUBSTR) alpar@1: { pos = eval_numeric(mpl, code->arg.arg.y); alpar@1: if (pos != floor(pos)) alpar@1: error(mpl, "substr('...', %.*g); non-integer secon" alpar@1: "d argument", DBL_DIG, pos); alpar@1: if (pos < 1 || pos > strlen(str) + 1) alpar@1: error(mpl, "substr('...', %.*g); substring out of " alpar@1: "range", DBL_DIG, pos); alpar@1: } alpar@1: else alpar@1: { pos = eval_numeric(mpl, code->arg.arg.y); alpar@1: len = eval_numeric(mpl, code->arg.arg.z); alpar@1: if (pos != floor(pos) || len != floor(len)) alpar@1: error(mpl, "substr('...', %.*g, %.*g); non-integer" alpar@1: " second and/or third argument", DBL_DIG, pos, alpar@1: DBL_DIG, len); alpar@1: if (pos < 1 || len < 0 || pos + len > strlen(str) + 1) alpar@1: error(mpl, "substr('...', %.*g, %.*g); substring o" alpar@1: "ut of range", DBL_DIG, pos, DBL_DIG, len); alpar@1: str[(int)pos + (int)len - 1] = '\0'; alpar@1: } alpar@1: value = create_symbol_str(mpl, create_string(mpl, str + alpar@1: (int)pos - 1)); alpar@1: } alpar@1: break; alpar@1: case O_TIME2STR: alpar@1: { double num; alpar@1: SYMBOL *sym; alpar@1: char str[MAX_LENGTH+1], fmt[MAX_LENGTH+1]; alpar@1: num = eval_numeric(mpl, code->arg.arg.x); alpar@1: sym = eval_symbolic(mpl, code->arg.arg.y); alpar@1: if (sym->str == NULL) alpar@1: sprintf(fmt, "%.*g", DBL_DIG, sym->num); alpar@1: else alpar@1: fetch_string(mpl, sym->str, fmt); alpar@1: delete_symbol(mpl, sym); alpar@1: fn_time2str(mpl, str, num, fmt); alpar@1: value = create_symbol_str(mpl, create_string(mpl, str)); alpar@1: } alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: /* save resultant value */ alpar@1: xassert(!code->valid); alpar@1: code->valid = 1; alpar@1: code->value.sym = copy_symbol(mpl, value); alpar@1: done: return value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_logical - evaluate pseudo-code to determine logical value. alpar@1: -- alpar@1: -- This routine evaluates specified pseudo-code to determine resultant alpar@1: -- logical value, which is returned on exit. */ alpar@1: alpar@1: struct iter_log_info alpar@1: { /* working info used by the routine iter_log_func */ alpar@1: CODE *code; alpar@1: /* pseudo-code for iterated operation to be performed */ alpar@1: int value; alpar@1: /* resultant value */ alpar@1: }; alpar@1: alpar@1: static int iter_log_func(MPL *mpl, void *_info) alpar@1: { /* this is auxiliary routine used to perform iterated operation alpar@1: on logical "integrand" within domain scope */ alpar@1: struct iter_log_info *info = _info; alpar@1: int ret = 0; alpar@1: switch (info->code->op) alpar@1: { case O_FORALL: alpar@1: /* conjunction over domain */ alpar@1: info->value &= eval_logical(mpl, info->code->arg.loop.x); alpar@1: if (!info->value) ret = 1; alpar@1: break; alpar@1: case O_EXISTS: alpar@1: /* disjunction over domain */ alpar@1: info->value |= eval_logical(mpl, info->code->arg.loop.x); alpar@1: if (info->value) ret = 1; alpar@1: break; alpar@1: default: alpar@1: xassert(info != info); alpar@1: } alpar@1: return ret; alpar@1: } alpar@1: alpar@1: int eval_logical(MPL *mpl, CODE *code) alpar@1: { int value; alpar@1: xassert(code->type == A_LOGICAL); alpar@1: xassert(code->dim == 0); alpar@1: /* if the operation has a side effect, invalidate and delete the alpar@1: resultant value */ alpar@1: if (code->vflag && code->valid) alpar@1: { code->valid = 0; alpar@1: delete_value(mpl, code->type, &code->value); alpar@1: } alpar@1: /* if resultant value is valid, no evaluation is needed */ alpar@1: if (code->valid) alpar@1: { value = code->value.bit; alpar@1: goto done; alpar@1: } alpar@1: /* evaluate pseudo-code recursively */ alpar@1: switch (code->op) alpar@1: { case O_CVTLOG: alpar@1: /* conversion to logical */ alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) != 0.0); alpar@1: break; alpar@1: case O_NOT: alpar@1: /* negation (logical "not") */ alpar@1: value = !eval_logical(mpl, code->arg.arg.x); alpar@1: break; alpar@1: case O_LT: alpar@1: /* comparison on 'less than' */ alpar@1: #if 0 /* 02/VIII-2008 */ alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) < alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: #else alpar@1: xassert(code->arg.arg.x != NULL); alpar@1: if (code->arg.arg.x->type == A_NUMERIC) alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) < alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: else alpar@1: { SYMBOL *sym1 = eval_symbolic(mpl, code->arg.arg.x); alpar@1: SYMBOL *sym2 = eval_symbolic(mpl, code->arg.arg.y); alpar@1: value = (compare_symbols(mpl, sym1, sym2) < 0); alpar@1: delete_symbol(mpl, sym1); alpar@1: delete_symbol(mpl, sym2); alpar@1: } alpar@1: #endif alpar@1: break; alpar@1: case O_LE: alpar@1: /* comparison on 'not greater than' */ alpar@1: #if 0 /* 02/VIII-2008 */ alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) <= alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: #else alpar@1: xassert(code->arg.arg.x != NULL); alpar@1: if (code->arg.arg.x->type == A_NUMERIC) alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) <= alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: else alpar@1: { SYMBOL *sym1 = eval_symbolic(mpl, code->arg.arg.x); alpar@1: SYMBOL *sym2 = eval_symbolic(mpl, code->arg.arg.y); alpar@1: value = (compare_symbols(mpl, sym1, sym2) <= 0); alpar@1: delete_symbol(mpl, sym1); alpar@1: delete_symbol(mpl, sym2); alpar@1: } alpar@1: #endif alpar@1: break; alpar@1: case O_EQ: alpar@1: /* comparison on 'equal to' */ alpar@1: xassert(code->arg.arg.x != NULL); alpar@1: if (code->arg.arg.x->type == A_NUMERIC) alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) == alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: else alpar@1: { SYMBOL *sym1 = eval_symbolic(mpl, code->arg.arg.x); alpar@1: SYMBOL *sym2 = eval_symbolic(mpl, code->arg.arg.y); alpar@1: value = (compare_symbols(mpl, sym1, sym2) == 0); alpar@1: delete_symbol(mpl, sym1); alpar@1: delete_symbol(mpl, sym2); alpar@1: } alpar@1: break; alpar@1: case O_GE: alpar@1: /* comparison on 'not less than' */ alpar@1: #if 0 /* 02/VIII-2008 */ alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) >= alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: #else alpar@1: xassert(code->arg.arg.x != NULL); alpar@1: if (code->arg.arg.x->type == A_NUMERIC) alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) >= alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: else alpar@1: { SYMBOL *sym1 = eval_symbolic(mpl, code->arg.arg.x); alpar@1: SYMBOL *sym2 = eval_symbolic(mpl, code->arg.arg.y); alpar@1: value = (compare_symbols(mpl, sym1, sym2) >= 0); alpar@1: delete_symbol(mpl, sym1); alpar@1: delete_symbol(mpl, sym2); alpar@1: } alpar@1: #endif alpar@1: break; alpar@1: case O_GT: alpar@1: /* comparison on 'greater than' */ alpar@1: #if 0 /* 02/VIII-2008 */ alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) > alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: #else alpar@1: xassert(code->arg.arg.x != NULL); alpar@1: if (code->arg.arg.x->type == A_NUMERIC) alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) > alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: else alpar@1: { SYMBOL *sym1 = eval_symbolic(mpl, code->arg.arg.x); alpar@1: SYMBOL *sym2 = eval_symbolic(mpl, code->arg.arg.y); alpar@1: value = (compare_symbols(mpl, sym1, sym2) > 0); alpar@1: delete_symbol(mpl, sym1); alpar@1: delete_symbol(mpl, sym2); alpar@1: } alpar@1: #endif alpar@1: break; alpar@1: case O_NE: alpar@1: /* comparison on 'not equal to' */ alpar@1: xassert(code->arg.arg.x != NULL); alpar@1: if (code->arg.arg.x->type == A_NUMERIC) alpar@1: value = (eval_numeric(mpl, code->arg.arg.x) != alpar@1: eval_numeric(mpl, code->arg.arg.y)); alpar@1: else alpar@1: { SYMBOL *sym1 = eval_symbolic(mpl, code->arg.arg.x); alpar@1: SYMBOL *sym2 = eval_symbolic(mpl, code->arg.arg.y); alpar@1: value = (compare_symbols(mpl, sym1, sym2) != 0); alpar@1: delete_symbol(mpl, sym1); alpar@1: delete_symbol(mpl, sym2); alpar@1: } alpar@1: break; alpar@1: case O_AND: alpar@1: /* conjunction (logical "and") */ alpar@1: value = eval_logical(mpl, code->arg.arg.x) && alpar@1: eval_logical(mpl, code->arg.arg.y); alpar@1: break; alpar@1: case O_OR: alpar@1: /* disjunction (logical "or") */ alpar@1: value = eval_logical(mpl, code->arg.arg.x) || alpar@1: eval_logical(mpl, code->arg.arg.y); alpar@1: break; alpar@1: case O_IN: alpar@1: /* test on 'x in Y' */ alpar@1: { TUPLE *tuple; alpar@1: tuple = eval_tuple(mpl, code->arg.arg.x); alpar@1: value = is_member(mpl, code->arg.arg.y, tuple); alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: break; alpar@1: case O_NOTIN: alpar@1: /* test on 'x not in Y' */ alpar@1: { TUPLE *tuple; alpar@1: tuple = eval_tuple(mpl, code->arg.arg.x); alpar@1: value = !is_member(mpl, code->arg.arg.y, tuple); alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: break; alpar@1: case O_WITHIN: alpar@1: /* test on 'X within Y' */ alpar@1: { ELEMSET *set; alpar@1: MEMBER *memb; alpar@1: set = eval_elemset(mpl, code->arg.arg.x); alpar@1: value = 1; alpar@1: for (memb = set->head; memb != NULL; memb = memb->next) alpar@1: { if (!is_member(mpl, code->arg.arg.y, memb->tuple)) alpar@1: { value = 0; alpar@1: break; alpar@1: } alpar@1: } alpar@1: delete_elemset(mpl, set); alpar@1: } alpar@1: break; alpar@1: case O_NOTWITHIN: alpar@1: /* test on 'X not within Y' */ alpar@1: { ELEMSET *set; alpar@1: MEMBER *memb; alpar@1: set = eval_elemset(mpl, code->arg.arg.x); alpar@1: value = 1; alpar@1: for (memb = set->head; memb != NULL; memb = memb->next) alpar@1: { if (is_member(mpl, code->arg.arg.y, memb->tuple)) alpar@1: { value = 0; alpar@1: break; alpar@1: } alpar@1: } alpar@1: delete_elemset(mpl, set); alpar@1: } alpar@1: break; alpar@1: case O_FORALL: alpar@1: /* conjunction (A-quantification) */ alpar@1: { struct iter_log_info _info, *info = &_info; alpar@1: info->code = code; alpar@1: info->value = 1; alpar@1: loop_within_domain(mpl, code->arg.loop.domain, info, alpar@1: iter_log_func); alpar@1: value = info->value; alpar@1: } alpar@1: break; alpar@1: case O_EXISTS: alpar@1: /* disjunction (E-quantification) */ alpar@1: { struct iter_log_info _info, *info = &_info; alpar@1: info->code = code; alpar@1: info->value = 0; alpar@1: loop_within_domain(mpl, code->arg.loop.domain, info, alpar@1: iter_log_func); alpar@1: value = info->value; alpar@1: } alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: /* save resultant value */ alpar@1: xassert(!code->valid); alpar@1: code->valid = 1; alpar@1: code->value.bit = value; alpar@1: done: return value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_tuple - evaluate pseudo-code to construct n-tuple. alpar@1: -- alpar@1: -- This routine evaluates specified pseudo-code to construct resultant alpar@1: -- n-tuple, which is returned on exit. */ alpar@1: alpar@1: TUPLE *eval_tuple(MPL *mpl, CODE *code) alpar@1: { TUPLE *value; alpar@1: xassert(code != NULL); alpar@1: xassert(code->type == A_TUPLE); alpar@1: xassert(code->dim > 0); alpar@1: /* if the operation has a side effect, invalidate and delete the alpar@1: resultant value */ alpar@1: if (code->vflag && code->valid) alpar@1: { code->valid = 0; alpar@1: delete_value(mpl, code->type, &code->value); alpar@1: } alpar@1: /* if resultant value is valid, no evaluation is needed */ alpar@1: if (code->valid) alpar@1: { value = copy_tuple(mpl, code->value.tuple); alpar@1: goto done; alpar@1: } alpar@1: /* evaluate pseudo-code recursively */ alpar@1: switch (code->op) alpar@1: { case O_TUPLE: alpar@1: /* make n-tuple */ alpar@1: { ARG_LIST *e; alpar@1: value = create_tuple(mpl); alpar@1: for (e = code->arg.list; e != NULL; e = e->next) alpar@1: value = expand_tuple(mpl, value, eval_symbolic(mpl, alpar@1: e->x)); alpar@1: } alpar@1: break; alpar@1: case O_CVTTUP: alpar@1: /* convert to 1-tuple */ alpar@1: value = expand_tuple(mpl, create_tuple(mpl), alpar@1: eval_symbolic(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: /* save resultant value */ alpar@1: xassert(!code->valid); alpar@1: code->valid = 1; alpar@1: code->value.tuple = copy_tuple(mpl, value); alpar@1: done: return value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_elemset - evaluate pseudo-code to construct elemental set. alpar@1: -- alpar@1: -- This routine evaluates specified pseudo-code to construct resultant alpar@1: -- elemental set, which is returned on exit. */ alpar@1: alpar@1: struct iter_set_info alpar@1: { /* working info used by the routine iter_set_func */ alpar@1: CODE *code; alpar@1: /* pseudo-code for iterated operation to be performed */ alpar@1: ELEMSET *value; alpar@1: /* resultant value */ alpar@1: }; alpar@1: alpar@1: static int iter_set_func(MPL *mpl, void *_info) alpar@1: { /* this is auxiliary routine used to perform iterated operation alpar@1: on n-tuple "integrand" within domain scope */ alpar@1: struct iter_set_info *info = _info; alpar@1: TUPLE *tuple; alpar@1: switch (info->code->op) alpar@1: { case O_SETOF: alpar@1: /* compute next n-tuple and add it to the set; in this case alpar@1: duplicate n-tuples are silently ignored */ alpar@1: tuple = eval_tuple(mpl, info->code->arg.loop.x); alpar@1: if (find_tuple(mpl, info->value, tuple) == NULL) alpar@1: add_tuple(mpl, info->value, tuple); alpar@1: else alpar@1: delete_tuple(mpl, tuple); alpar@1: break; alpar@1: case O_BUILD: alpar@1: /* construct next n-tuple using current values assigned to alpar@1: *free* dummy indices as its components and add it to the alpar@1: set; in this case duplicate n-tuples cannot appear */ alpar@1: add_tuple(mpl, info->value, get_domain_tuple(mpl, alpar@1: info->code->arg.loop.domain)); alpar@1: break; alpar@1: default: alpar@1: xassert(info != info); alpar@1: } alpar@1: return 0; alpar@1: } alpar@1: alpar@1: ELEMSET *eval_elemset(MPL *mpl, CODE *code) alpar@1: { ELEMSET *value; alpar@1: xassert(code != NULL); alpar@1: xassert(code->type == A_ELEMSET); alpar@1: xassert(code->dim > 0); alpar@1: /* if the operation has a side effect, invalidate and delete the alpar@1: resultant value */ alpar@1: if (code->vflag && code->valid) alpar@1: { code->valid = 0; alpar@1: delete_value(mpl, code->type, &code->value); alpar@1: } alpar@1: /* if resultant value is valid, no evaluation is needed */ alpar@1: if (code->valid) alpar@1: { value = copy_elemset(mpl, code->value.set); alpar@1: goto done; alpar@1: } alpar@1: /* evaluate pseudo-code recursively */ alpar@1: switch (code->op) alpar@1: { case O_MEMSET: alpar@1: /* take member of set */ alpar@1: { TUPLE *tuple; alpar@1: ARG_LIST *e; alpar@1: tuple = create_tuple(mpl); alpar@1: for (e = code->arg.set.list; e != NULL; e = e->next) alpar@1: tuple = expand_tuple(mpl, tuple, eval_symbolic(mpl, alpar@1: e->x)); alpar@1: value = copy_elemset(mpl, alpar@1: eval_member_set(mpl, code->arg.set.set, tuple)); alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: break; alpar@1: case O_MAKE: alpar@1: /* make elemental set of n-tuples */ alpar@1: { ARG_LIST *e; alpar@1: value = create_elemset(mpl, code->dim); alpar@1: for (e = code->arg.list; e != NULL; e = e->next) alpar@1: check_then_add(mpl, value, eval_tuple(mpl, e->x)); alpar@1: } alpar@1: break; alpar@1: case O_UNION: alpar@1: /* union of two elemental sets */ alpar@1: value = set_union(mpl, alpar@1: eval_elemset(mpl, code->arg.arg.x), alpar@1: eval_elemset(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_DIFF: alpar@1: /* difference between two elemental sets */ alpar@1: value = set_diff(mpl, alpar@1: eval_elemset(mpl, code->arg.arg.x), alpar@1: eval_elemset(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_SYMDIFF: alpar@1: /* symmetric difference between two elemental sets */ alpar@1: value = set_symdiff(mpl, alpar@1: eval_elemset(mpl, code->arg.arg.x), alpar@1: eval_elemset(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_INTER: alpar@1: /* intersection of two elemental sets */ alpar@1: value = set_inter(mpl, alpar@1: eval_elemset(mpl, code->arg.arg.x), alpar@1: eval_elemset(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_CROSS: alpar@1: /* cross (Cartesian) product of two elemental sets */ alpar@1: value = set_cross(mpl, alpar@1: eval_elemset(mpl, code->arg.arg.x), alpar@1: eval_elemset(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_DOTS: alpar@1: /* build "arithmetic" elemental set */ alpar@1: value = create_arelset(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_numeric(mpl, code->arg.arg.y), alpar@1: code->arg.arg.z == NULL ? 1.0 : eval_numeric(mpl, alpar@1: code->arg.arg.z)); alpar@1: break; alpar@1: case O_FORK: alpar@1: /* if-then-else */ alpar@1: if (eval_logical(mpl, code->arg.arg.x)) alpar@1: value = eval_elemset(mpl, code->arg.arg.y); alpar@1: else alpar@1: value = eval_elemset(mpl, code->arg.arg.z); alpar@1: break; alpar@1: case O_SETOF: alpar@1: /* compute elemental set */ alpar@1: { struct iter_set_info _info, *info = &_info; alpar@1: info->code = code; alpar@1: info->value = create_elemset(mpl, code->dim); alpar@1: loop_within_domain(mpl, code->arg.loop.domain, info, alpar@1: iter_set_func); alpar@1: value = info->value; alpar@1: } alpar@1: break; alpar@1: case O_BUILD: alpar@1: /* build elemental set identical to domain set */ alpar@1: { struct iter_set_info _info, *info = &_info; alpar@1: info->code = code; alpar@1: info->value = create_elemset(mpl, code->dim); alpar@1: loop_within_domain(mpl, code->arg.loop.domain, info, alpar@1: iter_set_func); alpar@1: value = info->value; alpar@1: } alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: /* save resultant value */ alpar@1: xassert(!code->valid); alpar@1: code->valid = 1; alpar@1: code->value.set = copy_elemset(mpl, value); alpar@1: done: return value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- is_member - check if n-tuple is in set specified by pseudo-code. alpar@1: -- alpar@1: -- This routine checks if given n-tuple is a member of elemental set alpar@1: -- specified in the form of pseudo-code (i.e. by expression). alpar@1: -- alpar@1: -- The n-tuple may have more components that dimension of the elemental alpar@1: -- set, in which case the extra components are ignored. */ alpar@1: alpar@1: static void null_func(MPL *mpl, void *info) alpar@1: { /* this is dummy routine used to enter the domain scope */ alpar@1: xassert(mpl == mpl); alpar@1: xassert(info == NULL); alpar@1: return; alpar@1: } alpar@1: alpar@1: int is_member(MPL *mpl, CODE *code, TUPLE *tuple) alpar@1: { int value; alpar@1: xassert(code != NULL); alpar@1: xassert(code->type == A_ELEMSET); alpar@1: xassert(code->dim > 0); alpar@1: xassert(tuple != NULL); alpar@1: switch (code->op) alpar@1: { case O_MEMSET: alpar@1: /* check if given n-tuple is member of elemental set, which alpar@1: is assigned to member of model set */ alpar@1: { ARG_LIST *e; alpar@1: TUPLE *temp; alpar@1: ELEMSET *set; alpar@1: /* evaluate reference to elemental set */ alpar@1: temp = create_tuple(mpl); alpar@1: for (e = code->arg.set.list; e != NULL; e = e->next) alpar@1: temp = expand_tuple(mpl, temp, eval_symbolic(mpl, alpar@1: e->x)); alpar@1: set = eval_member_set(mpl, code->arg.set.set, temp); alpar@1: delete_tuple(mpl, temp); alpar@1: /* check if the n-tuple is contained in the set array */ alpar@1: temp = build_subtuple(mpl, tuple, set->dim); alpar@1: value = (find_tuple(mpl, set, temp) != NULL); alpar@1: delete_tuple(mpl, temp); alpar@1: } alpar@1: break; alpar@1: case O_MAKE: alpar@1: /* check if given n-tuple is member of literal set */ alpar@1: { ARG_LIST *e; alpar@1: TUPLE *temp, *that; alpar@1: value = 0; alpar@1: temp = build_subtuple(mpl, tuple, code->dim); alpar@1: for (e = code->arg.list; e != NULL; e = e->next) alpar@1: { that = eval_tuple(mpl, e->x); alpar@1: value = (compare_tuples(mpl, temp, that) == 0); alpar@1: delete_tuple(mpl, that); alpar@1: if (value) break; alpar@1: } alpar@1: delete_tuple(mpl, temp); alpar@1: } alpar@1: break; alpar@1: case O_UNION: alpar@1: value = is_member(mpl, code->arg.arg.x, tuple) || alpar@1: is_member(mpl, code->arg.arg.y, tuple); alpar@1: break; alpar@1: case O_DIFF: alpar@1: value = is_member(mpl, code->arg.arg.x, tuple) && alpar@1: !is_member(mpl, code->arg.arg.y, tuple); alpar@1: break; alpar@1: case O_SYMDIFF: alpar@1: { int in1 = is_member(mpl, code->arg.arg.x, tuple); alpar@1: int in2 = is_member(mpl, code->arg.arg.y, tuple); alpar@1: value = (in1 && !in2) || (!in1 && in2); alpar@1: } alpar@1: break; alpar@1: case O_INTER: alpar@1: value = is_member(mpl, code->arg.arg.x, tuple) && alpar@1: is_member(mpl, code->arg.arg.y, tuple); alpar@1: break; alpar@1: case O_CROSS: alpar@1: { int j; alpar@1: value = is_member(mpl, code->arg.arg.x, tuple); alpar@1: if (value) alpar@1: { for (j = 1; j <= code->arg.arg.x->dim; j++) alpar@1: { xassert(tuple != NULL); alpar@1: tuple = tuple->next; alpar@1: } alpar@1: value = is_member(mpl, code->arg.arg.y, tuple); alpar@1: } alpar@1: } alpar@1: break; alpar@1: case O_DOTS: alpar@1: /* check if given 1-tuple is member of "arithmetic" set */ alpar@1: { int j; alpar@1: double x, t0, tf, dt; alpar@1: xassert(code->dim == 1); alpar@1: /* compute "parameters" of the "arithmetic" set */ alpar@1: t0 = eval_numeric(mpl, code->arg.arg.x); alpar@1: tf = eval_numeric(mpl, code->arg.arg.y); alpar@1: if (code->arg.arg.z == NULL) alpar@1: dt = 1.0; alpar@1: else alpar@1: dt = eval_numeric(mpl, code->arg.arg.z); alpar@1: /* make sure the parameters are correct */ alpar@1: arelset_size(mpl, t0, tf, dt); alpar@1: /* if component of 1-tuple is symbolic, not numeric, the alpar@1: 1-tuple cannot be member of "arithmetic" set */ alpar@1: xassert(tuple->sym != NULL); alpar@1: if (tuple->sym->str != NULL) alpar@1: { value = 0; alpar@1: break; alpar@1: } alpar@1: /* determine numeric value of the component */ alpar@1: x = tuple->sym->num; alpar@1: /* if the component value is out of the set range, the alpar@1: 1-tuple is not in the set */ alpar@1: if (dt > 0.0 && !(t0 <= x && x <= tf) || alpar@1: dt < 0.0 && !(tf <= x && x <= t0)) alpar@1: { value = 0; alpar@1: break; alpar@1: } alpar@1: /* estimate ordinal number of the 1-tuple in the set */ alpar@1: j = (int)(((x - t0) / dt) + 0.5) + 1; alpar@1: /* perform the main check */ alpar@1: value = (arelset_member(mpl, t0, tf, dt, j) == x); alpar@1: } alpar@1: break; alpar@1: case O_FORK: alpar@1: /* check if given n-tuple is member of conditional set */ alpar@1: if (eval_logical(mpl, code->arg.arg.x)) alpar@1: value = is_member(mpl, code->arg.arg.y, tuple); alpar@1: else alpar@1: value = is_member(mpl, code->arg.arg.z, tuple); alpar@1: break; alpar@1: case O_SETOF: alpar@1: /* check if given n-tuple is member of computed set */ alpar@1: /* it is not clear how to efficiently perform the check not alpar@1: computing the entire elemental set :+( */ alpar@1: error(mpl, "implementation restriction; in/within setof{} n" alpar@1: "ot allowed"); alpar@1: break; alpar@1: case O_BUILD: alpar@1: /* check if given n-tuple is member of domain set */ alpar@1: { TUPLE *temp; alpar@1: temp = build_subtuple(mpl, tuple, code->dim); alpar@1: /* try to enter the domain scope; if it is successful, alpar@1: the n-tuple is in the domain set */ alpar@1: value = (eval_within_domain(mpl, code->arg.loop.domain, alpar@1: temp, NULL, null_func) == 0); alpar@1: delete_tuple(mpl, temp); alpar@1: } alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: return value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- eval_formula - evaluate pseudo-code to construct linear form. alpar@1: -- alpar@1: -- This routine evaluates specified pseudo-code to construct resultant alpar@1: -- linear form, which is returned on exit. */ alpar@1: alpar@1: struct iter_form_info alpar@1: { /* working info used by the routine iter_form_func */ alpar@1: CODE *code; alpar@1: /* pseudo-code for iterated operation to be performed */ alpar@1: FORMULA *value; alpar@1: /* resultant value */ alpar@1: FORMULA *tail; alpar@1: /* pointer to the last term */ alpar@1: }; alpar@1: alpar@1: static int iter_form_func(MPL *mpl, void *_info) alpar@1: { /* this is auxiliary routine used to perform iterated operation alpar@1: on linear form "integrand" within domain scope */ alpar@1: struct iter_form_info *info = _info; alpar@1: switch (info->code->op) alpar@1: { case O_SUM: alpar@1: /* summation over domain */ alpar@1: #if 0 alpar@1: info->value = alpar@1: linear_comb(mpl, alpar@1: +1.0, info->value, alpar@1: +1.0, eval_formula(mpl, info->code->arg.loop.x)); alpar@1: #else alpar@1: /* the routine linear_comb needs to look through all terms alpar@1: of both linear forms to reduce identical terms, so using alpar@1: it here is not a good idea (for example, evaluation of alpar@1: sum{i in 1..n} x[i] required quadratic time); the better alpar@1: idea is to gather all terms of the integrand in one list alpar@1: and reduce identical terms only once after all terms of alpar@1: the resultant linear form have been evaluated */ alpar@1: { FORMULA *form, *term; alpar@1: form = eval_formula(mpl, info->code->arg.loop.x); alpar@1: if (info->value == NULL) alpar@1: { xassert(info->tail == NULL); alpar@1: info->value = form; alpar@1: } alpar@1: else alpar@1: { xassert(info->tail != NULL); alpar@1: info->tail->next = form; alpar@1: } alpar@1: for (term = form; term != NULL; term = term->next) alpar@1: info->tail = term; alpar@1: } alpar@1: #endif alpar@1: break; alpar@1: default: alpar@1: xassert(info != info); alpar@1: } alpar@1: return 0; alpar@1: } alpar@1: alpar@1: FORMULA *eval_formula(MPL *mpl, CODE *code) alpar@1: { FORMULA *value; alpar@1: xassert(code != NULL); alpar@1: xassert(code->type == A_FORMULA); alpar@1: xassert(code->dim == 0); alpar@1: /* if the operation has a side effect, invalidate and delete the alpar@1: resultant value */ alpar@1: if (code->vflag && code->valid) alpar@1: { code->valid = 0; alpar@1: delete_value(mpl, code->type, &code->value); alpar@1: } alpar@1: /* if resultant value is valid, no evaluation is needed */ alpar@1: if (code->valid) alpar@1: { value = copy_formula(mpl, code->value.form); alpar@1: goto done; alpar@1: } alpar@1: /* evaluate pseudo-code recursively */ alpar@1: switch (code->op) alpar@1: { case O_MEMVAR: alpar@1: /* take member of variable */ alpar@1: { TUPLE *tuple; alpar@1: ARG_LIST *e; alpar@1: tuple = create_tuple(mpl); alpar@1: for (e = code->arg.var.list; e != NULL; e = e->next) alpar@1: tuple = expand_tuple(mpl, tuple, eval_symbolic(mpl, alpar@1: e->x)); alpar@1: #if 1 /* 15/V-2010 */ alpar@1: xassert(code->arg.var.suff == DOT_NONE); alpar@1: #endif alpar@1: value = single_variable(mpl, alpar@1: eval_member_var(mpl, code->arg.var.var, tuple)); alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: break; alpar@1: case O_CVTLFM: alpar@1: /* convert to linear form */ alpar@1: value = constant_term(mpl, eval_numeric(mpl, alpar@1: code->arg.arg.x)); alpar@1: break; alpar@1: case O_PLUS: alpar@1: /* unary plus */ alpar@1: value = linear_comb(mpl, alpar@1: 0.0, constant_term(mpl, 0.0), alpar@1: +1.0, eval_formula(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_MINUS: alpar@1: /* unary minus */ alpar@1: value = linear_comb(mpl, alpar@1: 0.0, constant_term(mpl, 0.0), alpar@1: -1.0, eval_formula(mpl, code->arg.arg.x)); alpar@1: break; alpar@1: case O_ADD: alpar@1: /* addition */ alpar@1: value = linear_comb(mpl, alpar@1: +1.0, eval_formula(mpl, code->arg.arg.x), alpar@1: +1.0, eval_formula(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_SUB: alpar@1: /* subtraction */ alpar@1: value = linear_comb(mpl, alpar@1: +1.0, eval_formula(mpl, code->arg.arg.x), alpar@1: -1.0, eval_formula(mpl, code->arg.arg.y)); alpar@1: break; alpar@1: case O_MUL: alpar@1: /* multiplication */ alpar@1: xassert(code->arg.arg.x != NULL); alpar@1: xassert(code->arg.arg.y != NULL); alpar@1: if (code->arg.arg.x->type == A_NUMERIC) alpar@1: { xassert(code->arg.arg.y->type == A_FORMULA); alpar@1: value = linear_comb(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.x), alpar@1: eval_formula(mpl, code->arg.arg.y), alpar@1: 0.0, constant_term(mpl, 0.0)); alpar@1: } alpar@1: else alpar@1: { xassert(code->arg.arg.x->type == A_FORMULA); alpar@1: xassert(code->arg.arg.y->type == A_NUMERIC); alpar@1: value = linear_comb(mpl, alpar@1: eval_numeric(mpl, code->arg.arg.y), alpar@1: eval_formula(mpl, code->arg.arg.x), alpar@1: 0.0, constant_term(mpl, 0.0)); alpar@1: } alpar@1: break; alpar@1: case O_DIV: alpar@1: /* division */ alpar@1: value = linear_comb(mpl, alpar@1: fp_div(mpl, 1.0, eval_numeric(mpl, code->arg.arg.y)), alpar@1: eval_formula(mpl, code->arg.arg.x), alpar@1: 0.0, constant_term(mpl, 0.0)); alpar@1: break; alpar@1: case O_FORK: alpar@1: /* if-then-else */ alpar@1: if (eval_logical(mpl, code->arg.arg.x)) alpar@1: value = eval_formula(mpl, code->arg.arg.y); alpar@1: else if (code->arg.arg.z == NULL) alpar@1: value = constant_term(mpl, 0.0); alpar@1: else alpar@1: value = eval_formula(mpl, code->arg.arg.z); alpar@1: break; alpar@1: case O_SUM: alpar@1: /* summation over domain */ alpar@1: { struct iter_form_info _info, *info = &_info; alpar@1: info->code = code; alpar@1: info->value = constant_term(mpl, 0.0); alpar@1: info->tail = NULL; alpar@1: loop_within_domain(mpl, code->arg.loop.domain, info, alpar@1: iter_form_func); alpar@1: value = reduce_terms(mpl, info->value); alpar@1: } alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: /* save resultant value */ alpar@1: xassert(!code->valid); alpar@1: code->valid = 1; alpar@1: code->value.form = copy_formula(mpl, value); alpar@1: done: return value; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean_code - clean pseudo-code. alpar@1: -- alpar@1: -- This routine recursively cleans specified pseudo-code that assumes alpar@1: -- deleting all temporary resultant values. */ alpar@1: alpar@1: void clean_code(MPL *mpl, CODE *code) alpar@1: { ARG_LIST *e; alpar@1: /* if no pseudo-code is specified, do nothing */ alpar@1: if (code == NULL) goto done; alpar@1: /* if resultant value is valid (exists), delete it */ alpar@1: if (code->valid) alpar@1: { code->valid = 0; alpar@1: delete_value(mpl, code->type, &code->value); alpar@1: } alpar@1: /* recursively clean pseudo-code for operands */ alpar@1: switch (code->op) alpar@1: { case O_NUMBER: alpar@1: case O_STRING: alpar@1: case O_INDEX: alpar@1: break; alpar@1: case O_MEMNUM: alpar@1: case O_MEMSYM: alpar@1: for (e = code->arg.par.list; e != NULL; e = e->next) alpar@1: clean_code(mpl, e->x); alpar@1: break; alpar@1: case O_MEMSET: alpar@1: for (e = code->arg.set.list; e != NULL; e = e->next) alpar@1: clean_code(mpl, e->x); alpar@1: break; alpar@1: case O_MEMVAR: alpar@1: for (e = code->arg.var.list; e != NULL; e = e->next) alpar@1: clean_code(mpl, e->x); alpar@1: break; alpar@1: #if 1 /* 15/V-2010 */ alpar@1: case O_MEMCON: alpar@1: for (e = code->arg.con.list; e != NULL; e = e->next) alpar@1: clean_code(mpl, e->x); alpar@1: break; alpar@1: #endif alpar@1: case O_TUPLE: alpar@1: case O_MAKE: alpar@1: for (e = code->arg.list; e != NULL; e = e->next) alpar@1: clean_code(mpl, e->x); alpar@1: break; alpar@1: case O_SLICE: alpar@1: xassert(code != code); alpar@1: case O_IRAND224: alpar@1: case O_UNIFORM01: alpar@1: case O_NORMAL01: alpar@1: case O_GMTIME: alpar@1: break; alpar@1: case O_CVTNUM: alpar@1: case O_CVTSYM: alpar@1: case O_CVTLOG: alpar@1: case O_CVTTUP: alpar@1: case O_CVTLFM: alpar@1: case O_PLUS: alpar@1: case O_MINUS: alpar@1: case O_NOT: alpar@1: case O_ABS: alpar@1: case O_CEIL: alpar@1: case O_FLOOR: alpar@1: case O_EXP: alpar@1: case O_LOG: alpar@1: case O_LOG10: alpar@1: case O_SQRT: alpar@1: case O_SIN: alpar@1: case O_COS: alpar@1: case O_ATAN: alpar@1: case O_ROUND: alpar@1: case O_TRUNC: alpar@1: case O_CARD: alpar@1: case O_LENGTH: alpar@1: /* unary operation */ alpar@1: clean_code(mpl, code->arg.arg.x); alpar@1: break; alpar@1: case O_ADD: alpar@1: case O_SUB: alpar@1: case O_LESS: alpar@1: case O_MUL: alpar@1: case O_DIV: alpar@1: case O_IDIV: alpar@1: case O_MOD: alpar@1: case O_POWER: alpar@1: case O_ATAN2: alpar@1: case O_ROUND2: alpar@1: case O_TRUNC2: alpar@1: case O_UNIFORM: alpar@1: case O_NORMAL: alpar@1: case O_CONCAT: alpar@1: case O_LT: alpar@1: case O_LE: alpar@1: case O_EQ: alpar@1: case O_GE: alpar@1: case O_GT: alpar@1: case O_NE: alpar@1: case O_AND: alpar@1: case O_OR: alpar@1: case O_UNION: alpar@1: case O_DIFF: alpar@1: case O_SYMDIFF: alpar@1: case O_INTER: alpar@1: case O_CROSS: alpar@1: case O_IN: alpar@1: case O_NOTIN: alpar@1: case O_WITHIN: alpar@1: case O_NOTWITHIN: alpar@1: case O_SUBSTR: alpar@1: case O_STR2TIME: alpar@1: case O_TIME2STR: alpar@1: /* binary operation */ alpar@1: clean_code(mpl, code->arg.arg.x); alpar@1: clean_code(mpl, code->arg.arg.y); alpar@1: break; alpar@1: case O_DOTS: alpar@1: case O_FORK: alpar@1: case O_SUBSTR3: alpar@1: /* ternary operation */ alpar@1: clean_code(mpl, code->arg.arg.x); alpar@1: clean_code(mpl, code->arg.arg.y); alpar@1: clean_code(mpl, code->arg.arg.z); alpar@1: break; alpar@1: case O_MIN: alpar@1: case O_MAX: alpar@1: /* n-ary operation */ alpar@1: for (e = code->arg.list; e != NULL; e = e->next) alpar@1: clean_code(mpl, e->x); alpar@1: break; alpar@1: case O_SUM: alpar@1: case O_PROD: alpar@1: case O_MINIMUM: alpar@1: case O_MAXIMUM: alpar@1: case O_FORALL: alpar@1: case O_EXISTS: alpar@1: case O_SETOF: alpar@1: case O_BUILD: alpar@1: /* iterated operation */ alpar@1: clean_domain(mpl, code->arg.loop.domain); alpar@1: clean_code(mpl, code->arg.loop.x); alpar@1: break; alpar@1: default: alpar@1: xassert(code->op != code->op); alpar@1: } alpar@1: done: return; alpar@1: } alpar@1: alpar@1: #if 1 /* 11/II-2008 */ alpar@1: /**********************************************************************/ alpar@1: /* * * DATA TABLES * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: int mpl_tab_num_args(TABDCA *dca) alpar@1: { /* returns the number of arguments */ alpar@1: return dca->na; alpar@1: } alpar@1: alpar@1: const char *mpl_tab_get_arg(TABDCA *dca, int k) alpar@1: { /* returns pointer to k-th argument */ alpar@1: xassert(1 <= k && k <= dca->na); alpar@1: return dca->arg[k]; alpar@1: } alpar@1: alpar@1: int mpl_tab_num_flds(TABDCA *dca) alpar@1: { /* returns the number of fields */ alpar@1: return dca->nf; alpar@1: } alpar@1: alpar@1: const char *mpl_tab_get_name(TABDCA *dca, int k) alpar@1: { /* returns pointer to name of k-th field */ alpar@1: xassert(1 <= k && k <= dca->nf); alpar@1: return dca->name[k]; alpar@1: } alpar@1: alpar@1: int mpl_tab_get_type(TABDCA *dca, int k) alpar@1: { /* returns type of k-th field */ alpar@1: xassert(1 <= k && k <= dca->nf); alpar@1: return dca->type[k]; alpar@1: } alpar@1: alpar@1: double mpl_tab_get_num(TABDCA *dca, int k) alpar@1: { /* returns numeric value of k-th field */ alpar@1: xassert(1 <= k && k <= dca->nf); alpar@1: xassert(dca->type[k] == 'N'); alpar@1: return dca->num[k]; alpar@1: } alpar@1: alpar@1: const char *mpl_tab_get_str(TABDCA *dca, int k) alpar@1: { /* returns pointer to string value of k-th field */ alpar@1: xassert(1 <= k && k <= dca->nf); alpar@1: xassert(dca->type[k] == 'S'); alpar@1: xassert(dca->str[k] != NULL); alpar@1: return dca->str[k]; alpar@1: } alpar@1: alpar@1: void mpl_tab_set_num(TABDCA *dca, int k, double num) alpar@1: { /* assign numeric value to k-th field */ alpar@1: xassert(1 <= k && k <= dca->nf); alpar@1: xassert(dca->type[k] == '?'); alpar@1: dca->type[k] = 'N'; alpar@1: dca->num[k] = num; alpar@1: return; alpar@1: } alpar@1: alpar@1: void mpl_tab_set_str(TABDCA *dca, int k, const char *str) alpar@1: { /* assign string value to k-th field */ alpar@1: xassert(1 <= k && k <= dca->nf); alpar@1: xassert(dca->type[k] == '?'); alpar@1: xassert(strlen(str) <= MAX_LENGTH); alpar@1: xassert(dca->str[k] != NULL); alpar@1: dca->type[k] = 'S'; alpar@1: strcpy(dca->str[k], str); alpar@1: return; alpar@1: } alpar@1: alpar@1: static int write_func(MPL *mpl, void *info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: TABLE *tab = info; alpar@1: TABDCA *dca = mpl->dca; alpar@1: TABOUT *out; alpar@1: SYMBOL *sym; alpar@1: int k; alpar@1: char buf[MAX_LENGTH+1]; alpar@1: /* evaluate field values */ alpar@1: k = 0; alpar@1: for (out = tab->u.out.list; out != NULL; out = out->next) alpar@1: { k++; alpar@1: switch (out->code->type) alpar@1: { case A_NUMERIC: alpar@1: dca->type[k] = 'N'; alpar@1: dca->num[k] = eval_numeric(mpl, out->code); alpar@1: dca->str[k][0] = '\0'; alpar@1: break; alpar@1: case A_SYMBOLIC: alpar@1: sym = eval_symbolic(mpl, out->code); alpar@1: if (sym->str == NULL) alpar@1: { dca->type[k] = 'N'; alpar@1: dca->num[k] = sym->num; alpar@1: dca->str[k][0] = '\0'; alpar@1: } alpar@1: else alpar@1: { dca->type[k] = 'S'; alpar@1: dca->num[k] = 0.0; alpar@1: fetch_string(mpl, sym->str, buf); alpar@1: strcpy(dca->str[k], buf); alpar@1: } alpar@1: delete_symbol(mpl, sym); alpar@1: break; alpar@1: default: alpar@1: xassert(out != out); alpar@1: } alpar@1: } alpar@1: /* write record to output table */ alpar@1: mpl_tab_drv_write(mpl); alpar@1: return 0; alpar@1: } alpar@1: alpar@1: void execute_table(MPL *mpl, TABLE *tab) alpar@1: { /* execute table statement */ alpar@1: TABARG *arg; alpar@1: TABFLD *fld; alpar@1: TABIN *in; alpar@1: TABOUT *out; alpar@1: TABDCA *dca; alpar@1: SET *set; alpar@1: int k; alpar@1: char buf[MAX_LENGTH+1]; alpar@1: /* allocate table driver communication area */ alpar@1: xassert(mpl->dca == NULL); alpar@1: mpl->dca = dca = xmalloc(sizeof(TABDCA)); alpar@1: dca->id = 0; alpar@1: dca->link = NULL; alpar@1: dca->na = 0; alpar@1: dca->arg = NULL; alpar@1: dca->nf = 0; alpar@1: dca->name = NULL; alpar@1: dca->type = NULL; alpar@1: dca->num = NULL; alpar@1: dca->str = NULL; alpar@1: /* allocate arguments */ alpar@1: xassert(dca->na == 0); alpar@1: for (arg = tab->arg; arg != NULL; arg = arg->next) alpar@1: dca->na++; alpar@1: dca->arg = xcalloc(1+dca->na, sizeof(char *)); alpar@1: #if 1 /* 28/IX-2008 */ alpar@1: for (k = 1; k <= dca->na; k++) dca->arg[k] = NULL; alpar@1: #endif alpar@1: /* evaluate argument values */ alpar@1: k = 0; alpar@1: for (arg = tab->arg; arg != NULL; arg = arg->next) alpar@1: { SYMBOL *sym; alpar@1: k++; alpar@1: xassert(arg->code->type == A_SYMBOLIC); alpar@1: sym = eval_symbolic(mpl, arg->code); alpar@1: if (sym->str == NULL) alpar@1: sprintf(buf, "%.*g", DBL_DIG, sym->num); alpar@1: else alpar@1: fetch_string(mpl, sym->str, buf); alpar@1: delete_symbol(mpl, sym); alpar@1: dca->arg[k] = xmalloc(strlen(buf)+1); alpar@1: strcpy(dca->arg[k], buf); alpar@1: } alpar@1: /* perform table input/output */ alpar@1: switch (tab->type) alpar@1: { case A_INPUT: goto read_table; alpar@1: case A_OUTPUT: goto write_table; alpar@1: default: xassert(tab != tab); alpar@1: } alpar@1: read_table: alpar@1: /* read data from input table */ alpar@1: /* add the only member to the control set and assign it empty alpar@1: elemental set */ alpar@1: set = tab->u.in.set; alpar@1: if (set != NULL) alpar@1: { if (set->data) alpar@1: error(mpl, "%s already provided with data", set->name); alpar@1: xassert(set->array->head == NULL); alpar@1: add_member(mpl, set->array, NULL)->value.set = alpar@1: create_elemset(mpl, set->dimen); alpar@1: set->data = 1; alpar@1: } alpar@1: /* check parameters specified in the input list */ alpar@1: for (in = tab->u.in.list; in != NULL; in = in->next) alpar@1: { if (in->par->data) alpar@1: error(mpl, "%s already provided with data", in->par->name); alpar@1: in->par->data = 1; alpar@1: } alpar@1: /* allocate and initialize fields */ alpar@1: xassert(dca->nf == 0); alpar@1: for (fld = tab->u.in.fld; fld != NULL; fld = fld->next) alpar@1: dca->nf++; alpar@1: for (in = tab->u.in.list; in != NULL; in = in->next) alpar@1: dca->nf++; alpar@1: dca->name = xcalloc(1+dca->nf, sizeof(char *)); alpar@1: dca->type = xcalloc(1+dca->nf, sizeof(int)); alpar@1: dca->num = xcalloc(1+dca->nf, sizeof(double)); alpar@1: dca->str = xcalloc(1+dca->nf, sizeof(char *)); alpar@1: k = 0; alpar@1: for (fld = tab->u.in.fld; fld != NULL; fld = fld->next) alpar@1: { k++; alpar@1: dca->name[k] = fld->name; alpar@1: dca->type[k] = '?'; alpar@1: dca->num[k] = 0.0; alpar@1: dca->str[k] = xmalloc(MAX_LENGTH+1); alpar@1: dca->str[k][0] = '\0'; alpar@1: } alpar@1: for (in = tab->u.in.list; in != NULL; in = in->next) alpar@1: { k++; alpar@1: dca->name[k] = in->name; alpar@1: dca->type[k] = '?'; alpar@1: dca->num[k] = 0.0; alpar@1: dca->str[k] = xmalloc(MAX_LENGTH+1); alpar@1: dca->str[k][0] = '\0'; alpar@1: } alpar@1: /* open input table */ alpar@1: mpl_tab_drv_open(mpl, 'R'); alpar@1: /* read and process records */ alpar@1: for (;;) alpar@1: { TUPLE *tup; alpar@1: /* reset field types */ alpar@1: for (k = 1; k <= dca->nf; k++) alpar@1: dca->type[k] = '?'; alpar@1: /* read next record */ alpar@1: if (mpl_tab_drv_read(mpl)) break; alpar@1: /* all fields must be set by the driver */ alpar@1: for (k = 1; k <= dca->nf; k++) alpar@1: { if (dca->type[k] == '?') alpar@1: error(mpl, "field %s missing in input table", alpar@1: dca->name[k]); alpar@1: } alpar@1: /* construct n-tuple */ alpar@1: tup = create_tuple(mpl); alpar@1: k = 0; alpar@1: for (fld = tab->u.in.fld; fld != NULL; fld = fld->next) alpar@1: { k++; alpar@1: xassert(k <= dca->nf); alpar@1: switch (dca->type[k]) alpar@1: { case 'N': alpar@1: tup = expand_tuple(mpl, tup, create_symbol_num(mpl, alpar@1: dca->num[k])); alpar@1: break; alpar@1: case 'S': alpar@1: xassert(strlen(dca->str[k]) <= MAX_LENGTH); alpar@1: tup = expand_tuple(mpl, tup, create_symbol_str(mpl, alpar@1: create_string(mpl, dca->str[k]))); alpar@1: break; alpar@1: default: alpar@1: xassert(dca != dca); alpar@1: } alpar@1: } alpar@1: /* add n-tuple just read to the control set */ alpar@1: if (tab->u.in.set != NULL) alpar@1: check_then_add(mpl, tab->u.in.set->array->head->value.set, alpar@1: copy_tuple(mpl, tup)); alpar@1: /* assign values to the parameters in the input list */ alpar@1: for (in = tab->u.in.list; in != NULL; in = in->next) alpar@1: { MEMBER *memb; alpar@1: k++; alpar@1: xassert(k <= dca->nf); alpar@1: /* there must be no member with the same n-tuple */ alpar@1: if (find_member(mpl, in->par->array, tup) != NULL) alpar@1: error(mpl, "%s%s already defined", in->par->name, alpar@1: format_tuple(mpl, '[', tup)); alpar@1: /* create new parameter member with given n-tuple */ alpar@1: memb = add_member(mpl, in->par->array, copy_tuple(mpl, tup)) alpar@1: ; alpar@1: /* assign value to the parameter member */ alpar@1: switch (in->par->type) alpar@1: { case A_NUMERIC: alpar@1: case A_INTEGER: alpar@1: case A_BINARY: alpar@1: if (dca->type[k] != 'N') alpar@1: error(mpl, "%s requires numeric data", alpar@1: in->par->name); alpar@1: memb->value.num = dca->num[k]; alpar@1: break; alpar@1: case A_SYMBOLIC: alpar@1: switch (dca->type[k]) alpar@1: { case 'N': alpar@1: memb->value.sym = create_symbol_num(mpl, alpar@1: dca->num[k]); alpar@1: break; alpar@1: case 'S': alpar@1: xassert(strlen(dca->str[k]) <= MAX_LENGTH); alpar@1: memb->value.sym = create_symbol_str(mpl, alpar@1: create_string(mpl,dca->str[k])); alpar@1: break; alpar@1: default: alpar@1: xassert(dca != dca); alpar@1: } alpar@1: break; alpar@1: default: alpar@1: xassert(in != in); alpar@1: } alpar@1: } alpar@1: /* n-tuple is no more needed */ alpar@1: delete_tuple(mpl, tup); alpar@1: } alpar@1: /* close input table */ alpar@1: mpl_tab_drv_close(mpl); alpar@1: goto done; alpar@1: write_table: alpar@1: /* write data to output table */ alpar@1: /* allocate and initialize fields */ alpar@1: xassert(dca->nf == 0); alpar@1: for (out = tab->u.out.list; out != NULL; out = out->next) alpar@1: dca->nf++; alpar@1: dca->name = xcalloc(1+dca->nf, sizeof(char *)); alpar@1: dca->type = xcalloc(1+dca->nf, sizeof(int)); alpar@1: dca->num = xcalloc(1+dca->nf, sizeof(double)); alpar@1: dca->str = xcalloc(1+dca->nf, sizeof(char *)); alpar@1: k = 0; alpar@1: for (out = tab->u.out.list; out != NULL; out = out->next) alpar@1: { k++; alpar@1: dca->name[k] = out->name; alpar@1: dca->type[k] = '?'; alpar@1: dca->num[k] = 0.0; alpar@1: dca->str[k] = xmalloc(MAX_LENGTH+1); alpar@1: dca->str[k][0] = '\0'; alpar@1: } alpar@1: /* open output table */ alpar@1: mpl_tab_drv_open(mpl, 'W'); alpar@1: /* evaluate fields and write records */ alpar@1: loop_within_domain(mpl, tab->u.out.domain, tab, write_func); alpar@1: /* close output table */ alpar@1: mpl_tab_drv_close(mpl); alpar@1: done: /* free table driver communication area */ alpar@1: free_dca(mpl); alpar@1: return; alpar@1: } alpar@1: alpar@1: void free_dca(MPL *mpl) alpar@1: { /* free table driver communucation area */ alpar@1: TABDCA *dca = mpl->dca; alpar@1: int k; alpar@1: if (dca != NULL) alpar@1: { if (dca->link != NULL) alpar@1: mpl_tab_drv_close(mpl); alpar@1: if (dca->arg != NULL) alpar@1: { for (k = 1; k <= dca->na; k++) alpar@1: #if 1 /* 28/IX-2008 */ alpar@1: if (dca->arg[k] != NULL) alpar@1: #endif alpar@1: xfree(dca->arg[k]); alpar@1: xfree(dca->arg); alpar@1: } alpar@1: if (dca->name != NULL) xfree(dca->name); alpar@1: if (dca->type != NULL) xfree(dca->type); alpar@1: if (dca->num != NULL) xfree(dca->num); alpar@1: if (dca->str != NULL) alpar@1: { for (k = 1; k <= dca->nf; k++) alpar@1: xfree(dca->str[k]); alpar@1: xfree(dca->str); alpar@1: } alpar@1: xfree(dca), mpl->dca = NULL; alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: void clean_table(MPL *mpl, TABLE *tab) alpar@1: { /* clean table statement */ alpar@1: TABARG *arg; alpar@1: TABOUT *out; alpar@1: /* clean string list */ alpar@1: for (arg = tab->arg; arg != NULL; arg = arg->next) alpar@1: clean_code(mpl, arg->code); alpar@1: switch (tab->type) alpar@1: { case A_INPUT: alpar@1: break; alpar@1: case A_OUTPUT: alpar@1: /* clean subscript domain */ alpar@1: clean_domain(mpl, tab->u.out.domain); alpar@1: /* clean output list */ alpar@1: for (out = tab->u.out.list; out != NULL; out = out->next) alpar@1: clean_code(mpl, out->code); alpar@1: break; alpar@1: default: alpar@1: xassert(tab != tab); alpar@1: } alpar@1: return; alpar@1: } alpar@1: #endif alpar@1: alpar@1: /**********************************************************************/ alpar@1: /* * * MODEL STATEMENTS * * */ alpar@1: /**********************************************************************/ alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- execute_check - execute check statement. alpar@1: -- alpar@1: -- This routine executes specified check statement. */ alpar@1: alpar@1: static int check_func(MPL *mpl, void *info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: CHECK *chk = (CHECK *)info; alpar@1: if (!eval_logical(mpl, chk->code)) alpar@1: error(mpl, "check%s failed", format_tuple(mpl, '[', alpar@1: get_domain_tuple(mpl, chk->domain))); alpar@1: return 0; alpar@1: } alpar@1: alpar@1: void execute_check(MPL *mpl, CHECK *chk) alpar@1: { loop_within_domain(mpl, chk->domain, chk, check_func); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean_check - clean check statement. alpar@1: -- alpar@1: -- This routine cleans specified check statement that assumes deleting alpar@1: -- all stuff dynamically allocated on generating/postsolving phase. */ alpar@1: alpar@1: void clean_check(MPL *mpl, CHECK *chk) alpar@1: { /* clean subscript domain */ alpar@1: clean_domain(mpl, chk->domain); alpar@1: /* clean pseudo-code for computing predicate */ alpar@1: clean_code(mpl, chk->code); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- execute_display - execute display statement. alpar@1: -- alpar@1: -- This routine executes specified display statement. */ alpar@1: alpar@1: static void display_set(MPL *mpl, SET *set, MEMBER *memb) alpar@1: { /* display member of model set */ alpar@1: ELEMSET *s = memb->value.set; alpar@1: MEMBER *m; alpar@1: write_text(mpl, "%s%s%s\n", set->name, alpar@1: format_tuple(mpl, '[', memb->tuple), alpar@1: s->head == NULL ? " is empty" : ":"); alpar@1: for (m = s->head; m != NULL; m = m->next) alpar@1: write_text(mpl, " %s\n", format_tuple(mpl, '(', m->tuple)); alpar@1: return; alpar@1: } alpar@1: alpar@1: static void display_par(MPL *mpl, PARAMETER *par, MEMBER *memb) alpar@1: { /* display member of model parameter */ alpar@1: switch (par->type) alpar@1: { case A_NUMERIC: alpar@1: case A_INTEGER: alpar@1: case A_BINARY: alpar@1: write_text(mpl, "%s%s = %.*g\n", par->name, alpar@1: format_tuple(mpl, '[', memb->tuple), alpar@1: DBL_DIG, memb->value.num); alpar@1: break; alpar@1: case A_SYMBOLIC: alpar@1: write_text(mpl, "%s%s = %s\n", par->name, alpar@1: format_tuple(mpl, '[', memb->tuple), alpar@1: format_symbol(mpl, memb->value.sym)); alpar@1: break; alpar@1: default: alpar@1: xassert(par != par); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: #if 1 /* 15/V-2010 */ alpar@1: static void display_var(MPL *mpl, VARIABLE *var, MEMBER *memb, alpar@1: int suff) alpar@1: { /* display member of model variable */ alpar@1: if (suff == DOT_NONE || suff == DOT_VAL) alpar@1: write_text(mpl, "%s%s.val = %.*g\n", var->name, alpar@1: format_tuple(mpl, '[', memb->tuple), DBL_DIG, alpar@1: memb->value.var->prim); alpar@1: else if (suff == DOT_LB) alpar@1: write_text(mpl, "%s%s.lb = %.*g\n", var->name, alpar@1: format_tuple(mpl, '[', memb->tuple), DBL_DIG, alpar@1: memb->value.var->var->lbnd == NULL ? -DBL_MAX : alpar@1: memb->value.var->lbnd); alpar@1: else if (suff == DOT_UB) alpar@1: write_text(mpl, "%s%s.ub = %.*g\n", var->name, alpar@1: format_tuple(mpl, '[', memb->tuple), DBL_DIG, alpar@1: memb->value.var->var->ubnd == NULL ? +DBL_MAX : alpar@1: memb->value.var->ubnd); alpar@1: else if (suff == DOT_STATUS) alpar@1: write_text(mpl, "%s%s.status = %d\n", var->name, format_tuple alpar@1: (mpl, '[', memb->tuple), memb->value.var->stat); alpar@1: else if (suff == DOT_DUAL) alpar@1: write_text(mpl, "%s%s.dual = %.*g\n", var->name, alpar@1: format_tuple(mpl, '[', memb->tuple), DBL_DIG, alpar@1: memb->value.var->dual); alpar@1: else alpar@1: xassert(suff != suff); alpar@1: return; alpar@1: } alpar@1: #endif alpar@1: alpar@1: #if 1 /* 15/V-2010 */ alpar@1: static void display_con(MPL *mpl, CONSTRAINT *con, MEMBER *memb, alpar@1: int suff) alpar@1: { /* display member of model constraint */ alpar@1: if (suff == DOT_NONE || suff == DOT_VAL) alpar@1: write_text(mpl, "%s%s.val = %.*g\n", con->name, alpar@1: format_tuple(mpl, '[', memb->tuple), DBL_DIG, alpar@1: memb->value.con->prim); alpar@1: else if (suff == DOT_LB) alpar@1: write_text(mpl, "%s%s.lb = %.*g\n", con->name, alpar@1: format_tuple(mpl, '[', memb->tuple), DBL_DIG, alpar@1: memb->value.con->con->lbnd == NULL ? -DBL_MAX : alpar@1: memb->value.con->lbnd); alpar@1: else if (suff == DOT_UB) alpar@1: write_text(mpl, "%s%s.ub = %.*g\n", con->name, alpar@1: format_tuple(mpl, '[', memb->tuple), DBL_DIG, alpar@1: memb->value.con->con->ubnd == NULL ? +DBL_MAX : alpar@1: memb->value.con->ubnd); alpar@1: else if (suff == DOT_STATUS) alpar@1: write_text(mpl, "%s%s.status = %d\n", con->name, format_tuple alpar@1: (mpl, '[', memb->tuple), memb->value.con->stat); alpar@1: else if (suff == DOT_DUAL) alpar@1: write_text(mpl, "%s%s.dual = %.*g\n", con->name, alpar@1: format_tuple(mpl, '[', memb->tuple), DBL_DIG, alpar@1: memb->value.con->dual); alpar@1: else alpar@1: xassert(suff != suff); alpar@1: return; alpar@1: } alpar@1: #endif alpar@1: alpar@1: static void display_memb(MPL *mpl, CODE *code) alpar@1: { /* display member specified by pseudo-code */ alpar@1: MEMBER memb; alpar@1: ARG_LIST *e; alpar@1: xassert(code->op == O_MEMNUM || code->op == O_MEMSYM alpar@1: || code->op == O_MEMSET || code->op == O_MEMVAR alpar@1: || code->op == O_MEMCON); alpar@1: memb.tuple = create_tuple(mpl); alpar@1: for (e = code->arg.par.list; e != NULL; e = e->next) alpar@1: memb.tuple = expand_tuple(mpl, memb.tuple, eval_symbolic(mpl, alpar@1: e->x)); alpar@1: switch (code->op) alpar@1: { case O_MEMNUM: alpar@1: memb.value.num = eval_member_num(mpl, code->arg.par.par, alpar@1: memb.tuple); alpar@1: display_par(mpl, code->arg.par.par, &memb); alpar@1: break; alpar@1: case O_MEMSYM: alpar@1: memb.value.sym = eval_member_sym(mpl, code->arg.par.par, alpar@1: memb.tuple); alpar@1: display_par(mpl, code->arg.par.par, &memb); alpar@1: delete_symbol(mpl, memb.value.sym); alpar@1: break; alpar@1: case O_MEMSET: alpar@1: memb.value.set = eval_member_set(mpl, code->arg.set.set, alpar@1: memb.tuple); alpar@1: display_set(mpl, code->arg.set.set, &memb); alpar@1: break; alpar@1: case O_MEMVAR: alpar@1: memb.value.var = eval_member_var(mpl, code->arg.var.var, alpar@1: memb.tuple); alpar@1: display_var alpar@1: (mpl, code->arg.var.var, &memb, code->arg.var.suff); alpar@1: break; alpar@1: case O_MEMCON: alpar@1: memb.value.con = eval_member_con(mpl, code->arg.con.con, alpar@1: memb.tuple); alpar@1: display_con alpar@1: (mpl, code->arg.con.con, &memb, code->arg.con.suff); alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: delete_tuple(mpl, memb.tuple); alpar@1: return; alpar@1: } alpar@1: alpar@1: static void display_code(MPL *mpl, CODE *code) alpar@1: { /* display value of expression */ alpar@1: switch (code->type) alpar@1: { case A_NUMERIC: alpar@1: /* numeric value */ alpar@1: { double num; alpar@1: num = eval_numeric(mpl, code); alpar@1: write_text(mpl, "%.*g\n", DBL_DIG, num); alpar@1: } alpar@1: break; alpar@1: case A_SYMBOLIC: alpar@1: /* symbolic value */ alpar@1: { SYMBOL *sym; alpar@1: sym = eval_symbolic(mpl, code); alpar@1: write_text(mpl, "%s\n", format_symbol(mpl, sym)); alpar@1: delete_symbol(mpl, sym); alpar@1: } alpar@1: break; alpar@1: case A_LOGICAL: alpar@1: /* logical value */ alpar@1: { int bit; alpar@1: bit = eval_logical(mpl, code); alpar@1: write_text(mpl, "%s\n", bit ? "true" : "false"); alpar@1: } alpar@1: break; alpar@1: case A_TUPLE: alpar@1: /* n-tuple */ alpar@1: { TUPLE *tuple; alpar@1: tuple = eval_tuple(mpl, code); alpar@1: write_text(mpl, "%s\n", format_tuple(mpl, '(', tuple)); alpar@1: delete_tuple(mpl, tuple); alpar@1: } alpar@1: break; alpar@1: case A_ELEMSET: alpar@1: /* elemental set */ alpar@1: { ELEMSET *set; alpar@1: MEMBER *memb; alpar@1: set = eval_elemset(mpl, code); alpar@1: if (set->head == 0) alpar@1: write_text(mpl, "set is empty\n"); alpar@1: for (memb = set->head; memb != NULL; memb = memb->next) alpar@1: write_text(mpl, " %s\n", format_tuple(mpl, '(', alpar@1: memb->tuple)); alpar@1: delete_elemset(mpl, set); alpar@1: } alpar@1: break; alpar@1: case A_FORMULA: alpar@1: /* linear form */ alpar@1: { FORMULA *form, *term; alpar@1: form = eval_formula(mpl, code); alpar@1: if (form == NULL) alpar@1: write_text(mpl, "linear form is empty\n"); alpar@1: for (term = form; term != NULL; term = term->next) alpar@1: { if (term->var == NULL) alpar@1: write_text(mpl, " %.*g\n", term->coef); alpar@1: else alpar@1: write_text(mpl, " %.*g %s%s\n", DBL_DIG, alpar@1: term->coef, term->var->var->name, alpar@1: format_tuple(mpl, '[', term->var->memb->tuple)); alpar@1: } alpar@1: delete_formula(mpl, form); alpar@1: } alpar@1: break; alpar@1: default: alpar@1: xassert(code != code); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: static int display_func(MPL *mpl, void *info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: DISPLAY *dpy = (DISPLAY *)info; alpar@1: DISPLAY1 *entry; alpar@1: for (entry = dpy->list; entry != NULL; entry = entry->next) alpar@1: { if (entry->type == A_INDEX) alpar@1: { /* dummy index */ alpar@1: DOMAIN_SLOT *slot = entry->u.slot; alpar@1: write_text(mpl, "%s = %s\n", slot->name, alpar@1: format_symbol(mpl, slot->value)); alpar@1: } alpar@1: else if (entry->type == A_SET) alpar@1: { /* model set */ alpar@1: SET *set = entry->u.set; alpar@1: MEMBER *memb; alpar@1: if (set->assign != NULL) alpar@1: { /* the set has assignment expression; evaluate all its alpar@1: members over entire domain */ alpar@1: eval_whole_set(mpl, set); alpar@1: } alpar@1: else alpar@1: { /* the set has no assignment expression; refer to its alpar@1: any existing member ignoring resultant value to check alpar@1: the data provided the data section */ alpar@1: #if 1 /* 12/XII-2008 */ alpar@1: if (set->gadget != NULL && set->data == 0) alpar@1: { /* initialize the set with data from a plain set */ alpar@1: saturate_set(mpl, set); alpar@1: } alpar@1: #endif alpar@1: if (set->array->head != NULL) alpar@1: eval_member_set(mpl, set, set->array->head->tuple); alpar@1: } alpar@1: /* display all members of the set array */ alpar@1: if (set->array->head == NULL) alpar@1: write_text(mpl, "%s has empty content\n", set->name); alpar@1: for (memb = set->array->head; memb != NULL; memb = alpar@1: memb->next) display_set(mpl, set, memb); alpar@1: } alpar@1: else if (entry->type == A_PARAMETER) alpar@1: { /* model parameter */ alpar@1: PARAMETER *par = entry->u.par; alpar@1: MEMBER *memb; alpar@1: if (par->assign != NULL) alpar@1: { /* the parameter has an assignment expression; evaluate alpar@1: all its member over entire domain */ alpar@1: eval_whole_par(mpl, par); alpar@1: } alpar@1: else alpar@1: { /* the parameter has no assignment expression; refer to alpar@1: its any existing member ignoring resultant value to alpar@1: check the data provided in the data section */ alpar@1: if (par->array->head != NULL) alpar@1: { if (par->type != A_SYMBOLIC) alpar@1: eval_member_num(mpl, par, par->array->head->tuple); alpar@1: else alpar@1: delete_symbol(mpl, eval_member_sym(mpl, par, alpar@1: par->array->head->tuple)); alpar@1: } alpar@1: } alpar@1: /* display all members of the parameter array */ alpar@1: if (par->array->head == NULL) alpar@1: write_text(mpl, "%s has empty content\n", par->name); alpar@1: for (memb = par->array->head; memb != NULL; memb = alpar@1: memb->next) display_par(mpl, par, memb); alpar@1: } alpar@1: else if (entry->type == A_VARIABLE) alpar@1: { /* model variable */ alpar@1: VARIABLE *var = entry->u.var; alpar@1: MEMBER *memb; alpar@1: xassert(mpl->flag_p); alpar@1: /* display all members of the variable array */ alpar@1: if (var->array->head == NULL) alpar@1: write_text(mpl, "%s has empty content\n", var->name); alpar@1: for (memb = var->array->head; memb != NULL; memb = alpar@1: memb->next) display_var(mpl, var, memb, DOT_NONE); alpar@1: } alpar@1: else if (entry->type == A_CONSTRAINT) alpar@1: { /* model constraint */ alpar@1: CONSTRAINT *con = entry->u.con; alpar@1: MEMBER *memb; alpar@1: xassert(mpl->flag_p); alpar@1: /* display all members of the constraint array */ alpar@1: if (con->array->head == NULL) alpar@1: write_text(mpl, "%s has empty content\n", con->name); alpar@1: for (memb = con->array->head; memb != NULL; memb = alpar@1: memb->next) display_con(mpl, con, memb, DOT_NONE); alpar@1: } alpar@1: else if (entry->type == A_EXPRESSION) alpar@1: { /* expression */ alpar@1: CODE *code = entry->u.code; alpar@1: if (code->op == O_MEMNUM || code->op == O_MEMSYM || alpar@1: code->op == O_MEMSET || code->op == O_MEMVAR || alpar@1: code->op == O_MEMCON) alpar@1: display_memb(mpl, code); alpar@1: else alpar@1: display_code(mpl, code); alpar@1: } alpar@1: else alpar@1: xassert(entry != entry); alpar@1: } alpar@1: return 0; alpar@1: } alpar@1: alpar@1: void execute_display(MPL *mpl, DISPLAY *dpy) alpar@1: { loop_within_domain(mpl, dpy->domain, dpy, display_func); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean_display - clean display statement. alpar@1: -- alpar@1: -- This routine cleans specified display statement that assumes deleting alpar@1: -- all stuff dynamically allocated on generating/postsolving phase. */ alpar@1: alpar@1: void clean_display(MPL *mpl, DISPLAY *dpy) alpar@1: { DISPLAY1 *d; alpar@1: #if 0 /* 15/V-2010 */ alpar@1: ARG_LIST *e; alpar@1: #endif alpar@1: /* clean subscript domain */ alpar@1: clean_domain(mpl, dpy->domain); alpar@1: /* clean display list */ alpar@1: for (d = dpy->list; d != NULL; d = d->next) alpar@1: { /* clean pseudo-code for computing expression */ alpar@1: if (d->type == A_EXPRESSION) alpar@1: clean_code(mpl, d->u.code); alpar@1: #if 0 /* 15/V-2010 */ alpar@1: /* clean pseudo-code for computing subscripts */ alpar@1: for (e = d->list; e != NULL; e = e->next) alpar@1: clean_code(mpl, e->x); alpar@1: #endif alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- execute_printf - execute printf statement. alpar@1: -- alpar@1: -- This routine executes specified printf statement. */ alpar@1: alpar@1: #if 1 /* 14/VII-2006 */ alpar@1: static void print_char(MPL *mpl, int c) alpar@1: { if (mpl->prt_fp == NULL) alpar@1: write_char(mpl, c); alpar@1: else alpar@1: xfputc(c, mpl->prt_fp); alpar@1: return; alpar@1: } alpar@1: alpar@1: static void print_text(MPL *mpl, char *fmt, ...) alpar@1: { va_list arg; alpar@1: char buf[OUTBUF_SIZE], *c; alpar@1: va_start(arg, fmt); alpar@1: vsprintf(buf, fmt, arg); alpar@1: xassert(strlen(buf) < sizeof(buf)); alpar@1: va_end(arg); alpar@1: for (c = buf; *c != '\0'; c++) print_char(mpl, *c); alpar@1: return; alpar@1: } alpar@1: #endif alpar@1: alpar@1: static int printf_func(MPL *mpl, void *info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: PRINTF *prt = (PRINTF *)info; alpar@1: PRINTF1 *entry; alpar@1: SYMBOL *sym; alpar@1: char fmt[MAX_LENGTH+1], *c, *from, save; alpar@1: /* evaluate format control string */ alpar@1: sym = eval_symbolic(mpl, prt->fmt); alpar@1: if (sym->str == NULL) alpar@1: sprintf(fmt, "%.*g", DBL_DIG, sym->num); alpar@1: else alpar@1: fetch_string(mpl, sym->str, fmt); alpar@1: delete_symbol(mpl, sym); alpar@1: /* scan format control string and perform formatting output */ alpar@1: entry = prt->list; alpar@1: for (c = fmt; *c != '\0'; c++) alpar@1: { if (*c == '%') alpar@1: { /* scan format specifier */ alpar@1: from = c++; alpar@1: if (*c == '%') alpar@1: { print_char(mpl, '%'); alpar@1: continue; alpar@1: } alpar@1: if (entry == NULL) break; alpar@1: /* scan optional flags */ alpar@1: while (*c == '-' || *c == '+' || *c == ' ' || *c == '#' || alpar@1: *c == '0') c++; alpar@1: /* scan optional minimum field width */ alpar@1: while (isdigit((unsigned char)*c)) c++; alpar@1: /* scan optional precision */ alpar@1: if (*c == '.') alpar@1: { c++; alpar@1: while (isdigit((unsigned char)*c)) c++; alpar@1: } alpar@1: /* scan conversion specifier and perform formatting */ alpar@1: save = *(c+1), *(c+1) = '\0'; alpar@1: if (*c == 'd' || *c == 'i' || *c == 'e' || *c == 'E' || alpar@1: *c == 'f' || *c == 'F' || *c == 'g' || *c == 'G') alpar@1: { /* the specifier requires numeric value */ alpar@1: double value; alpar@1: xassert(entry != NULL); alpar@1: switch (entry->code->type) alpar@1: { case A_NUMERIC: alpar@1: value = eval_numeric(mpl, entry->code); alpar@1: break; alpar@1: case A_SYMBOLIC: alpar@1: sym = eval_symbolic(mpl, entry->code); alpar@1: if (sym->str != NULL) alpar@1: error(mpl, "cannot convert %s to floating-point" alpar@1: " number", format_symbol(mpl, sym)); alpar@1: value = sym->num; alpar@1: delete_symbol(mpl, sym); alpar@1: break; alpar@1: case A_LOGICAL: alpar@1: if (eval_logical(mpl, entry->code)) alpar@1: value = 1.0; alpar@1: else alpar@1: value = 0.0; alpar@1: break; alpar@1: default: alpar@1: xassert(entry != entry); alpar@1: } alpar@1: if (*c == 'd' || *c == 'i') alpar@1: { double int_max = (double)INT_MAX; alpar@1: if (!(-int_max <= value && value <= +int_max)) alpar@1: error(mpl, "cannot convert %.*g to integer", alpar@1: DBL_DIG, value); alpar@1: print_text(mpl, from, (int)floor(value + 0.5)); alpar@1: } alpar@1: else alpar@1: print_text(mpl, from, value); alpar@1: } alpar@1: else if (*c == 's') alpar@1: { /* the specifier requires symbolic value */ alpar@1: char value[MAX_LENGTH+1]; alpar@1: switch (entry->code->type) alpar@1: { case A_NUMERIC: alpar@1: sprintf(value, "%.*g", DBL_DIG, eval_numeric(mpl, alpar@1: entry->code)); alpar@1: break; alpar@1: case A_LOGICAL: alpar@1: if (eval_logical(mpl, entry->code)) alpar@1: strcpy(value, "T"); alpar@1: else alpar@1: strcpy(value, "F"); alpar@1: break; alpar@1: case A_SYMBOLIC: alpar@1: sym = eval_symbolic(mpl, entry->code); alpar@1: if (sym->str == NULL) alpar@1: sprintf(value, "%.*g", DBL_DIG, sym->num); alpar@1: else alpar@1: fetch_string(mpl, sym->str, value); alpar@1: delete_symbol(mpl, sym); alpar@1: break; alpar@1: default: alpar@1: xassert(entry != entry); alpar@1: } alpar@1: print_text(mpl, from, value); alpar@1: } alpar@1: else alpar@1: error(mpl, "format specifier missing or invalid"); alpar@1: *(c+1) = save; alpar@1: entry = entry->next; alpar@1: } alpar@1: else if (*c == '\\') alpar@1: { /* write some control character */ alpar@1: c++; alpar@1: if (*c == 't') alpar@1: print_char(mpl, '\t'); alpar@1: else if (*c == 'n') alpar@1: print_char(mpl, '\n'); alpar@1: #if 1 /* 28/X-2010 */ alpar@1: else if (*c == '\0') alpar@1: { /* format string ends with backslash */ alpar@1: error(mpl, "invalid use of escape character \\ in format" alpar@1: " control string"); alpar@1: } alpar@1: #endif alpar@1: else alpar@1: print_char(mpl, *c); alpar@1: } alpar@1: else alpar@1: { /* write character without formatting */ alpar@1: print_char(mpl, *c); alpar@1: } alpar@1: } alpar@1: return 0; alpar@1: } alpar@1: alpar@1: #if 0 /* 14/VII-2006 */ alpar@1: void execute_printf(MPL *mpl, PRINTF *prt) alpar@1: { loop_within_domain(mpl, prt->domain, prt, printf_func); alpar@1: return; alpar@1: } alpar@1: #else alpar@1: void execute_printf(MPL *mpl, PRINTF *prt) alpar@1: { if (prt->fname == NULL) alpar@1: { /* switch to the standard output */ alpar@1: if (mpl->prt_fp != NULL) alpar@1: { xfclose(mpl->prt_fp), mpl->prt_fp = NULL; alpar@1: xfree(mpl->prt_file), mpl->prt_file = NULL; alpar@1: } alpar@1: } alpar@1: else alpar@1: { /* evaluate file name string */ alpar@1: SYMBOL *sym; alpar@1: char fname[MAX_LENGTH+1]; alpar@1: sym = eval_symbolic(mpl, prt->fname); alpar@1: if (sym->str == NULL) alpar@1: sprintf(fname, "%.*g", DBL_DIG, sym->num); alpar@1: else alpar@1: fetch_string(mpl, sym->str, fname); alpar@1: delete_symbol(mpl, sym); alpar@1: /* close the current print file, if necessary */ alpar@1: if (mpl->prt_fp != NULL && alpar@1: (!prt->app || strcmp(mpl->prt_file, fname) != 0)) alpar@1: { xfclose(mpl->prt_fp), mpl->prt_fp = NULL; alpar@1: xfree(mpl->prt_file), mpl->prt_file = NULL; alpar@1: } alpar@1: /* open the specified print file, if necessary */ alpar@1: if (mpl->prt_fp == NULL) alpar@1: { mpl->prt_fp = xfopen(fname, prt->app ? "a" : "w"); alpar@1: if (mpl->prt_fp == NULL) alpar@1: error(mpl, "unable to open `%s' for writing - %s", alpar@1: fname, xerrmsg()); alpar@1: mpl->prt_file = xmalloc(strlen(fname)+1); alpar@1: strcpy(mpl->prt_file, fname); alpar@1: } alpar@1: } alpar@1: loop_within_domain(mpl, prt->domain, prt, printf_func); alpar@1: if (mpl->prt_fp != NULL) alpar@1: { xfflush(mpl->prt_fp); alpar@1: if (xferror(mpl->prt_fp)) alpar@1: error(mpl, "writing error to `%s' - %s", mpl->prt_file, alpar@1: xerrmsg()); alpar@1: } alpar@1: return; alpar@1: } alpar@1: #endif alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean_printf - clean printf statement. alpar@1: -- alpar@1: -- This routine cleans specified printf statement that assumes deleting alpar@1: -- all stuff dynamically allocated on generating/postsolving phase. */ alpar@1: alpar@1: void clean_printf(MPL *mpl, PRINTF *prt) alpar@1: { PRINTF1 *p; alpar@1: /* clean subscript domain */ alpar@1: clean_domain(mpl, prt->domain); alpar@1: /* clean pseudo-code for computing format string */ alpar@1: clean_code(mpl, prt->fmt); alpar@1: /* clean printf list */ alpar@1: for (p = prt->list; p != NULL; p = p->next) alpar@1: { /* clean pseudo-code for computing value to be printed */ alpar@1: clean_code(mpl, p->code); alpar@1: } alpar@1: #if 1 /* 14/VII-2006 */ alpar@1: /* clean pseudo-code for computing file name string */ alpar@1: clean_code(mpl, prt->fname); alpar@1: #endif alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- execute_for - execute for statement. alpar@1: -- alpar@1: -- This routine executes specified for statement. */ alpar@1: alpar@1: static int for_func(MPL *mpl, void *info) alpar@1: { /* this is auxiliary routine to work within domain scope */ alpar@1: FOR *fur = (FOR *)info; alpar@1: STATEMENT *stmt, *save; alpar@1: save = mpl->stmt; alpar@1: for (stmt = fur->list; stmt != NULL; stmt = stmt->next) alpar@1: execute_statement(mpl, stmt); alpar@1: mpl->stmt = save; alpar@1: return 0; alpar@1: } alpar@1: alpar@1: void execute_for(MPL *mpl, FOR *fur) alpar@1: { loop_within_domain(mpl, fur->domain, fur, for_func); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean_for - clean for statement. alpar@1: -- alpar@1: -- This routine cleans specified for statement that assumes deleting all alpar@1: -- stuff dynamically allocated on generating/postsolving phase. */ alpar@1: alpar@1: void clean_for(MPL *mpl, FOR *fur) alpar@1: { STATEMENT *stmt; alpar@1: /* clean subscript domain */ alpar@1: clean_domain(mpl, fur->domain); alpar@1: /* clean all sub-statements */ alpar@1: for (stmt = fur->list; stmt != NULL; stmt = stmt->next) alpar@1: clean_statement(mpl, stmt); alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- execute_statement - execute specified model statement. alpar@1: -- alpar@1: -- This routine executes specified model statement. */ alpar@1: alpar@1: void execute_statement(MPL *mpl, STATEMENT *stmt) alpar@1: { mpl->stmt = stmt; alpar@1: switch (stmt->type) alpar@1: { case A_SET: alpar@1: case A_PARAMETER: alpar@1: case A_VARIABLE: alpar@1: break; alpar@1: case A_CONSTRAINT: alpar@1: xprintf("Generating %s...\n", stmt->u.con->name); alpar@1: eval_whole_con(mpl, stmt->u.con); alpar@1: break; alpar@1: case A_TABLE: alpar@1: switch (stmt->u.tab->type) alpar@1: { case A_INPUT: alpar@1: xprintf("Reading %s...\n", stmt->u.tab->name); alpar@1: break; alpar@1: case A_OUTPUT: alpar@1: xprintf("Writing %s...\n", stmt->u.tab->name); alpar@1: break; alpar@1: default: alpar@1: xassert(stmt != stmt); alpar@1: } alpar@1: execute_table(mpl, stmt->u.tab); alpar@1: break; alpar@1: case A_SOLVE: alpar@1: break; alpar@1: case A_CHECK: alpar@1: xprintf("Checking (line %d)...\n", stmt->line); alpar@1: execute_check(mpl, stmt->u.chk); alpar@1: break; alpar@1: case A_DISPLAY: alpar@1: write_text(mpl, "Display statement at line %d\n", alpar@1: stmt->line); alpar@1: execute_display(mpl, stmt->u.dpy); alpar@1: break; alpar@1: case A_PRINTF: alpar@1: execute_printf(mpl, stmt->u.prt); alpar@1: break; alpar@1: case A_FOR: alpar@1: execute_for(mpl, stmt->u.fur); alpar@1: break; alpar@1: default: alpar@1: xassert(stmt != stmt); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /*---------------------------------------------------------------------- alpar@1: -- clean_statement - clean specified model statement. alpar@1: -- alpar@1: -- This routine cleans specified model statement that assumes deleting alpar@1: -- all stuff dynamically allocated on generating/postsolving phase. */ alpar@1: alpar@1: void clean_statement(MPL *mpl, STATEMENT *stmt) alpar@1: { switch(stmt->type) alpar@1: { case A_SET: alpar@1: clean_set(mpl, stmt->u.set); break; alpar@1: case A_PARAMETER: alpar@1: clean_parameter(mpl, stmt->u.par); break; alpar@1: case A_VARIABLE: alpar@1: clean_variable(mpl, stmt->u.var); break; alpar@1: case A_CONSTRAINT: alpar@1: clean_constraint(mpl, stmt->u.con); break; alpar@1: #if 1 /* 11/II-2008 */ alpar@1: case A_TABLE: alpar@1: clean_table(mpl, stmt->u.tab); break; alpar@1: #endif alpar@1: case A_SOLVE: alpar@1: break; alpar@1: case A_CHECK: alpar@1: clean_check(mpl, stmt->u.chk); break; alpar@1: case A_DISPLAY: alpar@1: clean_display(mpl, stmt->u.dpy); break; alpar@1: case A_PRINTF: alpar@1: clean_printf(mpl, stmt->u.prt); break; alpar@1: case A_FOR: alpar@1: clean_for(mpl, stmt->u.fur); break; alpar@1: default: alpar@1: xassert(stmt != stmt); alpar@1: } alpar@1: return; alpar@1: } alpar@1: alpar@1: /* eof */