From aee256e211673c09a9c3e0140da50c749b2b175e Mon Sep 17 00:00:00 2001 From: Jim Monte Date: Sun, 12 Jan 2020 23:02:36 -0500 Subject: [PATCH] Fix of buffer overrun in interpolation at endpoint of interval. Made cfunc.mod for tables more modular. Prevented buffer overrun when building file name. Added error checking for allocation failures in many locations. Made binary search for interpolation more efficient. --- src/xspice/icm/dlmain.c | 105 ++- src/xspice/icm/table/mada/alloc.c | 20 +- src/xspice/icm/table/mada/eno.c | 57 +- src/xspice/icm/table/mada/eno2.c | 71 +- src/xspice/icm/table/mada/eno3.c | 85 +- src/xspice/icm/table/support/gettokens.c | 163 +++- src/xspice/icm/table/support/interp.c | 87 +- src/xspice/icm/table/table2D/cfunc.mod | 945 +++++++++---------- src/xspice/icm/table/table3D/cfunc.mod | 1055 +++++++++++----------- visualc/xspice/analog.vcxproj | 6 +- visualc/xspice/digital.vcxproj | 6 +- visualc/xspice/spice2poly.vcxproj | 6 +- visualc/xspice/table.vcxproj | 15 +- visualc/xspice/xtradev.vcxproj | 6 +- visualc/xspice/xtraevt.vcxproj | 6 +- 15 files changed, 1484 insertions(+), 1149 deletions(-) diff --git a/src/xspice/icm/dlmain.c b/src/xspice/icm/dlmain.c index d993ff374..93afe745e 100644 --- a/src/xspice/icm/dlmain.c +++ b/src/xspice/icm/dlmain.c @@ -7,17 +7,18 @@ // Author: Arpad Buermen ////////////////////////////////////////////////////////////////////////////// -#include "ngspice/inpdefs.h" -#include "ngspice/devdefs.h" -#include "ngspice/evtudn.h" -#include "ngspice/dllitf.h" -#include "cmextrn.h" -#include "udnextrn.h" - +#include #include #include -#include +#include "ngspice/devdefs.h" +#include "ngspice/dstring.h" +#include "ngspice/dllitf.h" +#include "ngspice/evtudn.h" +#include "ngspice/inpdefs.h" +#include "cmextrn.h" +#include "udnextrn.h" + ////////////////////////////////////////////////////////////////////////////// @@ -416,44 +417,76 @@ Infile_Path/ NGSPICE_INPUT_DIR/, where the path is given by the environmental variable , where the path is the current directory */ - -#define MAX_PATH_LEN 1024 - +#define DFLT_BUF_SIZE 256 FILE *fopen_with_path(const char *path, const char *mode) { - char buf[MAX_PATH_LEN+1]; FILE *fp; - if((path[0] != '/') && (path[1] != ':')) { + if((path[0] != '/') && (path[1] != ':')) { /* path absolue (probably) */ // const char *x = getenv("ngspice_vpath"); const char *x = cm_get_path(); - if(x) { - strncpy(buf, x, MAX_PATH_LEN); - strcat(buf, "/"); - strcat(buf, path); - fp = fopen(buf, mode); - if (fp) - return fp; - else { - char *y = getenv( "NGSPICE_INPUT_DIR" ); - if (y && *y) { - char *a; - strncpy(buf, y, MAX_PATH_LEN); - a = strrchr(buf, '/'); - if(a && a[1] == '\0') - a[0] = '\0'; - strcat(buf, "/"); - strcat(buf, path); - fp = fopen(buf, mode); - if (fp) - return fp; - } + if (x) { + DS_CREATE(ds, DFLT_BUF_SIZE); + + /* Build file /path> */ + if (ds_cat_printf(&ds, "%s/%s", x, path) != 0) { + cm_message_printf( + "Unable to build cm_get_path() path for opening file."); + ds_free(&ds); + return (FILE *) NULL; } + + /* Try opening file. If fail, try using NGSPICE_INPUT_DIR + * env variable location */ + if ((fp = fopen(ds_get_buf(&ds), mode)) == (FILE *) NULL) { + char *y = getenv("NGSPICE_INPUT_DIR"); + if (y && *y) { /* have env var and not "" */ + int rc_ds = 0; + /* Build /path and try opening. If the env var + * ends with a slash, do not add a second slash */ + ds_clear(&ds); + rc_ds |= ds_cat_str(&ds, y); + + /* Add slash if not present. Note that check for + * length > 0 is done on the remote chance that the + * ds_cat_str() failed. */ + const size_t len = ds_get_length(&ds); + if (len > 0 && ds_get_buf(&ds)[len - 1] != '/') { + rc_ds |= ds_cat_char(&ds, '/'); /* add dir sep */ + } + rc_ds |= ds_cat_str(&ds, path); /* add input path */ + + /* Ensure path built OK */ + if (rc_ds != 0) { + cm_message_printf( + "Unable to build NGSPICE_INPUT_DIR " + "path for opening file."); + ds_free(&ds); + return (FILE *) NULL; + } + + /* Try opening file name that was built */ + if ((fp = fopen(ds_get_buf(&ds), + mode)) != (FILE *) NULL) { + ds_free(&ds); + return fp; + } + } + } /* end of open using prefix from cm_get_path() failed */ + else { /* Opened OK */ + ds_free(&ds); + return fp; + } + ds_free(&ds); /* free dstring resources, if any */ } - } + } /* end of case that path is not absolute */ + + /* If not opened yet, try opening exactly as given */ fp = fopen(path, mode); + return fp; -} +} /* end of function fopen_with_path */ + int diff --git a/src/xspice/icm/table/mada/alloc.c b/src/xspice/icm/table/mada/alloc.c index cc0753581..b45d611aa 100644 --- a/src/xspice/icm/table/mada/alloc.c +++ b/src/xspice/icm/table/mada/alloc.c @@ -33,15 +33,19 @@ sf_alloc(int n /* number of elements */, { void *ptr; - size *= (size_t) n; + if (n < 0) { + cm_message_printf("%s: illegal allocation(%d X %zd bytes)", + __FILE__, n, size); + return NULL; + } - if (0 >= size) - cm_message_printf("%s: illegal allocation(%d bytes)", __FILE__, (int) size); - - ptr = malloc(size); - - if (NULL == ptr) - cm_message_printf("%s: cannot allocate %d bytes : ", __FILE__, (int) size); + /* Use calloc so that any internal allocations will be set to NULL to + * facilitate error recovery */ + if ((ptr = calloc(n, size)) == NULL) { + cm_message_printf("%s: cannot allocate %zd bytes : ", + __FILE__, n * size); + return NULL; + } return ptr; } diff --git a/src/xspice/icm/table/mada/eno.c b/src/xspice/icm/table/mada/eno.c index e6fc11126..36acbe338 100644 --- a/src/xspice/icm/table/mada/eno.c +++ b/src/xspice/icm/table/mada/eno.c @@ -21,8 +21,9 @@ #include #include -#include "eno.h" +#include "../xspice/icm/dlmain.h" #include "alloc.h" +#include "eno.h" #define SF_MAX(a,b) ((a) < (b) ? (b) : (a)) #define SF_MIN(a,b) ((a) < (b) ? (a) : (b)) @@ -39,16 +40,49 @@ sf_eno_init (int order, /* interpolation order */ int n /* data size */) /*< Initialize interpolation object. >*/ { - sf_eno ent; - int i; + int xrc = 0; + sf_eno ent = (sf_eno) NULL; + + if ((ent = (sf_eno) sf_alloc( + 1, sizeof(*ent))) == (sf_eno) NULL) { + cm_message_printf("Unable to allocate sf_eno structure " + "in sf_eno_init"); + xrc = -1; + goto EXITPOINT; + } - ent = (sf_eno) sf_alloc(1, sizeof(*ent)); ent->order = order; ent->n = n; - ent->diff = (double**) sf_alloc(order, sizeof(double*)); - for (i = 0; i < order; i++) - ent->diff[i] = sf_doublealloc(n - i); + + if ((ent->diff = (double **) sf_alloc( + order, sizeof(double *))) == (double **) NULL) { + cm_message_printf("Unable to allocate diff field " + "in sf_eno_init"); + xrc = -1; + goto EXITPOINT; + } + { + int i; + for (i = 0; i < order; i++) { + if ((ent->diff[i] = sf_doublealloc( + n - i)) == (double *) NULL) { + cm_message_printf("Unable to allocate in sf_eno_init " + "at index %d.", + i); + xrc = -1; + goto EXITPOINT; + } + } + } + +EXITPOINT: + if (xrc != 0) { + if (ent != (sf_eno) NULL) { + sf_eno_close(ent); + ent = (sf_eno) NULL; + } + } return ent; } @@ -56,10 +90,15 @@ void sf_eno_close (sf_eno ent) /*< Free internal storage >*/ { - int i; + if (ent == (sf_eno) NULL) { + return; + } - for (i = 0; i < ent->order; i++) + int i; + const int n = ent->order; + for (i = 0; i < n; i++) { free (ent->diff[i]); + } free (ent->diff); free (ent); } diff --git a/src/xspice/icm/table/mada/eno2.c b/src/xspice/icm/table/mada/eno2.c index 1adb003e0..e8a0cc867 100644 --- a/src/xspice/icm/table/mada/eno2.c +++ b/src/xspice/icm/table/mada/eno2.c @@ -36,23 +36,74 @@ sf_eno2_init (int order, /* interpolation order */ int n1, int n2 /* data dimensions */) /*< Initialize interpolation object >*/ { - sf_eno2 pnt; - int i2; + int xrc = 0; + sf_eno2 pnt = (sf_eno2) NULL; + + if ((pnt = (sf_eno2) sf_alloc( + 1, sizeof(*pnt))) == (sf_eno2) NULL) { + cm_message_printf("Unable to allocate sf_eno2 structure " + "in sf_eno2_init"); + xrc = -1; + goto EXITPOINT; + } - pnt = (sf_eno2) sf_alloc(1, sizeof(*pnt)); pnt->order = order; pnt->n1 = n1; pnt->n2 = n2; pnt->ng = 2 * order - 2; - if (pnt->ng > pnt->n2) + if (pnt->ng > pnt->n2) { cm_message_printf("%s: ng=%d is too big", __FILE__, pnt->ng); - pnt->jnt = sf_eno_init (order, pnt->ng); - pnt->f = sf_doublealloc (pnt->ng); - pnt->f1 = sf_doublealloc (pnt->ng); - pnt->ent = (sf_eno*) sf_alloc(n2, sizeof(sf_eno)); - for (i2 = 0; i2 < n2; i2++) - pnt->ent[i2] = sf_eno_init (order, n1); + xrc = -1; + goto EXITPOINT; + } + if ((pnt->jnt = sf_eno_init(order, pnt->ng)) == (sf_eno) NULL) { + cm_message_printf("Unable to initialize field jnt " + "in sf_eno2_init"); + xrc = -1; + goto EXITPOINT; + } + + if ((pnt->f = sf_doublealloc (pnt->ng)) == (double *) NULL) { + cm_message_printf("Unable to allocate field f in sf_eno2_init()"); + xrc = -1; + goto EXITPOINT; + } + + if ((pnt->f1 = sf_doublealloc (pnt->ng)) == (double *) NULL) { + cm_message_printf("Unable to allocate field f1 in sf_eno2_init()"); + xrc = -1; + goto EXITPOINT; + } + + if ((pnt->ent = (sf_eno *) sf_alloc( + n2, sizeof(sf_eno))) == (sf_eno *) NULL) { + cm_message_printf("Unable to allocate field ent in sf_eno2_init()"); + xrc = -1; + goto EXITPOINT; + } + + { + int i2; + for (i2 = 0; i2 < n2; i2++) { + if ((pnt->ent[i2] = sf_eno_init( + order, n1)) == (sf_eno) NULL) { + cm_message_printf("Unable to initialize field ent[%d] " + "in sf_eno3_init()", + i2); + xrc = -1; + goto EXITPOINT; + } + } + } + +EXITPOINT: + if (xrc != 0) { + if (pnt != (sf_eno2) NULL) { + sf_eno2_close(pnt); + pnt = (sf_eno2) NULL; + } + } return pnt; } diff --git a/src/xspice/icm/table/mada/eno3.c b/src/xspice/icm/table/mada/eno3.c index 2942f142c..a5cb76f64 100644 --- a/src/xspice/icm/table/mada/eno3.c +++ b/src/xspice/icm/table/mada/eno3.c @@ -37,27 +37,88 @@ sf_eno3_init(int order, /* interpolation order */ int n1, int n2, int n3 /* data dimensions */) /*< Initialize interpolation object >*/ { - sf_eno3 pnt; - int i2, i3; + int xrc = 0; + sf_eno3 pnt = (sf_eno3) NULL; + + /* Allocate base structrue */ + if ((pnt = (sf_eno3) sf_alloc (1, sizeof *pnt)) == (sf_eno3) NULL) { + cm_message_printf("Unable to allocate sf_eno3 structure " + "in sf_eno3_init"); + xrc = -1; + goto EXITPOINT; + } - pnt = (sf_eno3) sf_alloc (1, sizeof(*pnt)); pnt->order = order; pnt->n1 = n1; pnt->n2 = n2; pnt->n3 = n3; pnt->ng = 2 * order - 2; - if (pnt->ng > n2 || pnt->ng > n3) + + if (pnt->ng > n2 || pnt->ng > n3) { cm_message_printf("%s: ng=%d is too big", __FILE__, pnt->ng); - pnt->jnt = sf_eno2_init (order, pnt->ng, pnt->ng); - pnt->f = sf_doublealloc2 (pnt->ng, pnt->ng); - pnt->f1 = sf_doublealloc2 (pnt->ng, pnt->ng); - pnt->ent = (sf_eno**) sf_alloc (n3, sizeof(sf_eno*)); - for (i3 = 0; i3 < n3; i3++) { - pnt->ent[i3] = (sf_eno*) sf_alloc (n2, sizeof(sf_eno)); - for (i2 = 0; i2 < n2; i2++) - pnt->ent[i3][i2] = sf_eno_init (order, n1); + xrc = -1; + goto EXITPOINT; } + if ((pnt->jnt = sf_eno2_init( + order, pnt->ng, pnt->ng)) == (sf_eno2) NULL) { + cm_message_printf("Unable to initialize field jnt " + "in sf_eno3_init"); + xrc = -1; + goto EXITPOINT; + } + + if ((pnt->f = sf_doublealloc2(pnt->ng, pnt->ng)) == (double **) NULL) { + cm_message_printf("Unable to allocate field f in sf_eno3_init()"); + xrc = -1; + goto EXITPOINT; + } + + if ((pnt->f1 = sf_doublealloc2(pnt->ng, pnt->ng)) == (double **) NULL) { + cm_message_printf("Unable to allocate field f1 in sf_eno3_init()"); + xrc = -1; + goto EXITPOINT; + } + + if ((pnt->ent = (sf_eno **) sf_alloc( + n3, sizeof(sf_eno*))) == (sf_eno **) NULL) { + cm_message_printf("Unable to allocate field ent in sf_eno3_init()"); + xrc = -1; + goto EXITPOINT; + } + + { + int i3; + for (i3 = 0; i3 < n3; i3++) { + if ((pnt->ent[i3] = (sf_eno*) sf_alloc( + n2, sizeof(sf_eno))) == (sf_eno *) NULL) { + cm_message_printf("Unable to allocate field ent[%d] " + "in sf_eno3_init()", i3); + xrc = -1; + goto EXITPOINT; + } + int i2; + for (i2 = 0; i2 < n2; i2++) { + if ((pnt->ent[i3][i2] = sf_eno_init( + order, n1)) == (sf_eno) NULL) { + cm_message_printf("Unable to initialize field " + "ent[%d][%d] in sf_eno3_init()", + i2, i3); + xrc = -1; + goto EXITPOINT; + } + } + } + } + +EXITPOINT: + if (xrc != 0) { + if (pnt != (sf_eno3) NULL) { + sf_eno3_close(pnt); + free(pnt); + pnt = (sf_eno3) NULL; + } + } return pnt; } diff --git a/src/xspice/icm/table/support/gettokens.c b/src/xspice/icm/table/support/gettokens.c index 4566ad932..4c936c017 100644 --- a/src/xspice/icm/table/support/gettokens.c +++ b/src/xspice/icm/table/support/gettokens.c @@ -12,33 +12,11 @@ is returned. The original input string is undisturbed. #include #include -/*=== CONSTANTS ========================*/ - -#define OK 0 -#define FAIL 1 - -/* Type definition for each possible token returned. */ -typedef enum token_type_s { CNV_NO_TOK, CNV_STRING_TOK } Cnv_Token_Type_t; - -extern char *CNVget_token(char **s, Cnv_Token_Type_t *type); - -/*=== MACROS ===========================*/ - -#if defined(__MINGW32__) || defined(_MSC_VER) -#define DIR_PATHSEP "\\" -#else -#define DIR_PATHSEP "/" -#endif - -#if defined(_MSC_VER) -#define strdup _strdup -#define snprintf _snprintf -#endif +#include "gettokens.h" -char * -CNVgettok(char **s) +char *CNVgettok(char **s) { char *buf; /* temporary storage to copy token into */ /*char *temp;*/ /* temporary storage to copy token into */ @@ -52,7 +30,7 @@ CNVgettok(char **s) /* skip over any white space */ - while (isspace_c(**s) || (**s == '=') || + while (isspace(**s) || (**s == '=') || (**s == '(') || (**s == ')') || (**s == ',')) (*s)++; @@ -71,7 +49,7 @@ CNVgettok(char **s) /* or a mess o' characters. */ i = 0; while ( (**s != '\0') && - (! ( isspace_c(**s) || (**s == '=') || + (! ( isspace(**s) || (**s == '=') || (**s == '(') || (**s == ')') || (**s == ',') ) ) ) { @@ -85,7 +63,7 @@ CNVgettok(char **s) /* skip over white space up to next token */ - while (isspace_c(**s) || (**s == '=') || + while (isspace(**s) || (**s == '=') || (**s == '(') || (**s == ')') || (**s == ',')) (*s)++; @@ -98,11 +76,10 @@ CNVgettok(char **s) if (buf) free(buf); return ret_str; -} +} /* end of function CNVgettok */ -/*=== Static CNVget_token ROUTINE =============*/ /* Get the next token from the input string together with its type. The input string pointer @@ -110,9 +87,7 @@ is advanced to the following token and the token from the input string is copied to malloced storage and a pointer to that storage is returned. The original input string is undisturbed. */ - -char * -CNVget_token(char **s, Cnv_Token_Type_t *type) +char * CNVget_token(char **s, Cnv_Token_Type_t *type) { char *ret_str; /* storage for returned string */ @@ -132,4 +107,126 @@ CNVget_token(char **s, Cnv_Token_Type_t *type) break; } return ret_str; -} +} /* end of function CNVget_token */ + + + +/* + Function takes as input a string token from a SPICE + deck and returns a floating point equivalent value. +*/ +int cnv_get_spice_value(char *str, /* IN - The value text e.g. 1.2K */ + double *p_value) /* OUT - The numerical value */ +{ + /* the following were "int4" devices - jpm */ + size_t len; + size_t i; + int n_matched; + + /* A SPICE size line. <= 80 characters plus '\n\0' */ + typedef char line_t[82]; + line_t val_str; + + char c = ' '; + char c1; + + double scale_factor; + double value; + + /* Scan the input string looking for an alpha character that is not */ + /* 'e' or 'E'. Such a character is assumed to be an engineering */ + /* suffix as defined in the Spice 2G.6 user's manual. */ + + len = strlen(str); + if (len > sizeof(val_str) - 1) + len = sizeof(val_str) - 1; + + for (i = 0; i < len; i++) { + c = str[i]; + if (isalpha(c) && (c != 'E') && (c != 'e')) + break; + else if (isspace(c)) + break; + else + val_str[i] = c; + } + val_str[i] = '\0'; + + /* Determine the scale factor */ + + if ((i >= len) || (! isalpha(c))) + scale_factor = 1.0; + else { + c = (char) tolower(c); + + switch (c) { + + case 't': + scale_factor = 1.0e12; + break; + + case 'g': + scale_factor = 1.0e9; + break; + + case 'k': + scale_factor = 1.0e3; + break; + + case 'u': + scale_factor = 1.0e-6; + break; + + case 'n': + scale_factor = 1.0e-9; + break; + + case 'p': + scale_factor = 1.0e-12; + break; + + case 'f': + scale_factor = 1.0e-15; + break; + + case 'm': + i++; + if (i >= len) { + scale_factor = 1.0e-3; + break; + } + c1 = str[i]; + if (!isalpha(c1)) { + scale_factor = 1.0e-3; + break; + } + c1 = (char) toupper(c1); + if (c1 == 'E') + scale_factor = 1.0e6; + else if (c1 == 'I') + scale_factor = 25.4e-6; + else + scale_factor = 1.0e-3; + break; + + default: + scale_factor = 1.0; + } + } + + /* Convert the numeric portion to a float and multiply by the */ + /* scale factor. */ + + n_matched = sscanf(val_str, "%le", &value); + + if (n_matched < 1) { + *p_value = 0.0; + return -1; + } + + *p_value = value * scale_factor; + return 0; +} /* end of function cnv_get_spice_value */ + + + diff --git a/src/xspice/icm/table/support/interp.c b/src/xspice/icm/table/support/interp.c index 3ace153b5..c3c5d2327 100644 --- a/src/xspice/icm/table/support/interp.c +++ b/src/xspice/icm/table/support/interp.c @@ -13,39 +13,33 @@ typedef struct Point3Struct { /* 3d point */ } Point3; typedef Point3 Vector3; -//FIXME -double BilinearInterpolation(double q11, double q12, double q21, double q22, double x1, double x2, double y1, double y2, double x, double y); - /* Function to find the cross over point (the point before which elements are smaller than or equal to x and after which greater than x) - http://www.geeksforgeeks.org/find-k-closest-elements-given-value/ */ -int -findCrossOver(double arr[], int low, int high, double x) + Returns the highest index of an element in arr whose value is less than + or equal to x or 0 if all elements are greater than x. It is assumed that + arr is sorted in order of increasing values. */ +int findCrossOver(double arr[], int n, double x) { - int mid; - // Base cases - if (arr[high] <= x) // x is greater than all - return high; - if (arr[low] > x) // x is smaller than all - return low; + int low = 0; + int high = n; /* 1 more than highest index */ + while (high - low > 1) { + const int mid = (low + high) / 2; + if (arr[mid] > x) { /* search lower */ + high = mid; + } + else { /* search higher */ + low = mid; + } + } /* end of bisecting loop */ - // Find the middle point - mid = (low + high)/2; /* low + (high - low)/2 */ + return low; +} /* end of function findCrossOver */ - /* If x is same as middle element, then return mid */ - if (arr[mid] <= x && arr[mid+1] > x) - return mid; - /* If x is greater than arr[mid], then either arr[mid + 1] - is ceiling of x or ceiling lies in arr[mid+1...high] */ - if (arr[mid] < x) - return findCrossOver(arr, mid+1, high, x); - - return findCrossOver(arr, low, mid - 1, x); -} +#if 0 /* https://helloacm.com/cc-function-to-compute-the-bilinear-interpolation/ */ double BilinearInterpolation(double q11, double q12, double q21, double q22, double x1, double x2, double y1, double y2, double x, double y) @@ -74,7 +68,6 @@ BilinearInterpolation(double q11, double q12, double q21, double q22, double x1, * */ -#if 0 double trilinear(Point3 *p, double *d, int xsize, int ysize, int zsize, double def) @@ -152,13 +145,31 @@ trilinear(Point3 *p, double *d, int xsize, int ysize, int zsize, double def) #endif +double BilinearInterpolation(double x, double y, + int xind, int yind, double **td) +{ + double V00, V10, V01, V11, Vxyz; + + V00 = td[yind][xind]; + V10 = td[yind][xind+1]; + V01 = td[yind+1][xind]; + V11 = td[yind+1][xind+1]; + + Vxyz = V00 * (1 - x) * (1 - y) + + V10 * x * (1 - y) + + V01 * (1 - x) * y + + V11 * x * y; + return Vxyz; +} /* end of function BilinearInterpolation */ + + + /* trilinear interpolation Paul Bourke July 1997 http://paulbourke.net/miscellaneous/interpolation/ */ - -double -TrilinearInterpolation(double x, double y, double z, int xind, int yind, int zind, double ***td) +double TrilinearInterpolation(double x, double y, double z, + int xind, int yind, int zind, double ***td) { double V000, V100, V010, V001, V101, V011, V110, V111, Vxyz; @@ -172,12 +183,18 @@ TrilinearInterpolation(double x, double y, double z, int xind, int yind, int zin V111 = td[zind+1][yind+1][xind+1]; Vxyz = V000 * (1 - x) * (1 - y) * (1 - z) + - V100 * x * (1 - y) * (1 - z) + - V010 * (1 - x) * y * (1 - z) + - V001 * (1 - x) * (1 - y) * z + - V101 * x * (1 - y) * z + - V011 * (1 - x) * y * z + - V110 * x * y * (1 - z) + - V111 * x * y * z; + V100 * x * (1 - y) * (1 - z) + + V010 * (1 - x) * y * (1 - z) + + V001 * (1 - x) * (1 - y) * z + + V101 * x * (1 - y) * z + + V011 * (1 - x) * y * z + + V110 * x * y * (1 - z) + + V111 * x * y * z; return Vxyz; } + + + + + + diff --git a/src/xspice/icm/table/table2D/cfunc.mod b/src/xspice/icm/table/table2D/cfunc.mod index df8d363a9..77760dae6 100644 --- a/src/xspice/icm/table/table2D/cfunc.mod +++ b/src/xspice/icm/table/table2D/cfunc.mod @@ -62,7 +62,7 @@ NON-STANDARD FEATURES ===============================================================================*/ /*=== INCLUDE FILES ====================*/ - +#include #include #include #include @@ -74,10 +74,21 @@ NON-STANDARD FEATURES #include "mada/eno2.h" -/*=== CONSTANTS ========================*/ +typedef struct { + int ix; /* size of array in x */ + int iy; /* size of array in y */ -#define OK 0 -#define FAIL 1 + sf_eno2 newtable; /* the table, code borrowed from madagascar project */ + + /* Input values corresponding to each index. They define the value + * in the domain at each index value */ + double *xcol; /* array of doubles in x */ + double *ycol; /* array of doubles in y */ + + double **table; /* f(xi, yj) */ +} Table2_Data_t; + +typedef Table2_Data_t Local_Data_t; /*=== MACROS ===========================*/ @@ -100,35 +111,22 @@ struct filesource_state { }; -typedef struct { - int ix; /* size of array in x */ - int iy; /* size of array in y */ - struct filesource_state *state; /* the storage array for the - filesource status. */ - int init_err; - - sf_eno2 newtable; /* the table, code borrowed from madagascar project */ - - double *xcol; /* array of floats in x */ - double *ycol; /* array of floats in y */ - - double **table; - -} Local_Data_t; - -/* Type definition for each possible token returned. */ -typedef enum token_type_s {CNV_NO_TOK, CNV_STRING_TOK} Cnv_Token_Type_t; - -typedef char line_t[82]; /* A SPICE size line. <= 80 characters plus '\n\0' */ /*=== FUNCTION PROTOTYPE DEFINITIONS ===*/ -extern int findCrossOver(double arr[], int low, int high, double x); -extern double BilinearInterpolation(double q11, double q12, double q21, double q22, double x1, double x2, double y1, double y2, double x, double y); - +double BilinearInterpolation(double x, double y, + int xind, int yind, double **td); extern char *CNVgettok(char **s); +int cnv_get_spice_value(char *str, double *p_value); +extern int findCrossOver(double arr[], int n, double x); + +static void free_local_data(Table2_Data_t *loc); +static inline double get_local_diff(int n, double *col, int ind); +static Table2_Data_t *init_local_data(const char *filename, int order); + + /*============================================================================== @@ -168,151 +166,23 @@ NON-STANDARD FEATURES ==============================================================================*/ -/*=== Static CNV_get_spice_value ROUTINE =============*/ - -/* - Function takes as input a string token from a SPICE - deck and returns a floating point equivalent value. -*/ -static int -cnv_get_spice_value(char *str, /* IN - The value text e.g. 1.2K */ - double *p_value) /* OUT - The numerical value */ -{ - /* the following were "int4" devices - jpm */ - size_t len; - size_t i; - int n_matched; - line_t val_str; - - char c = ' '; - char c1; - - double scale_factor; - double value; - - /* Scan the input string looking for an alpha character that is not */ - /* 'e' or 'E'. Such a character is assumed to be an engineering */ - /* suffix as defined in the Spice 2G.6 user's manual. */ - - len = strlen(str); - if (len > sizeof(val_str) - 1) - len = sizeof(val_str) - 1; - - for (i = 0; i < len; i++) { - c = str[i]; - if (isalpha_c(c) && (c != 'E') && (c != 'e')) - break; - else if (isspace_c(c)) - break; - else - val_str[i] = c; - } - val_str[i] = '\0'; - - /* Determine the scale factor */ - - if ((i >= len) || !isalpha_c(c)) - scale_factor = 1.0; - else { - - if (isupper_c(c)) - c = tolower_c(c); - - switch (c) { - - case 't': - scale_factor = 1.0e12; - break; - - case 'g': - scale_factor = 1.0e9; - break; - - case 'k': - scale_factor = 1.0e3; - break; - - case 'u': - scale_factor = 1.0e-6; - break; - - case 'n': - scale_factor = 1.0e-9; - break; - - case 'p': - scale_factor = 1.0e-12; - break; - - case 'f': - scale_factor = 1.0e-15; - break; - - case 'm': - i++; - if (i >= len) { - scale_factor = 1.0e-3; - break; - } - c1 = str[i]; - if (!isalpha_c(c1)) { - scale_factor = 1.0e-3; - break; - } - if (islower_c(c1)) - c1 = toupper_c(c1); - if (c1 == 'E') - scale_factor = 1.0e6; - else if (c1 == 'I') - scale_factor = 25.4e-6; - else - scale_factor = 1.0e-3; - break; - - default: - scale_factor = 1.0; - } - } - - /* Convert the numeric portion to a float and multiply by the */ - /* scale factor. */ - - n_matched = sscanf(val_str, "%le", &value); - - if (n_matched < 1) { - *p_value = 0.0; - return FAIL; - } - - *p_value = value * scale_factor; - return OK; -} - -static void -cm_table2D_callback(ARGS, Mif_Callback_Reason_t reason) +static void cm_table2D_callback(ARGS, + Mif_Callback_Reason_t reason) { switch (reason) { case MIF_CB_DESTROY: { - int i; - Local_Data_t *loc = STATIC_VAR (locdata); - if (!loc) - break; - free(loc->state); - for (i = 0; i < loc->iy; i++) - free(loc->table[i]); - free(loc->table); - free(loc->xcol); - free(loc->ycol); - sf_eno2_close (loc->newtable); - free(loc); - STATIC_VAR (locdata) = loc = NULL; + Table2_Data_t *loc = STATIC_VAR(locdata); + if (loc) { + free_local_data(loc); + STATIC_VAR(locdata) = loc = NULL; + } break; - } - } -} + } /* end of case MIF_CB_DESTROY */ + } /* end of switch over reason being called */ +} /* end of function cm_table2D_callback */ @@ -359,7 +229,7 @@ ix iy * x row independent variables x0 x1 x2 x3 ... xix-1 -y column independent variables +* y column independent variables y0 y1 y2 y3 ... yiy-1 * table x0y0 x1y0 x2y0 ... xix-1y0 @@ -374,339 +244,89 @@ x0yiy-1 x1yiy-1 x2yiy-1 ... xix-1yiy-1 /*=== CM_table2D ROUTINE ===*/ -void -cm_table2D(ARGS) /* structure holding parms, inputs, outputs, etc. */ +void cm_table2D(ARGS) /* structure holding parms, inputs, outputs, etc. */ { int size, xind, yind; double xval, yval, xoff, yoff, xdiff, ydiff; double derivval[2], outval; - double q11, q12, q21, q22, x1, x2, y1, y2; - Local_Data_t *loc; /* Pointer to local static data, not to be included + Table2_Data_t *loc; /* Pointer to local static data, not to be included in the state vector */ - Mif_Complex_t ac_gain; size = PORT_SIZE(out); - if (INIT == 1) { - - int i; - int ix = 0, /* elements in a row */ - iy = 0; /* number of rows */ - double **table_data; - double tmp; - char *cFile, *cThisPtr, *cThisLine, *cThisLinePtr; - int isNewline; /* Boolean indicating we've read a CR or LF */ - size_t lFileLen; /* Length of file */ - size_t lFileRead; /* Length of file read in */ - long lIndex; /* Index into cThisLine array */ - int lLineCount; /* Current line number */ - size_t lStartPos; /* Offset of start of current line */ - size_t lTotalChars; /* Total characters read */ - int interporder; /* order of interpolation for eno */ - + if (INIT == 1) { /* Must do initializations */ + STATIC_VAR(locdata) = init_local_data( + PARAM(file), + PARAM(order)); CALLBACK = cm_table2D_callback; - - /* allocate static storage for *loc */ - STATIC_VAR (locdata) = calloc(1, sizeof(Local_Data_t)); - loc = STATIC_VAR (locdata); - - /* Allocate storage for internal state */ - loc->state = (struct filesource_state*) malloc(sizeof(struct filesource_state)); - loc->ix = loc->iy = 0; - - /* open file */ - loc->state->fp = fopen_with_path(PARAM(file), "r"); - loc->state->pos = 0; - loc->state->atend = 0; - if (!loc->state->fp) { - char *lbuffer, *p; - lbuffer = getenv("NGSPICE_INPUT_DIR"); - if (lbuffer && *lbuffer) { - p = (char*) malloc(strlen(lbuffer) + strlen(DIR_PATHSEP) + strlen(PARAM(file)) + 1); - sprintf(p, "%s%s%s", lbuffer, DIR_PATHSEP, PARAM(file)); - loc->state->fp = fopen(p, "r"); - free(p); - } - } - struct stat st; - if (!loc->state->fp || fstat(fileno(loc->state->fp), &st)) { - cm_message_printf("cannot open file %s", PARAM(file)); - free(loc->state); - free(loc); - STATIC_VAR (locdata) = loc = NULL; - return; - } - /* get file length */ - lFileLen = (size_t) st.st_size; - - /* create string to hold the whole file */ - cFile = calloc(lFileLen + 1, sizeof(char)); - /* create another string long enough for file manipulation */ - cThisLine = calloc(lFileLen + 1, sizeof(char)); - if (cFile == NULL || cThisLine == NULL) { - cm_message_printf("Insufficient memory to read file %s", PARAM(file)); - free(loc->state); - free(loc); - STATIC_VAR (locdata) = loc = NULL; - if(cFile) free(cFile); - if(cThisLine) free(cThisLine); - return; - } - /* read whole file into cFile */ - lFileRead = fread(cFile, sizeof(char), lFileLen, loc->state->fp); - fclose(loc->state->fp); - /* Number of chars read may be less than lFileLen, because /r are skipt by 'fread' */ - cFile[lFileRead] = '\0'; - - cThisPtr = cFile; - cThisLinePtr = cThisLine; - lLineCount = 0L; - lTotalChars = 0L; - - while (*cThisPtr) { /* Read until reaching null char */ - - lIndex = 0L; /* Reset counters and flags */ - isNewline = 0; - lStartPos = lTotalChars; - - while (*cThisPtr) { - if (!isNewline) { /* Haven't read a LF yet */ - if (*cThisPtr == '\n') /* This char is LF */ - isNewline = 1; /* Set flag */ - } - else if (*cThisPtr != '\n') /* Already found LF */ - break; /* Done with line */ - - cThisLinePtr[lIndex++] = *cThisPtr++; /* Add char to output and increment */ - lTotalChars++; - } - - cThisLinePtr[lIndex] = '\0'; /* Terminate the string */ - lLineCount++; /* Increment the line counter */ - /* continue if comment or empty */ - if (cThisLinePtr[0] == '*' || cThisLinePtr[0] == '\n') { - lLineCount--; /* we count only real lines */ - continue; - } - - if (lLineCount == 1) { - cnv_get_spice_value(cThisLinePtr, &tmp); - loc->ix = ix = (int) tmp; - /* generate row data structure (x) */ - loc->xcol = (double*) calloc((size_t) ix, sizeof(double)); - } - else if (lLineCount == 2) { - cnv_get_spice_value(cThisLinePtr, &tmp); - loc->iy = iy = (int) tmp; - /* generate column data structure (y) */ - loc->ycol = (double*) calloc((size_t) iy, sizeof(double)); - } - else if (lLineCount == 3) { - char *token = CNVgettok(&cThisLinePtr); - i = 0; - while (token) { - if (i == ix) { - cm_message_printf("Too many numbers in x row."); - loc->init_err = 1; - return; - } - cnv_get_spice_value(token, &loc->xcol[i++]); - free(token); - token = CNVgettok(&cThisLinePtr); - } - if (i < ix) { - cm_message_printf("Not enough numbers in x row."); - loc->init_err = 1; - return; - } - } - else if (lLineCount == 4) { - char *token = CNVgettok(&cThisLinePtr); - i = 0; - while (token) { - if (i == iy) { - cm_message_printf("Too many numbers in y row."); - loc->init_err = 1; - return; - } - cnv_get_spice_value(token, &loc->ycol[i++]); - free(token); - token = CNVgettok(&cThisLinePtr); - } - if (i < iy) { - cm_message_printf("Not enough numbers in y row."); - loc->init_err = 1; - return; - } - /* jump out of while loop to read in the table */ - break; - } - } - - /* generate table core */ - interporder = PARAM(order); - /* boundary limits set to param 'order' aren't recognized, - so limit them here */ - if (interporder < 2) { - cm_message_printf("Parameter Order=%d not possible, set to minimum value 2", interporder); - interporder = 2; - } - /* int order : interpolation order, - int n1, int n2 : data dimensions */ - loc->newtable = sf_eno2_init(interporder, ix, iy); - - /* create table_data in memory */ - /* data [n2][n1] */ - table_data = calloc((size_t) iy, sizeof(double *)); - for (i = 0; i < iy; i++) - table_data[i] = calloc((size_t) ix, sizeof(double)); - - loc->table = table_data; - - /* continue reading from cFile */ - lLineCount = 0; - while (*cThisPtr) { /* Read until reaching null char */ - char *token; - - lIndex = 0L; /* Reset counters and flags */ - isNewline = 0; - lStartPos = lTotalChars; - - while (*cThisPtr) { /* Read until reaching null char */ - - if (!isNewline) { /* Haven't read a LF yet */ - if (*cThisPtr == '\n') /* This char is a LF */ - isNewline = 1; /* Set flag */ - } - - else if (*cThisPtr != '\n') /* Already found LF */ - break; /* Done with line */ - - cThisLinePtr[lIndex++] = *cThisPtr++; /* Add char to output and increment */ - lTotalChars++; - } - - cThisLinePtr[lIndex] = '\0'; /* Terminate the string */ - lLineCount++; /* Increment the line counter */ - /* continue if comment or empty */ - if (cThisLinePtr[0] == '*' || cThisLinePtr[0] == '\0') { - if (lTotalChars >= lFileLen) { - cm_message_printf("Not enough data in file %s", PARAM(file)); - loc->init_err = 1; - return; - } - lLineCount--; /* we count only real lines */ - continue; - } - token = CNVgettok(&cThisLinePtr); - i = 0; - while (token) { - double tmpval; - if (i == ix) { - cm_message_printf("Too many numbers in y row no. %d.", lLineCount); - loc->init_err = 1; - return; - } - - /* read table core from cFile, fill local static table structure table_data */ - cnv_get_spice_value(token, &tmpval); - table_data[lLineCount - 1][i++] = tmpval; - free(token); - token = CNVgettok(&cThisLinePtr); - } - if (i < ix) { - cm_message_printf("Not enough numbers in y row no. %d.", lLineCount); - loc->init_err = 1; - return; - } - } - - /* fill table data into eno2 structure */ - sf_eno2_set (loc->newtable, table_data /* data [n2][n1] */); - - /* free the file memory allocated */ - free(cFile); - free(cThisLine); - } /* end of initialization "if (INIT == 1)" */ - - loc = STATIC_VAR (locdata); + } /* return immediately if there was an initialization error */ - if (!loc || loc->init_err == 1) + if ((loc = STATIC_VAR(locdata)) == (Table2_Data_t *) NULL) { return; + } - /* get input x, y, - find corresponding indices, - get x and y offsets, - call interpolation function with value and derivative */ + /* get input x, y; + find corresponding indices; + get x and y offsets; + call interpolation functions with value and derivative */ xval = INPUT(inx); yval = INPUT(iny); - /* find index */ + /* check table ranges */ if (xval < loc->xcol[0] || xval > loc->xcol[loc->ix - 1]) { - if (PARAM(verbose) > 0) + if (PARAM(verbose) > 0) { cm_message_printf("x value %g exceeds table limits,\n" - " please enlarge range of your table", - xval); + " please enlarge range of your table", + xval); + } return; } if (yval < loc->ycol[0] || yval > loc->ycol[loc->iy - 1]) { - if (PARAM(verbose) > 0) + if (PARAM(verbose) > 0) { cm_message_printf("y value %g exceeds table limits,\n" - " please enlarge range of your table", - yval); + " please enlarge range of your table", + yval); + } return; } + + /*** find indices where interpolation will be done ***/ /* something like binary search to get the index */ - xind = findCrossOver(loc->xcol, 0, loc->ix - 1, xval); - /* find index with minimum distance between xval and row value */ - if (fabs(loc->xcol[xind + 1] - xval) < fabs(xval - loc->xcol[xind])) - xind++; + xind = findCrossOver(loc->xcol, loc->ix, xval); xoff = xval - loc->xcol[xind]; - yind = findCrossOver(loc->ycol, 0, loc->iy - 1, yval); - /* find index with minimum distance between yval and column value */ - if (fabs(loc->ycol[yind + 1] - yval) < fabs(yval - loc->ycol[yind])) - yind++; + yind = findCrossOver(loc->ycol, loc->iy, yval); yoff = yval - loc->ycol[yind]; - /* find local difference around index of independent row and column values */ - if (xind == loc->ix - 1) - xdiff = loc->xcol[xind] - loc->xcol[xind - 1]; - else if (xind == 0) - xdiff = loc->xcol[xind + 1] - loc->xcol[xind]; - else - xdiff = 0.5 * (loc->xcol[xind + 1] - loc->xcol[xind - 1]); - - if (yind == loc->iy - 1) - ydiff = loc->ycol[yind] - loc->ycol[yind - 1]; - else if (yind == 0) - ydiff = loc->ycol[yind + 1] - loc->ycol[yind]; - else - ydiff = 0.5 * (loc->ycol[yind + 1] - loc->ycol[yind - 1]); + /* Find local difference around index of independent row and + * column values */ + xdiff = get_local_diff(loc->ix, loc->xcol, xind); + ydiff = get_local_diff(loc->iy, loc->ycol, yind); /* Essentially non-oscillatory (ENO) interpolation to obtain the derivatives only. Using outval for now yields ngspice op non-convergence */ - sf_eno2_apply (loc->newtable, - xind, yind, /* grid location */ - xoff, yoff, /* offset from grid */ - &outval, /* output data value */ - derivval, /* output derivatives [2] */ - DER /* what to compute [FUNC, DER, BOTH] */ - ); + sf_eno2_apply(loc->newtable, + xind, yind, /* grid location */ + xoff, yoff, /* offset from grid */ + &outval, /* output data value */ + derivval, /* output derivatives [2] */ + DER /* what to compute [FUNC, DER, BOTH] */ + ); - /* bilinear interpolation to obtain the output value */ - xind = findCrossOver(loc->xcol, 0, loc->ix - 1, xval); - yind = findCrossOver(loc->ycol, 0, loc->iy - 1, yval); - x1 = loc->xcol[xind]; - x2 = loc->xcol[xind + 1]; - y1 = loc->ycol[yind]; - y2 = loc->ycol[yind + 1]; - q11 = loc->table[yind][xind]; - q12 = loc->table[yind + 1][xind]; - q21 = loc->table[yind][xind + 1]; - q22 = loc->table[yind + 1][xind + 1]; - outval = BilinearInterpolation(q11, q12, q21, q22, x1, x2, y1, y2, xval, yval); + /* xind and yind may become too large */ + if (xind == loc->ix - 1) { + --xind; + } + if (yind == loc->iy - 1) { + --yind; + } + + /* Overwrite outval from sf_eno2_apply by bilinear interpolation */ + outval = BilinearInterpolation( + xoff / (loc->xcol[xind + 1] - loc->xcol[xind]), + yoff / (loc->ycol[yind + 1] - loc->ycol[yind]), + xind, yind, loc->table); if (ANALYSIS != MIF_AC) { double xderiv, yderiv, outv; @@ -717,10 +337,14 @@ cm_table2D(ARGS) /* structure holding parms, inputs, outputs, etc. */ yderiv = PARAM(gain) * derivval[1] / ydiff; PARTIAL(out, iny) = yderiv; - if (PARAM(verbose) > 1) - cm_message_printf("\nI: %g, xval: %g, yval: %g, xderiv: %g, yderiv: %g", outv, xval, yval, xderiv, yderiv); + if (PARAM(verbose) > 1) { + cm_message_printf("\nI: %g, xval: %g, yval: %g, " + "xderiv: %g, yderiv: %g", + outv, xval, yval, xderiv, yderiv); + } } else { + Mif_Complex_t ac_gain; ac_gain.real = PARAM(gain) * derivval[0] / xdiff; ac_gain.imag= 0.0; AC_GAIN(out, inx) = ac_gain; @@ -728,5 +352,388 @@ cm_table2D(ARGS) /* structure holding parms, inputs, outputs, etc. */ ac_gain.imag= 0.0; AC_GAIN(out, iny) = ac_gain; } +} /* end of function cm_table2D */ + + + +/* This function initializes local data */ +static Table2_Data_t *init_local_data(const char *filename, int order) +{ + int xrc = 0; + int ix = 0, /* elements in a row */ + iy = 0; /* number of rows */ + + double **table_data; + double tmp; + FILE *fp = (FILE *) NULL; /* Handle to file */ + char *cFile = (char *) NULL; + char *cThisLine = (char *) NULL; + char *cThisPtr, *cThisLinePtr; + size_t lFileLen; /* Length of file */ + size_t lFileRead; /* Length of file read in */ + int lLineCount; /* Current line number */ + size_t lTotalChar; /* Total characters read */ + int interporder; /* order of interpolation for eno */ + Table2_Data_t *loc = (Table2_Data_t *) NULL; /* local data */ + + + /* Allocate static storage for *loc */ + if ((loc = (Table2_Data_t *) calloc(1, + sizeof(Table2_Data_t))) == (Table2_Data_t *) NULL) { + cm_message_printf("cannot allocate memory for lookup table."); + xrc = -1; + goto EXITPOINT; + } + + /* Init row and column counts to 0 (actually already were due + * to calloc) */ + loc->ix = loc->iy = 0; + + /* open file */ + fp = fopen_with_path(filename, "r"); + if (!fp) { /* Standard open attempt failed */ + const char * const lbuffer = getenv("NGSPICE_INPUT_DIR"); + if (lbuffer && *lbuffer) { + char * const p = (char *) malloc(strlen(lbuffer) + + strlen(DIR_PATHSEP) + + strlen(filename) + 1); + if (p == (char *) NULL) { + cm_message_printf("cannot allocate buffer to " + "attempt alternate file open."); + xrc = -1; + goto EXITPOINT; + } + (void) sprintf(p, "%s%s%s", + lbuffer, DIR_PATHSEP, filename); + fp = fopen(p, "r"); + free(p); + } + } + + /* Test for valid file pointer */ + if (!fp) { + cm_message_printf("cannot open file %s", filename); + xrc = -1; + goto EXITPOINT; + } + + /* Find the size of the data file */ + { + struct stat st; + if (fstat(fileno(fp), &st)) { + cm_message_printf("cannot get length of file %s", + filename); + xrc = -1; + goto EXITPOINT; + } + /* Copy file length */ + lFileLen = (size_t) st.st_size; + } + + /* create string to hold the whole file */ + cFile = calloc(lFileLen + 1, sizeof(char)); + /* create another string long enough for file manipulation */ + cThisLine = calloc(lFileLen + 1, sizeof(char)); + if (cFile == NULL || cThisLine == NULL) { + cm_message_printf("Insufficient memory to read file %s", + filename); + xrc = -1; + goto EXITPOINT; + } + + /* read whole file into cFile */ + { + /* Number of chars read may be less than lFileLen, because /r are + * skipped by 'fread' when file opened in text mode */ + lFileRead = fread(cFile, sizeof(char), lFileLen, fp); + const int file_error = ferror(fp); + fclose(fp); /* done with file */ + if (file_error) { + cm_message_printf("Error reading data file %s", filename); + xrc = -1; + goto EXITPOINT; + } + } + /* Number of chars read may be less than lFileLen, because /r are + * skipped by 'fread' when file opened in text mode */ + cFile[lFileRead] = '\0'; + + cThisPtr = cFile; + cThisLinePtr = cThisLine; + lLineCount = 0L; + lTotalChar = 0L; + + while (*cThisPtr) { /* Read until reaching null char */ + long lIndex = 0L; /* Index into cThisLine array */ + bool isNewline = false; /* Boolean indicating read a CR or LF */ + + while (*cThisPtr) { /* Read until reaching null char */ + if (!isNewline) { /* Haven't read a LF yet */ + if (*cThisPtr == '\n') { /* This char is a LF */ + isNewline = true; /* Set flag */ + } + } + else if (*cThisPtr != '\n') { /* Already found LF */ + break; /* Done with line */ + } + + /* Add char to output and increment */ + cThisLinePtr[lIndex++] = *cThisPtr++; + lTotalChar++; + } + + cThisLinePtr[lIndex] = '\0'; /* Terminate the string */ + lLineCount++; /* Increment the line counter */ + /* continue if comment or empty */ + if (cThisLinePtr[0] == '*' || cThisLinePtr[0] == '\n') { + lLineCount--; /* we count only real lines */ + continue; + } + + if (lLineCount == 1) { + cnv_get_spice_value(cThisLinePtr, &tmp); + loc->ix = ix = (int) tmp; + /* generate row data structure (x) */ + if ((loc->xcol = (double *) calloc((size_t) ix, + sizeof(double))) == (double *) NULL) { + cm_message_printf("Unable to allocate row structure."); + xrc = -1; + goto EXITPOINT; + } + } + else if (lLineCount == 2) { + cnv_get_spice_value(cThisLinePtr, &tmp); + loc->iy = iy = (int) tmp; + /* generate column data structure (y) */ + if ((loc->ycol = (double *) calloc((size_t) iy, + sizeof(double))) == (double *) NULL) { + cm_message_printf("Unable to allocate colum structure."); + xrc = -1; + goto EXITPOINT; + } + } + else if (lLineCount == 3) { + char *token = CNVgettok(&cThisLinePtr); + int i = 0; + while (token) { + if (i == ix) { + cm_message_printf("Too many numbers in x row."); + xrc = -1; + goto EXITPOINT; + } + cnv_get_spice_value(token, &loc->xcol[i++]); + free(token); + token = CNVgettok(&cThisLinePtr); + } + if (i < ix) { + cm_message_printf("Not enough numbers in x row."); + xrc = -1; + goto EXITPOINT; + } + } + else if (lLineCount == 4) { + char *token = CNVgettok(&cThisLinePtr); + int i = 0; + while (token) { + if (i == iy) { + cm_message_printf("Too many numbers in y row."); + xrc = -1; + goto EXITPOINT; + } + cnv_get_spice_value(token, &loc->ycol[i++]); + free(token); + token = CNVgettok(&cThisLinePtr); + } + if (i < iy) { + cm_message_printf("Not enough numbers in y row."); + xrc = -1; + goto EXITPOINT; + + } + /* jump out of while loop to read in the table */ + break; + } + } + + /* generate table core */ + interporder = order; + /* boundary limits set to param 'order' aren't recognized, + so limit them here */ + if (interporder < 2) { + cm_message_printf("Parameter Order=%d not possible, " + "set to minimum value 2", + interporder); + interporder = 2; + } + /* int order : interpolation order, + int n1, int n2 : data dimensions */ + if ((loc->newtable = sf_eno2_init( + interporder, ix, iy)) == (sf_eno2) NULL) { + cm_message_printf("eno2 initialization failure."); + xrc = -1; + goto EXITPOINT; + + } + + /* create table_data in memory */ + /* data [n2][n1] */ + if ((loc->table = table_data = (double **) calloc((size_t) iy, + sizeof(double *))) == (double **) NULL) { + cm_message_printf("Unable to allocate data table."); + free(cFile); + free(cThisLine); + free_local_data(loc); + return (Local_Data_t *) NULL; + } + + { + int i; + for (i = 0; i < iy; i++) { + if ((table_data[i] = (double *) calloc((size_t) ix, + sizeof(double))) == (double *) NULL) { + cm_message_printf("Unable to allocate data table " + "row %d", i + 1); + free(cFile); + free(cThisLine); + free_local_data(loc); + return (Local_Data_t *) NULL; + } + } + } + + loc->table = table_data; /* give to local data structure */ + + + /* continue reading f(x,y) values from cFile */ + lLineCount = 0; + while (*cThisPtr) { /* Read until reaching null char */ + char *token; + long int lIndex = 0; /* Index into cThisLine array */ + bool isNewline = 0; + + while (*cThisPtr) { /* Read until reaching null char */ + + if (!isNewline) { /* Haven't read a LF yet */ + if (*cThisPtr == '\n') /* This char is a LF */ + isNewline = 1; /* Set flag */ + } + + else if (*cThisPtr != '\n') { /* Already found LF */ + break; /* Done with line */ + } + + /* Add char to output and increment */ + cThisLinePtr[lIndex++] = *cThisPtr++; + lTotalChar++; + } + + cThisLinePtr[lIndex] = '\0'; /* Terminate the string */ + lLineCount++; /* Increment the line counter */ + /* continue if comment or empty */ + if (cThisLinePtr[0] == '*' || cThisLinePtr[0] == '\0') { + if (lTotalChar >= lFileLen) { + cm_message_printf("Not enough data in file %s", + filename); + free(cFile); + free(cThisLine); + free_local_data(loc); + return (Local_Data_t *) NULL; + } + lLineCount--; /* we count only real lines */ + continue; + } + token = CNVgettok(&cThisLinePtr); + { + int i = 0; + while (token) { + double tmpval; + if (i == ix) { + cm_message_printf("Too many numbers in y row no. %d.", + lLineCount); + xrc = -1; + goto EXITPOINT; + } + + /* read table core from cFile, fill local static table + * structure table_data */ + cnv_get_spice_value(token, &tmpval); + table_data[lLineCount - 1][i++] = tmpval; + free(token); + token = CNVgettok(&cThisLinePtr); + } + if (i < ix) { + cm_message_printf("Not enough numbers in y row no. %d.", + lLineCount); + xrc = -1; + goto EXITPOINT; + } + } + } /* end of loop over characters read from file */ + + /* fill table data into eno2 structure */ + sf_eno2_set(loc->newtable, table_data /* data [n2][n1] */); + +EXITPOINT: + /* free the file and memory allocated */ + if (cFile != (char *) NULL) { + free(cFile); + } + if (cThisLine != (char *) NULL) { + free(cThisLine); + } + if (fp != (FILE *) NULL) { + (void) fclose(fp); + } + + /* On error free any initialization that was started */ + if (xrc != 0) { + if (loc != (Table2_Data_t *) NULL) { + free_local_data(loc); + loc = (Table2_Data_t *) NULL; + } + } + return loc; +} /* end of function init_local_data */ + + + +/* Free memory allocations in Local_Data_t structure */ +static void free_local_data(Table2_Data_t *loc) +{ + if (loc == (Table2_Data_t *) NULL) { + return; + } + + /* Free data table and related values */ + if (loc->table) { + int i; + int n_y = loc->iy; + for (i = 0; i < n_y; i++) { + free(loc->table[i]); + } + + free(loc->table); + } + + free(loc->xcol); + free(loc->ycol); + sf_eno2_close(loc->newtable); + free(loc); +} /* end of function free_local_data */ + + + +/* Finds difference between column values */ +static inline double get_local_diff(int n, double *col, int ind) +{ + if (ind >= n - 1) { + return col[n - 1] - col[n - 2]; + } + if (ind <= 0) { + return col[1] - col[0]; + } + return 0.5 * (col[ind + 1] - col[ind - 1]); +} /* end of function get_local_diff */ + + -} diff --git a/src/xspice/icm/table/table3D/cfunc.mod b/src/xspice/icm/table/table3D/cfunc.mod index 9fe1cf5b1..d354cc9c3 100644 --- a/src/xspice/icm/table/table3D/cfunc.mod +++ b/src/xspice/icm/table/table3D/cfunc.mod @@ -62,7 +62,7 @@ NON-STANDARD FEATURES ===============================================================================*/ /*=== INCLUDE FILES ====================*/ - +#include #include #include #include @@ -72,13 +72,27 @@ NON-STANDARD FEATURES #include #include +#include "support/gettokens.h" #include "mada/eno2.h" #include "mada/eno3.h" -/*=== CONSTANTS ========================*/ +typedef struct { + int ix; /* size of array in x */ + int iy; /* size of array in y */ + int iz; /* size of array in z */ -#define OK 0 -#define FAIL 1 + sf_eno3 newtable; /* the table, code borrowed from madagascar project */ + + /* Input values corresponding to each index. They define the value + * in the domain at each index value */ + double *xcol; /* array of doubles in x */ + double *ycol; /* array of doubles in y */ + double *zcol; /* array of doubles in z */ + + double ***table; /* f(xi, yj, zk) */ +} Table3_Data_t; + +typedef Table3_Data_t Local_Data_t; /*=== MACROS ===========================*/ @@ -94,47 +108,21 @@ NON-STANDARD FEATURES /*=== LOCAL VARIABLES & TYPEDEFS =======*/ -struct filesource_state { - FILE *fp; - long pos; - unsigned char atend; -}; -typedef struct { - - int ix; /* size of array in x */ - int iy; /* size of array in y */ - int iz; /* size of array in z */ - - struct filesource_state *state; /* the storage array for the - filesource status. */ - - int init_err; - - sf_eno3 newtable; /* the table, code borrowed from madagascar project */ - - double *xcol; /* array of doubles in x */ - double *ycol; /* array of doubles in y */ - double *zcol; /* array of doubles in z */ - - double ***table; - -} Local_Data_t; - -/*********************/ -/* 3d geometry types */ -/*********************/ - -typedef char line_t[82]; /* A SPICE size line. <= 80 characters plus '\n\0' */ /*=== FUNCTION PROTOTYPE DEFINITIONS ===*/ -extern int findCrossOver(double arr[], int low, int high, double x); - -extern double TrilinearInterpolation(double x, double y, double z, int xind, int yind, int zind, double ***td); - extern char *CNVgettok(char **s); +extern double TrilinearInterpolation(double x, double y, double z, int xind, int yind, int zind, double ***td); +int cnv_get_spice_value(char *str, double *p_value); +extern int findCrossOver(double arr[], int n, double x); + +static void free_local_data(Table3_Data_t *loc); +static inline double get_local_diff(int n, double *col, int ind); +static Table3_Data_t *init_local_data(const char *filename, int order); + + /*============================================================================== @@ -174,156 +162,25 @@ NON-STANDARD FEATURES ==============================================================================*/ -/*=== Static CNV_get_spice_value ROUTINE =============*/ - -/* - Function takes as input a string token from a SPICE - deck and returns a floating point equivalent value. -*/ -static int -cnv_get_spice_value(char *str, /* IN - The value text e.g. 1.2K */ - double *p_value) /* OUT - The numerical value */ -{ - /* the following were "int4" devices - jpm */ - size_t len; - size_t i; - int n_matched; - line_t val_str; - char c = ' '; - char c1; - - double scale_factor; - double value; - - /* Scan the input string looking for an alpha character that is not */ - /* 'e' or 'E'. Such a character is assumed to be an engineering */ - /* suffix as defined in the Spice 2G.6 user's manual. */ - - len = strlen(str); - if (len > sizeof(val_str) - 1) - len = sizeof(val_str) - 1; - - for (i = 0; i < len; i++) { - c = str[i]; - if (isalpha_c(c) && (c != 'E') && (c != 'e')) - break; - else if (isspace_c(c)) - break; - else - val_str[i] = c; - } - val_str[i] = '\0'; - - /* Determine the scale factor */ - - if ((i >= len) || (! isalpha_c(c))) - scale_factor = 1.0; - else { - - if (isupper_c(c)) - c = tolower_c(c); - - switch (c) { - - case 't': - scale_factor = 1.0e12; - break; - - case 'g': - scale_factor = 1.0e9; - break; - - case 'k': - scale_factor = 1.0e3; - break; - - case 'u': - scale_factor = 1.0e-6; - break; - - case 'n': - scale_factor = 1.0e-9; - break; - - case 'p': - scale_factor = 1.0e-12; - break; - - case 'f': - scale_factor = 1.0e-15; - break; - - case 'm': - i++; - if (i >= len) { - scale_factor = 1.0e-3; - break; - } - c1 = str[i]; - if (!isalpha_c(c1)) { - scale_factor = 1.0e-3; - break; - } - if (islower_c(c1)) - c1 = toupper_c(c1); - if (c1 == 'E') - scale_factor = 1.0e6; - else if (c1 == 'I') - scale_factor = 25.4e-6; - else - scale_factor = 1.0e-3; - break; - - default: - scale_factor = 1.0; - } - } - - /* Convert the numeric portion to a float and multiply by the */ - /* scale factor. */ - - n_matched = sscanf(val_str, "%le", &value); - - if (n_matched < 1) { - *p_value = 0.0; - return FAIL; - } - - *p_value = value * scale_factor; - return OK; -} - -static void -cm_table3D_callback(ARGS, Mif_Callback_Reason_t reason) +static void cm_table3D_callback(ARGS, + Mif_Callback_Reason_t reason) { switch (reason) { case MIF_CB_DESTROY: { - int i, j; - Local_Data_t *loc = STATIC_VAR (locdata); - if (!loc) - break; - free(loc->state); - - for (i = 0; i < loc->iz; i++) { - for (j = 0; j < loc->iy; j++) - free(loc->table[i][j]); - free(loc->table[i]); + Table3_Data_t *loc = STATIC_VAR(locdata); + if (loc) { + free_local_data(loc); + STATIC_VAR(locdata) = loc = NULL; } - free(loc->table); - free(loc->xcol); - free(loc->ycol); - free(loc->zcol); - sf_eno3_close (loc->newtable); - free(loc); - STATIC_VAR (locdata) = loc = NULL; break; - } - } -} + } /* end of case MIF_CB_DESTROY */ + } /* end of switch over reason being called */ +} /* end of function cm_table2D_callback */ + /*============================================================================== @@ -369,8 +226,10 @@ ix iy * x row independent variables x0 x1 x2 x3 ... xix-1 -y column independent variables +* y column independent variables y0 y1 y2 y3 ... yiy-1 +* z column independent variables +z0 z1 z2 z3 ... ziz-1 * table x0y0 x1y0 x2y0 ... xix-1y0 ... @@ -384,391 +243,105 @@ x0yiy-1 x1yiy-1 x2yiy-1 ... xix-1yiy-1 /*=== CM_table3D ROUTINE ===*/ -void -cm_table3D(ARGS) /* structure holding parms, inputs, outputs, etc. */ +void cm_table3D(ARGS) /* structure holding parms, inputs, outputs, etc. */ { int size, xind, yind, zind; double xval, yval, zval, xoff, yoff, zoff, xdiff, ydiff, zdiff; double derivval[3], outval; - Local_Data_t *loc; /* Pointer to local static data, not to be included + Table3_Data_t *loc; /* Pointer to local static data, not to be included in the state vector */ - Mif_Complex_t ac_gain; size = PORT_SIZE(out); - if (INIT == 1) { - - int i, j; - int ix = 0, /* elements in a row */ - iy = 0, /* number of rows */ - iz = 0; /* number of 2D tables */ - - double ***table_data; - - double tmp; - char *cFile, *cThisPtr, *cThisLine, *cThisLinePtr; - int isNewline; /* Boolean indicating we've read a CR or LF */ - size_t lFileLen; /* Length of file */ - size_t lFileRead; /* Length of file read in */ - long lIndex; /* Index into cThisLine array */ - int lLineCount; /* Current line number */ - size_t lStartPos; /* Offset of start of current line */ - size_t lTotalChars; /* Total characters read */ - int lTableCount; /* Number of tables */ - int interporder; /* order of interpolation for eno */ - + if (INIT == 1) { /* Must do initializations */ + STATIC_VAR(locdata) = init_local_data( + PARAM(file), + PARAM(order)); CALLBACK = cm_table3D_callback; - - /* allocate static storage for *loc */ - STATIC_VAR (locdata) = calloc (1, sizeof(Local_Data_t)); - loc = STATIC_VAR (locdata); - - /* Allocate storage for internal state */ - loc->state = (struct filesource_state*) malloc(sizeof(struct filesource_state)); - loc->ix = loc->iy = loc->iz = 0; - loc->init_err = 0; - - /* open file */ - loc->state->fp = fopen_with_path(PARAM(file), "r"); - loc->state->pos = 0; - loc->state->atend = 0; - if (!loc->state->fp) { - char *lbuffer, *pp; - lbuffer = getenv("NGSPICE_INPUT_DIR"); - if (lbuffer && *lbuffer) { - pp = (char*) malloc(strlen(lbuffer) + strlen(DIR_PATHSEP) + strlen(PARAM(file)) + 1); - sprintf(pp, "%s%s%s", lbuffer, DIR_PATHSEP, PARAM(file)); - loc->state->fp = fopen(pp, "r"); - free(pp); - } - } - struct stat st; - if (!loc->state->fp || fstat(fileno(loc->state->fp), &st)) { - cm_message_printf("cannot open file %s", PARAM(file)); - free(loc->state); - free(loc); - STATIC_VAR (locdata) = loc = NULL; - return; - } - /* get file length */ - lFileLen = (size_t) st.st_size; - - /* create string to hold the whole file */ - cFile = calloc(lFileLen + 1, sizeof(char)); - /* create another string long enough for file manipulation */ - cThisLine = calloc(lFileLen + 1, sizeof(char)); - if (cFile == NULL || cThisLine == NULL) { - cm_message_printf("Insufficient memory to read file %s", PARAM(file)); - free(loc->state); - free(loc); - STATIC_VAR (locdata) = loc = NULL; - if(cFile) free(cFile); - if(cThisLine) free(cThisLine); - return; - } - /* read whole file into cFile */ - lFileRead = fread(cFile, sizeof(char), lFileLen, loc->state->fp); - fclose(loc->state->fp); - /* Number of chars read may be less than lFileLen, because /r are skipt by 'fread' */ - cFile[lFileRead] = '\0'; - - cThisPtr = cFile; - cThisLinePtr = cThisLine; - lLineCount = 0L; - lTotalChars = 0L; - - while (*cThisPtr) { /* Read until reaching null char */ - lIndex = 0L; /* Reset counters and flags */ - isNewline = 0; - lStartPos = lTotalChars; - - while (*cThisPtr) { /* Read until reaching null char */ - if (!isNewline) { /* Haven't read a LF yet */ - if (*cThisPtr == '\n') /* This char is a LF */ - isNewline = 1; /* Set flag */ - } - - else if (*cThisPtr != '\n') /* Already found LF */ - break; /* Done with line */ - - cThisLinePtr[lIndex++] = *cThisPtr++; /* Add char to output and increment */ - lTotalChars++; - } - - cThisLinePtr[lIndex] = '\0'; /* Terminate the string */ - lLineCount++; /* Increment the line counter */ - /* continue if comment or empty */ - if (cThisLinePtr[0] == '*' || cThisLinePtr[0] == '\n') { - lLineCount--; /* we count only real lines */ - continue; - } - if (lLineCount == 1) { - cnv_get_spice_value(cThisLinePtr, &tmp); - loc->ix = ix = (int) tmp; - /* generate row data structure (x) */ - loc->xcol = (double*) calloc((size_t) ix, sizeof(double)); - } else if (lLineCount == 2) { - cnv_get_spice_value(cThisLinePtr, &tmp); - loc->iy = iy = (int) tmp; - /* generate column data structure (y) */ - loc->ycol = (double*) calloc((size_t) iy, sizeof(double)); - } else if (lLineCount == 3) { - cnv_get_spice_value(cThisLinePtr, &tmp); - loc->iz = iz = (int) tmp; - /* generate column data structure (y) */ - loc->zcol = (double*) calloc((size_t) iz, sizeof(double)); - } else if (lLineCount == 4) { - char *token = CNVgettok(&cThisLinePtr); - i = 0; - while (token) { - if (i == ix) { - cm_message_printf("Too many numbers in x row."); - loc->init_err = 1; - return; - } - cnv_get_spice_value(token, &loc->xcol[i++]); - free(token); - token = CNVgettok(&cThisLinePtr); - } - if (i < ix) { - cm_message_printf("Not enough numbers in x row."); - loc->init_err = 1; - return; - } - } else if (lLineCount == 5) { - char *token = CNVgettok(&cThisLinePtr); - i = 0; - while (token) { - if (i == iy) { - cm_message_printf("Too many numbers in y row."); - loc->init_err = 1; - return; - } - cnv_get_spice_value(token, &loc->ycol[i++]); - free(token); - token = CNVgettok(&cThisLinePtr); - } - if (i < iy) { - cm_message_printf("Not enough numbers in y row."); - loc->init_err = 1; - return; - } - } else if (lLineCount == 6) { - char *token = CNVgettok(&cThisLinePtr); - i = 0; - while (token) { - if (i == iz) { - cm_message_printf("Too many numbers in z row."); - loc->init_err = 1; - return; - } - cnv_get_spice_value(token, &loc->zcol[i++]); - free(token); - token = CNVgettok(&cThisLinePtr); - } - if (i < iz) { - cm_message_printf("Not enough numbers in z row."); - loc->init_err = 1; - return; - } - /* jump out of while loop to read in the table */ - break; - } - } - - /* generate table core */ - interporder = PARAM(order); - /* boundary limits set to param 'order' aren't recognized, - so limit them here */ - if (interporder < 2) { - cm_message_printf("Parameter Order=%d not possible, set to minimum value 2", interporder); - interporder = 2; - } - /* int order : interpolation order, - int n1, int n2, int n3 : data dimensions */ - loc->newtable = sf_eno3_init(interporder, ix, iy, iz); - - /* create table_data in memory */ - /* data [n3][n2][n1] */ - table_data = calloc((size_t) iy, sizeof(double *)); - for (i = 0; i < iz; i++) { - table_data[i] = calloc((size_t) iy, sizeof(double *)); - for (j = 0; j < iy; j++) - table_data[i][j] = calloc((size_t) ix, sizeof(double)); - } - - loc->table = table_data; - - /* continue reading from cFile */ - for (lTableCount = 0; lTableCount < iz; lTableCount++) { - lLineCount = 0; - while (lLineCount < iy) { - char *token; - - lIndex = 0L; /* Reset counters and flags */ - isNewline = 0; - lStartPos = lTotalChars; - - /* read a line */ - while (*cThisPtr) { /* Read until reaching null char */ - if (!isNewline) { /* Haven't read a CR or LF yet */ - if (*cThisPtr == '\n') /* This char is LF */ - isNewline = 1; /* Set flag */ - } - - else if (*cThisPtr != '\n') /* Already found LF */ - break; /* Done with line */ - - cThisLinePtr[lIndex++] = *cThisPtr++; /* Add char to output and increment */ - lTotalChars++; - } - - cThisLinePtr[lIndex] = '\0'; /* Terminate the string */ - /* continue if comment or empty */ - if (cThisLinePtr[0] == '*' || cThisLinePtr[0] == '\0') { - if (lTotalChars >= lFileLen) { - cm_message_printf("Not enough data in file %s", PARAM(file)); - loc->init_err = 1; - return; - } - continue; - } - token = CNVgettok(&cThisLinePtr); - i = 0; - while (token) { - double tmpval; - - if (i == ix) { - cm_message_printf("Too many numbers in y row no. %d of table %d.", lLineCount, lTableCount); - loc->init_err = 1; - return; - } - - /* read table core from cFile, fill local static table structure table_data */ - cnv_get_spice_value(token, &tmpval); - - table_data[lTableCount][lLineCount][i++] = tmpval; - - free(token); - token = CNVgettok(&cThisLinePtr); - } - if (i < ix) { - cm_message_printf("Not enough numbers in y row no. %d of table %d.", lLineCount, lTableCount); - loc->init_err = 1; - return; - } - lLineCount++; - } - } - - /* fill table data into eno3 structure */ - - sf_eno3_set(loc->newtable, table_data /* data [n3][n2][n1] */); - - /* free file memory allocated */ - free(cFile); - free(cThisLine); - } /* end of initialization "if (INIT == 1)" */ - - loc = STATIC_VAR (locdata); + } /* return immediately if there was an initialization error */ - if (!loc || loc->init_err == 1) + if ((loc = STATIC_VAR(locdata)) == (Table3_Data_t *) NULL) { return; + } /* get input x, y, z; find corresponding indices; get x and y offsets; call interpolation functions with value and derivative */ - xval = INPUT(inx); yval = INPUT(iny); zval = INPUT(inz); /* check table ranges */ if (xval < loc->xcol[0] || xval > loc->xcol[loc->ix - 1]) { - if (PARAM(verbose) > 0) - cm_message_printf("x value %g exceeds table limits, \nplease enlarge range of your table", xval); + if (PARAM(verbose) > 0) { + cm_message_printf("x value %g exceeds table limits,\n" + " please enlarge range of your table", + xval); + } return; } if (yval < loc->ycol[0] || yval > loc->ycol[loc->iy - 1]) { - if (PARAM(verbose) > 0) - cm_message_printf("y value %g exceeds table limits, \nplease enlarge range of your table", yval); + if (PARAM(verbose) > 0) { + cm_message_printf("y value %g exceeds table limits,\n" + " please enlarge range of your table", + yval); + } return; } if (zval < loc->zcol[0] || zval > loc->zcol[loc->iz - 1]) { - if (PARAM(verbose) > 0) - cm_message_printf("z value %g exceeds table limits, \nplease enlarge range of your table", zval); + if (PARAM(verbose) > 0) { + cm_message_printf("z value %g exceeds table limits,\n" + " please enlarge range of your table", + zval); + } return; } - /* find index */ - /* something like binary search to get the index */ - xind = findCrossOver(loc->xcol, 0, loc->ix - 1, xval); - /* find index with minimum distance between xval and row value - if (fabs(loc->xcol[xind + 1] - xval) < fabs(xval - loc->xcol[xind])) - xind++; - */ + /*** find indices where interpolation will be done ***/ + /* something like binary search to get the index */ + xind = findCrossOver(loc->xcol, loc->ix, xval); xoff = xval - loc->xcol[xind]; - yind = findCrossOver(loc->ycol, 0, loc->iy - 1, yval); - /* find index with minimum distance between yval and column value - if (fabs(loc->ycol[yind + 1] - yval) < fabs(yval - loc->ycol[yind])) - yind++; - */ + yind = findCrossOver(loc->ycol, loc->iy, yval); yoff = yval - loc->ycol[yind]; - zind = findCrossOver(loc->zcol, 0, loc->iz - 1, zval); - /* find index with minimum distance between zval and table value - if (fabs(loc->zcol[zind + 1] - zval) < fabs(zval - loc->zcol[zind])) - zind++; - */ + zind = findCrossOver(loc->zcol, loc->iz, zval); zoff = zval - loc->zcol[zind]; - /* find local difference around index of independent row and column values */ - if (xind == loc->ix - 1) - xdiff = loc->xcol[xind] - loc->xcol[xind - 1]; - else if (xind == 0) - xdiff = loc->xcol[xind + 1] - loc->xcol[xind]; - else - xdiff = 0.5 * (loc->xcol[xind + 1] - loc->xcol[xind - 1]); - - if (yind == loc->iy - 1) - ydiff = loc->ycol[yind] - loc->ycol[yind - 1]; - else if (yind == 0) - ydiff = loc->ycol[yind + 1] - loc->ycol[yind]; - else - ydiff = 0.5 * (loc->ycol[yind + 1] - loc->ycol[yind - 1]); - - if (zind == loc->iz - 1) - zdiff = loc->zcol[zind] - loc->zcol[zind - 1]; - else if (zind == 0) - zdiff = loc->zcol[zind + 1] - loc->zcol[zind]; - else - zdiff = 0.5 * (loc->zcol[zind + 1] - loc->zcol[zind - 1]); + /* Find local difference around index of independent row and + * column values */ + xdiff = get_local_diff(loc->ix, loc->xcol, xind); + ydiff = get_local_diff(loc->iy, loc->ycol, yind); + zdiff = get_local_diff(loc->iz, loc->zcol, zind); /* Essentially non-oscillatory (ENO) interpolation to obtain the derivatives only. Using outval for now yields ngspice op non-convergence */ - sf_eno3_apply (loc->newtable, - xind, yind, zind, /* grid location */ - xoff, yoff, zoff, /* offset from grid */ - &outval, /* output data value */ - derivval, /* output derivatives [3] */ - DER /* what to compute [FUNC, DER, BOTH] */ - ); + sf_eno3_apply(loc->newtable, + xind, yind, zind, /* grid location */ + xoff, yoff, zoff, /* offset from grid */ + &outval, /* output data value */ + derivval, /* output derivatives [3] */ + DER /* what to compute [FUNC, DER, BOTH] */ + ); -/* xind yind zind may become too large */ - if (xind == loc->ix - 1) - xind--; - if (yind == loc->iy - 1) - yind--; - if (zind == loc->iz - 1) - zind--; + /* xind and yind may become too large */ + if (xind == loc->ix - 1) { + --xind; + } + if (yind == loc->iy - 1) { + --yind; + } + if (zind == loc->iz - 1) { + --zind; + } /* overwrite outval from sf_eno3_apply by trilinear interpolation */ - outval = TrilinearInterpolation(xoff / (loc->xcol[xind + 1] - loc->xcol[xind]), - yoff / (loc->ycol[yind + 1] - loc->ycol[yind]), - zoff / (loc->zcol[zind + 1] - loc->zcol[zind]), - xind, yind, zind, loc->table); + outval = TrilinearInterpolation( + xoff / (loc->xcol[xind + 1] - loc->xcol[xind]), + yoff / (loc->ycol[yind + 1] - loc->ycol[yind]), + zoff / (loc->zcol[zind + 1] - loc->zcol[zind]), + xind, yind, zind, loc->table); if (ANALYSIS != MIF_AC) { double xderiv, yderiv, zderiv, outv; @@ -781,17 +354,445 @@ cm_table3D(ARGS) /* structure holding parms, inputs, outputs, etc. */ zderiv = PARAM(gain) * derivval[2] / zdiff; PARTIAL(out, inz) = zderiv; - if (PARAM(verbose) > 1) - cm_message_printf("\nI: %g, xval: %g, yval: %g, zval: %g, xderiv: %g, yderiv: %g, zderiv: %g", outv, xval, yval, zval, xderiv, yderiv, zderiv); - } else { + if (PARAM(verbose) > 1) { + cm_message_printf("\nI: %g, xval: %g, yval: %g, zval: %g, " + "xderiv: %g, yderiv: %g, zderiv: %g", + outv, xval, yval, zval, xderiv, yderiv, zderiv); + } + } + else { + Mif_Complex_t ac_gain; ac_gain.real = PARAM(gain) * derivval[0] / xdiff; ac_gain.imag= 0.0; AC_GAIN(out, inx) = ac_gain; ac_gain.real = PARAM(gain) * derivval[1] / ydiff; ac_gain.imag= 0.0; AC_GAIN(out, iny) = ac_gain; + ac_gain.real = PARAM(gain) * derivval[2] / zdiff; + ac_gain.imag= 0.0; + AC_GAIN(out, iny) = ac_gain; } -} +} /* end of function cm_table3D */ + + + +/* This function initializes local data */ +static Table3_Data_t *init_local_data(const char *filename, int interporder) +{ + int xrc = 0; + int ix = 0, /* elements in a row */ + iy = 0, /* number of rows */ + iz = 0; /* number of 2D tables */ + double ***table_data; + double tmp; + FILE *fp = (FILE *) NULL; /* Handle to file */ + char *cFile = (char *) NULL; + char *cThisLine = (char *) NULL; + char *cThisPtr, *cThisLinePtr; + size_t lFileLen; /* Length of file */ + size_t lFileRead; /* Length of file read in */ + int lLineCount; /* Current line number */ + size_t lTotalChar; /* Total characters read */ + int lTableCount; /* Number of tables */ + Table3_Data_t *loc = (Table3_Data_t *) NULL; /* local data */ + + + /* Allocate static storage for *loc */ + if ((loc = (Table3_Data_t *) calloc(1, + sizeof(Table3_Data_t))) == (Table3_Data_t *) NULL) { + cm_message_printf("cannot allocate memory for lookup table."); + xrc = -1; + goto EXITPOINT; + } + + /* Init row and column counts to 0 (actually already were due + * to calloc) */ + loc->ix = loc->iy = loc->iz = 0; + + /* open file */ + fp = fopen_with_path(filename, "r"); + if (!fp) { /* Standard open attempt failed */ + const char * const lbuffer = getenv("NGSPICE_INPUT_DIR"); + if (lbuffer && *lbuffer) { + char * const p = (char *) malloc(strlen(lbuffer) + + strlen(DIR_PATHSEP) + + strlen(filename) + 1); + if (p == (char *) NULL) { + cm_message_printf("cannot allocate buffer to " + "attempt alternate file open."); + xrc = -1; + goto EXITPOINT; + } + (void) sprintf(p, "%s%s%s", + lbuffer, DIR_PATHSEP, filename); + fp = fopen(p, "r"); + free(p); + } + } + + /* Test for valid file pointer */ + if (!fp) { + cm_message_printf("cannot open file %s", filename); + xrc = -1; + goto EXITPOINT; + } + + /* Find the size of the data file */ + { + struct stat st; + if (fstat(fileno(fp), &st)) { + cm_message_printf("cannot get length of file %s", + filename); + xrc = -1; + goto EXITPOINT; + } + /* Copy file length */ + lFileLen = (size_t) st.st_size; + } + + /* create string to hold the whole file */ + cFile = calloc(lFileLen + 1, sizeof(char)); + /* create another string long enough for file manipulation */ + cThisLine = calloc(lFileLen + 1, sizeof(char)); + if (cFile == NULL || cThisLine == NULL) { + cm_message_printf("Insufficient memory to read file %s", + filename); + xrc = -1; + goto EXITPOINT; + } + + /* read whole file into cFile */ + { + /* Number of chars read may be less than lFileLen, because /r are + * skipped by 'fread' when file opened in text mode */ + lFileRead = fread(cFile, sizeof(char), lFileLen, fp); + const int file_error = ferror(fp); + fclose(fp); /* done with file */ + if (file_error) { + cm_message_printf("Error reading data file %s", filename); + xrc = -1; + goto EXITPOINT; + } + } + /* Number of chars read may be less than lFileLen, because /r are + * skipped by 'fread' when file opened in text mode */ + cFile[lFileRead] = '\0'; + + cThisPtr = cFile; + cThisLinePtr = cThisLine; + lLineCount = 0L; + lTotalChar = 0L; + + while (*cThisPtr) { /* Read until reaching null char */ + long lIndex = 0L; /* Index into cThisLine array */ + bool isNewline = false; /* Boolean indicating read a CR or LF */ + + while (*cThisPtr) { /* Read until reaching null char */ + if (!isNewline) { /* Haven't read a LF yet */ + if (*cThisPtr == '\n') { /* This char is a LF */ + isNewline = true; /* Set flag */ + } + } + else if (*cThisPtr != '\n') { /* Already found LF */ + break; /* Done with line */ + } + + /* Add char to output and increment */ + cThisLinePtr[lIndex++] = *cThisPtr++; + lTotalChar++; + } + + cThisLinePtr[lIndex] = '\0'; /* Terminate the string */ + lLineCount++; /* Increment the line counter */ + /* continue if comment or empty */ + if (cThisLinePtr[0] == '*' || cThisLinePtr[0] == '\n') { + lLineCount--; /* we count only real lines */ + continue; + } + + if (lLineCount == 1) { + cnv_get_spice_value(cThisLinePtr, &tmp); + loc->ix = ix = (int) tmp; + /* generate row data structure (x) */ + if ((loc->xcol = (double *) calloc((size_t) ix, + sizeof(double))) == (double *) NULL) { + cm_message_printf("Unable to allocate row structure."); + xrc = -1; + goto EXITPOINT; + } + } + else if (lLineCount == 2) { + cnv_get_spice_value(cThisLinePtr, &tmp); + loc->iy = iy = (int) tmp; + /* generate column data structure (y) */ + if ((loc->ycol = (double *) calloc((size_t) iy, + sizeof(double))) == (double *) NULL) { + cm_message_printf("Unable to allocate colum structure."); + xrc = -1; + goto EXITPOINT; + } + } + else if (lLineCount == 3) { + cnv_get_spice_value(cThisLinePtr, &tmp); + loc->iz = iz = (int) tmp; + /* generate column data structure (z) */ + if ((loc->zcol = (double *) calloc((size_t) iz, + sizeof(double))) == (double *) NULL) { + cm_message_printf("Unable to allocate \"z\" structure."); + xrc = -1; + goto EXITPOINT; + } + } + else if (lLineCount == 4) { + char *token = CNVgettok(&cThisLinePtr); + int i = 0; + while (token) { + if (i == ix) { + cm_message_printf("Too many numbers in x row."); + xrc = -1; + goto EXITPOINT; + } + cnv_get_spice_value(token, &loc->xcol[i++]); + free(token); + token = CNVgettok(&cThisLinePtr); + } + if (i < ix) { + cm_message_printf("Not enough numbers in x row."); + xrc = -1; + goto EXITPOINT; + } + } + else if (lLineCount == 5) { + char *token = CNVgettok(&cThisLinePtr); + int i = 0; + while (token) { + if (i == iy) { + cm_message_printf("Too many numbers in y row."); + xrc = -1; + goto EXITPOINT; + } + cnv_get_spice_value(token, &loc->ycol[i++]); + free(token); + token = CNVgettok(&cThisLinePtr); + } + if (i < iy) { + cm_message_printf("Not enough numbers in y row."); + xrc = -1; + goto EXITPOINT; + } + } + else if (lLineCount == 6) { + char *token = CNVgettok(&cThisLinePtr); + int i = 0; + while (token) { + if (i == iz) { + cm_message_printf("Too many numbers in z row."); + xrc = -1; + goto EXITPOINT; + } + cnv_get_spice_value(token, &loc->zcol[i++]); + free(token); + token = CNVgettok(&cThisLinePtr); + } + if (i < iz) { + cm_message_printf("Not enough numbers in z row."); + xrc = -1; + goto EXITPOINT; + } + + /* jump out of while loop to read in the table */ + break; + } + } + + /* generate table core */ + /* boundary limits set to param 'order' aren't recognized, + so limit them here */ + if (interporder < 2) { + cm_message_printf("Parameter Order=%d not possible, " + "set to minimum value 2", + interporder); + interporder = 2; + } + /* int order : interpolation order, + int n1, int n2, int n3 : data dimensions */ + if ((loc->newtable = sf_eno3_init( + interporder, ix, iy, iz)) == (sf_eno3) NULL) { + cm_message_printf("eno3 initialization failure."); + xrc = -1; + goto EXITPOINT; + } + + /* create table_data in memory */ + /* data [n3][n2][n1] */ + if ((loc->table = table_data = (double ***) calloc((size_t) iz, + sizeof(double **))) == (double ***) NULL) { + cm_message_printf("Unable to allocate data table."); + xrc = -1; + goto EXITPOINT; + } + + { + int i, j; + for (i = 0; i < iz; i++) { + if ((table_data[i] = (double **) calloc((size_t) iy, + sizeof(double *))) == (double **) NULL) { + cm_message_printf("Unable to allocate data table " + "z=%d", + i); + xrc = -1; + goto EXITPOINT; + } + for (j = 0; j < iy; j++) { + if ((table_data[i][j] = (double *) calloc((size_t) ix, + sizeof(double))) == (double *) NULL) { + cm_message_printf("Unable to allocate data table " + "z=%d y=%d", + i, j); + xrc = -1; + goto EXITPOINT; + } + } + } + } + + + /* continue reading f(x,y,z) values from cFile */ + for (lTableCount = 0; lTableCount < iz; lTableCount++) { + lLineCount = 0; + while (lLineCount < iy) { + char *token; + long int lIndex = 0; /* Index into cThisLine array */ + bool isNewline = 0; + + while (*cThisPtr) { /* Read until reaching null char */ + if (!isNewline) { /* Haven't read a CR or LF yet */ + if (*cThisPtr == '\n') /* This char is LF */ + isNewline = 1; /* Set flag */ + } + + else if (*cThisPtr != '\n') { /* Already found LF */ + break; /* Done with line */ + } + + cThisLinePtr[lIndex++] = *cThisPtr++; /* Add char to output and increment */ + lTotalChar++; + } + + cThisLinePtr[lIndex] = '\0'; /* Terminate the string */ + /* continue if comment or empty */ + if (cThisLinePtr[0] == '*' || cThisLinePtr[0] == '\0') { + if (lTotalChar >= lFileLen) { + cm_message_printf("Not enough data in file %s", + filename); + xrc = -1; + goto EXITPOINT; + } + continue; + } + token = CNVgettok(&cThisLinePtr); + { + int i = 0; + while (token) { + double tmpval; + + if (i == ix) { + cm_message_printf("Too many numbers in y row " + "no. %d of table %d.", + lLineCount, lTableCount); + xrc = -1; + goto EXITPOINT; + } + + /* read table core from cFile, fill local static table structure table_data */ + cnv_get_spice_value(token, &tmpval); + + table_data[lTableCount][lLineCount][i++] = tmpval; + + free(token); + token = CNVgettok(&cThisLinePtr); + } + if (i < ix) { + cm_message_printf("Not enough numbers in y row " + "no. %d of table %d.", + lLineCount, lTableCount); + xrc = -1; + goto EXITPOINT; + } + } + lLineCount++; + } + } /* end of loop over characters read from file */ + + /* fill table data into eno3 structure */ + sf_eno3_set(loc->newtable, table_data /* data [n3][n2][n1] */); + +EXITPOINT: + /* free the file and memory allocated */ + if (cFile != (char *) NULL) { + free(cFile); + } + if (cThisLine != (char *) NULL) { + free(cThisLine); + } + if (fp != (FILE *) NULL) { + (void) fclose(fp); + } + + /* On error free any initialization that was started */ + if (xrc != 0) { + if (loc != (Table3_Data_t *) NULL) { + free_local_data(loc); + loc = (Table3_Data_t *) NULL; + } + } + return loc; +} /* end of function init_local_data */ + + + +/* Free memory allocations in Local_Data_t structure */ +static void free_local_data(Table3_Data_t *loc) +{ + if (loc == (Table3_Data_t *) NULL) { + return; + } + + /* Free data table and related values */ + if (loc->table) { + int i, j; + int n_y = loc->iy; + int n_z = loc->iz; + for (i = 0; i < n_z; i++) { + for (j = 0; j < n_y; j++) { + free(loc->table[i][j]); + } + free(loc->table[i]); + } + + free(loc->table); + } + + free(loc->xcol); + free(loc->ycol); + free(loc->zcol); + sf_eno3_close(loc->newtable); + free(loc); +} /* end of function free_local_data */ + + + +/* Finds difference between column values */ +static inline double get_local_diff(int n, double *col, int ind) +{ + if (ind >= n - 1) { + return col[n - 1] - col[n - 2]; + } + if (ind <= 0) { + return col[1] - col[0]; + } + return 0.5 * (col[ind + 1] - col[ind - 1]); +} /* end of function get_local_diff */ /* These includes add functions from extra source code files, diff --git a/visualc/xspice/analog.vcxproj b/visualc/xspice/analog.vcxproj index 14576b480..4bbc0fd19 100644 --- a/visualc/xspice/analog.vcxproj +++ b/visualc/xspice/analog.vcxproj @@ -203,6 +203,7 @@ + @@ -291,7 +292,10 @@ + + + - \ No newline at end of file + diff --git a/visualc/xspice/digital.vcxproj b/visualc/xspice/digital.vcxproj index f48a5e999..b97b08a73 100644 --- a/visualc/xspice/digital.vcxproj +++ b/visualc/xspice/digital.vcxproj @@ -203,6 +203,7 @@ + ..\..\src\xspice\%(RelativeDir);%(AdditionalIncludeDirectories) @@ -371,7 +372,10 @@ + + + - \ No newline at end of file + diff --git a/visualc/xspice/spice2poly.vcxproj b/visualc/xspice/spice2poly.vcxproj index d72ac9cda..a784b6016 100644 --- a/visualc/xspice/spice2poly.vcxproj +++ b/visualc/xspice/spice2poly.vcxproj @@ -203,6 +203,7 @@ + @@ -210,7 +211,10 @@ + + + - \ No newline at end of file + diff --git a/visualc/xspice/table.vcxproj b/visualc/xspice/table.vcxproj index ec972ee53..d65d404c4 100644 --- a/visualc/xspice/table.vcxproj +++ b/visualc/xspice/table.vcxproj @@ -94,7 +94,7 @@ Disabled - icm\$(ProjectName);..\src\include;..\..\src\include;%(AdditionalIncludeDirectories) + icm\$(ProjectName);..\src\include;..\..\src\include;..\mada;%(AdditionalIncludeDirectories) _CRT_SECURE_NO_DEPRECATE;CIDER;%(PreprocessorDefinitions) false @@ -124,7 +124,7 @@ MaxSpeed true - icm\$(ProjectName);..\src\include;..\..\src\include;%(AdditionalIncludeDirectories) + icm\$(ProjectName);..\src\include;..\..\src\include;..\mada;%(AdditionalIncludeDirectories) _CRT_SECURE_NO_DEPRECATE;%(PreprocessorDefinitions) @@ -151,7 +151,7 @@ Disabled - icm\$(ProjectName);..\src\include;..\..\src\include;%(AdditionalIncludeDirectories) + icm\$(ProjectName);..\src\include;..\..\src\include;..\mada;%(AdditionalIncludeDirectories) _CRT_SECURE_NO_DEPRECATE;CIDER;%(PreprocessorDefinitions) false @@ -183,7 +183,7 @@ MaxSpeed true - icm\$(ProjectName);..\src\include;..\..\src\include;%(AdditionalIncludeDirectories) + icm\$(ProjectName);..\src\include;..\..\src\include;..\mada;%(AdditionalIncludeDirectories) _CRT_SECURE_NO_DEPRECATE;%(PreprocessorDefinitions) @@ -203,6 +203,7 @@ + ..\..\src\xspice\icm\$(ProjectName);..\..\src\xspice\%(RelativeDir);%(AdditionalIncludeDirectories) @@ -219,7 +220,11 @@ + + + + - \ No newline at end of file + diff --git a/visualc/xspice/xtradev.vcxproj b/visualc/xspice/xtradev.vcxproj index 31bf6689f..3ff05ff78 100644 --- a/visualc/xspice/xtradev.vcxproj +++ b/visualc/xspice/xtradev.vcxproj @@ -203,6 +203,7 @@ + ..\..\src\xspice\%(RelativeDir);%(AdditionalIncludeDirectories) @@ -271,7 +272,10 @@ + + + - \ No newline at end of file + diff --git a/visualc/xspice/xtraevt.vcxproj b/visualc/xspice/xtraevt.vcxproj index e83a4c5c6..ddba7fde6 100644 --- a/visualc/xspice/xtraevt.vcxproj +++ b/visualc/xspice/xtraevt.vcxproj @@ -204,6 +204,7 @@ call .\aux-udnfunc.bat $(ProjectName) + @@ -226,7 +227,10 @@ call .\aux-udnfunc.bat $(ProjectName) + + + - \ No newline at end of file +