inpptree.c

spice中支持多层次元件模型仿真的可单独运行的插件源码

/**********
Copyright 1990 Regents of the University of California.  All rights reserved.
Author: 1987 Wayne A. Christopher, U. C. Berkeley CAD Group 
**********/

#include "spice.h"
#include "misc.h"
#include <stdio.h>
#include <ctype.h>
#include "strext.h"
#include "ifsim.h"
#include "iferrmsg.h"
#include "inpdefs.h"
#include "inpptree.h"
#include "util.h"
#include "suffix.h"

static IFvalue *values = NULL;
static int *types;
static int numvalues;
static GENERIC *circuit;
static INPtables *tables;

static INPparseNode *PTdifferentiate();
static INPparseNode *mkcon(), *mkb(), *mkf();
static INPparseNode *PTparse();
static INPparseNode *makepnode(), *mkbnode(), *mkfnode();
static INPparseNode *mknnode(), *mksnode();
static PTelement *PTlexer();
static int PTcheck();

extern IFsimulator *ft_sim; /* XXX */

/* Some tables that the parser uses. */

static struct op {
    int number;
    char *name;
    double (*funcptr)();
} ops[] = { 
    { PT_COMMA, ",",    NULL } ,
    { PT_PLUS,  "+",    PTplus } ,
    { PT_MINUS, "-",    PTminus } ,
    { PT_TIMES, "*",    PTtimes } ,
    { PT_DIVIDE,    "/",    PTdivide } ,
    { PT_POWER, "^",    PTpower }
} ;

#define NUM_OPS (sizeof (ops) / sizeof (struct op))

static struct func {
    char *name;
    int number;
    double (*funcptr)();
} funcs[] = {
    { "abs",    PTF_ABS,    PTabs } ,
    { "acos",   PTF_ACOS,   PTacos } ,
    { "acosh",  PTF_ACOSH,  PTacosh } ,
    { "asin",   PTF_ASIN,   PTasin } ,
    { "asinh",  PTF_ASINH,  PTasinh } ,
    { "atan",   PTF_ATAN,   PTatan } ,
    { "atanh",  PTF_ATANH,  PTatanh } ,
    { "cos",    PTF_COS,    PTcos } ,
    { "cosh",   PTF_COSH,   PTcosh } ,
    { "exp",    PTF_EXP,    PTexp } ,
    { "ln",     PTF_LN,     PTln } ,
    { "log",    PTF_LOG,    PTlog } ,
    { "sgn",    PTF_SGN,    PTsgn } ,
    { "sin",    PTF_SIN,    PTsin } ,
    { "sinh",   PTF_SINH,   PTsinh } ,
    { "sqrt",   PTF_SQRT,   PTsqrt } ,
    { "tan",    PTF_TAN,    PTtan } ,
    { "tanh",   PTF_TANH,   PTtanh } ,
    { "u",   	PTF_USTEP,  PTustep } ,
    { "uramp",  PTF_URAMP,  PTuramp } ,
    { "-",      PTF_UMINUS, PTuminus }
} ;

#define NUM_FUNCS (sizeof (funcs) / sizeof (struct func))

/* These are all the constants any sane person needs. */

static struct constant {
    char *name;
    double value;
} constants[] = {
    { "e",  M_E },
    { "pi", M_PI }
} ;

#define NUM_CONSTANTS (sizeof (constants) / sizeof (struct constant))

/* Parse the expression in *line as far as possible, and return the parse
 * tree obtained.  If there is an error, *pt will be set to NULL and an error
 * message will be printed.
 */

void
INPgetTree(line, pt, ckt, tab)
    char **line;
    INPparseTree **pt;
    GENERIC *ckt;
    INPtables *tab;
{
    INPparseNode *p;
    int i;

    values = NULL;
    types = NULL;
    numvalues = 0;

    circuit = ckt;
    tables = tab;

    p = PTparse(line);

    if (!p || !PTcheck(p)) {
        *pt = NULL;
        return;
    }

    (*pt) = (INPparseTree *) MALLOC(sizeof (INPparseTree));

    (*pt)->p.numVars = numvalues;
    (*pt)->p.varTypes = types;
    (*pt)->p.vars = values;
    (*pt)->p.IFeval = IFeval;
    (*pt)->tree = p;

    (*pt)->derivs = (INPparseNode **) 
            MALLOC(numvalues * sizeof (INPparseNode *));

    for (i = 0; i < numvalues; i++)
        (*pt)->derivs[i] = PTdifferentiate(p, i);

    return;
}

/* This routine takes the partial derivative of the parse tree with respect to
 * the i'th variable.  We try to do optimizations like getting rid of 0-valued
 * terms.
 *
 *** Note that in the interests of simplicity we share some subtrees between
 *** the function and its derivatives.  This means that you can't free the
 *** trees.
 */

static INPparseNode *
PTdifferentiate(p, varnum)
    INPparseNode *p;
    int varnum;
{
    INPparseNode *arg1, *arg2, *newp;

/* printf("differentiating: "); printTree(p); printf(" wrt var %d\n", varnum);*/

    switch (p->type) {
    case PT_CONSTANT:
        newp = mkcon((double) 0);
        break;

    case PT_VAR:
        /* Is this the variable we're differentiating wrt? */
        if (p->valueIndex == varnum)
            newp = mkcon((double) 1);
        else
            newp = mkcon((double) 0);
        break;

    case PT_PLUS:
    case PT_MINUS:
        arg1 = PTdifferentiate(p->left, varnum);
        arg2 = PTdifferentiate(p->right, varnum);
            newp = mkb(p->type, arg1, arg2);
        break;

    case PT_TIMES:
        /* d(a * b) = d(a) * b + d(b) * a */
        arg1 = PTdifferentiate(p->left, varnum);
        arg2 = PTdifferentiate(p->right, varnum);

        newp = mkb(PT_PLUS, mkb(PT_TIMES, arg1, p->right),
                    mkb(PT_TIMES, p->left, arg2));
        break;

    case PT_DIVIDE:
        /* d(a / b) = (d(a) * b - d(b) * a) / b^2 */
        arg1 = PTdifferentiate(p->left, varnum);
        arg2 = PTdifferentiate(p->right, varnum);

        newp = mkb(PT_DIVIDE, mkb(PT_MINUS, mkb(PT_TIMES, arg1,
                p->right), mkb(PT_TIMES, p->left, arg2)),
                mkb(PT_POWER, p->right, mkcon((double) 2)));
        break;

    case PT_POWER:
        /* Two cases... If the power is a constant then we're cool.
         * Otherwise we have to be tricky.
         */
        if (p->right->type == PT_CONSTANT) {
            arg1 = PTdifferentiate(p->left, varnum);

            newp = mkb(PT_TIMES, mkb(PT_TIMES,
                    mkcon(p->right->constant),
                    mkb(PT_POWER, p->left,
                    mkcon(p->right->constant - 1))),
                    arg1);
        } else {
            /* This is complicated.  f(x) ^ g(x) ->
             * exp(y(x) * ln(f(x)) ...
             */
            arg1 = PTdifferentiate(p->left, varnum);
            arg2 = PTdifferentiate(p->right, varnum);
            newp = mkb(PT_TIMES, mkf(PTF_EXP, mkb(PT_TIMES,
                    p->right, mkf(PTF_LN, p->left))),
                    mkb(PT_PLUS, mkb(PT_TIMES, p->right,
                    mkb(PT_DIVIDE, arg1, p->left)),
                    mkb(PT_TIMES, arg2,
                    mkf(PTF_LN, arg1))));

        }
        break;

    case PT_FUNCTION:
        /* Many cases.  Set arg1 to the derivative of the function,
         * and arg2 to the derivative of the argument.
         */
        switch (p->funcnum) {
        case PTF_ABS:  /* sgn(u) */
            arg1 = mkf(PTF_SGN, p->left, 0);
            break;

        case PTF_SGN:
            arg1 = mkcon((double) 0.0);
            break;

        case PTF_ACOS:  /* - 1 / sqrt(1 - u^2) */
            arg1 = mkb(PT_DIVIDE, mkcon((double) -1), mkf(PTF_SQRT,
                    mkb(PT_MINUS, mkcon((double) 1),
                    mkb(PT_POWER, p->left, mkcon((double)
                    2)))));
            break;

        case PTF_ACOSH: /* 1 / sqrt(u^2 - 1) */
            arg1 = mkb(PT_DIVIDE, mkcon((double) 1), mkf(PTF_SQRT,
                    mkb(PT_MINUS, mkb(PT_POWER, p->left,
                    mkcon((double) 2)),
                    mkcon((double) 1))));

            break;

        case PTF_ASIN:  /* 1 / sqrt(1 - u^2) */
            arg1 = mkb(PT_DIVIDE, mkcon((double) 1), mkf(PTF_SQRT,
                    mkb(PT_MINUS, mkcon((double) 1),
                    mkb(PT_POWER, p->left, mkcon((double)
                    2)))));
            break;

        case PTF_ASINH: /* 1 / sqrt(u^2 + 1) */
            arg1 = mkb(PT_DIVIDE, mkcon((double) 1), mkf(PTF_SQRT,
                    mkb(PT_PLUS, mkb(PT_POWER, p->left,
                    mkcon((double) 2)),
                    mkcon((double) 1))));
            break;

        case PTF_ATAN:  /* 1 / (1 + u^2) */
            arg1 = mkb(PT_DIVIDE, mkcon((double) 1), mkb(PT_PLUS,
                    mkb(PT_POWER, p->left, mkcon((double)
                    2)), mkcon((double) 1)));
            break;

        case PTF_ATANH: /* 1 / (1 - u^2) */
            arg1 = mkb(PT_DIVIDE, mkcon((double) 1), mkb(PT_MINUS,
                    mkcon((double) 1), mkb(PT_POWER,
                    p->left, mkcon((double) 2))));
            break;

        case PTF_COS:   /* - sin(u) */
            arg1 = mkf(PTF_UMINUS, mkf(PTF_SIN, p->left));
            break;

        case PTF_COSH:  /* sinh(u) */
            arg1 = mkf(PTF_SINH, p->left);
            break;

        case PTF_EXP:   /* exp(u) */
            arg1 = mkf(PTF_EXP, p->left, 0);
            break;

        case PTF_LN:    /* 1 / u */
            arg1 = mkb(PT_DIVIDE, mkcon((double) 1), p->left);
            break;

        case PTF_LOG:   /* log(e) / u */
            arg1 = mkb(PT_DIVIDE, mkcon((double) M_LOG10E), p->left);
            break;

        case PTF_SIN:   /* cos(u) */
            arg1 = mkf(PTF_COS, p->left);
            break;

        case PTF_SINH:  /* cosh(u) */
            arg1 = mkf(PTF_COSH, p->left);
            break;

        case PTF_SQRT:  /* 1 / (2 * sqrt(u)) */
            arg1 = mkb(PT_DIVIDE, mkcon((double) 1), mkb(PT_TIMES,
                    mkcon((double) 2), mkf(PTF_SQRT,
                    p->left)));
            break;

        case PTF_TAN:   /* 1 / (cos(u) ^ 2) */
            arg1 = mkb(PT_DIVIDE, mkcon((double) 1), mkb(PT_POWER,
                    mkf(PTF_COS, p->left), mkcon((double)
                    2)));
            break;

        case PTF_TANH:  /* 1 / (cosh(u) ^ 2) */
            arg1 = mkb(PT_DIVIDE, mkcon((double) 1), mkb(PT_POWER,
                    mkf(PTF_COSH, p->left), mkcon((double)
                    2)));
            break;

        case PTF_USTEP:
            arg1 = mkcon((double) 0.0);
            break;

        case PTF_URAMP:
	    arg1 = mkf(PTF_USTEP, p->left);
            break;

        case PTF_UMINUS:    /* - 1 ; like a constant (was 0 !) */
            arg1 = mkcon((double) - 1.0);
            break;

	default:
            fprintf(stderr, "Internal Error: bad function # %d\n",
                    p->funcnum);
            newp = NULL;
            break;
        }

        arg2 = PTdifferentiate(p->left, varnum);

        newp = mkb(PT_TIMES, arg1, arg2);

        break;

    default:
        fprintf(stderr, "Internal error: bad node type %d\n", p->type);
        newp = NULL;
        break;
    }

/* printf("result is: "); printTree(newp); printf("\n"); */
    return (newp);
}

static INPparseNode *
mkcon(value)
    double value;
{
    INPparseNode *p = (INPparseNode *) MALLOC(sizeof (INPparseNode));

    p->type = PT_CONSTANT;
    p->constant = value;

    return (p);
}

static INPparseNode *
mkb(type, left, right)
    int type;
    INPparseNode *left, *right;
{
    INPparseNode *p = (INPparseNode *) MALLOC(sizeof (INPparseNode));
    int i;

    if ((right->type == PT_CONSTANT) && (left->type == PT_CONSTANT)) {
        switch (type) {
            case PT_TIMES:
            return (mkcon(left->constant * right->constant));

            case PT_DIVIDE:
            return (mkcon(left->constant / right->constant));

            case PT_PLUS:
            return (mkcon(left->constant + right->constant));

            case PT_MINUS:
            return (mkcon(left->constant - right->constant));

            case PT_POWER:
            return (mkcon(pow(left->constant, right->constant)));
        }
    }
    switch (type) {
        case PT_TIMES:
        if ((left->type == PT_CONSTANT) && (left->constant == 0))
            return (left);
        else if ((right->type == PT_CONSTANT) && (right->constant == 0))
            return (right);
        else if ((left->type == PT_CONSTANT) && (left->constant == 1))
            return (right);
        else if ((right->type == PT_CONSTANT) && (right->constant == 1))
            return (left);
        break;
    
        case PT_DIVIDE:
        if ((left->type == PT_CONSTANT) && (left->constant == 0))
            return (left);
        else if ((right->type == PT_CONSTANT) && (right->constant == 1))
            return (left);
        break;

        case PT_PLUS:
        if ((left->type == PT_CONSTANT) && (left->constant == 0))
            return (right);
        else if ((right->type == PT_CONSTANT) && (right->constant == 0))
            return (left);
        break;

        case PT_MINUS:
        if ((right->type == PT_CONSTANT) && (right->constant == 0))
            return (left);
        else if ((left->type == PT_CONSTANT) && (left->constant == 0))
            return (mkf(PTF_UMINUS, right));
        break;
    
        case PT_POWER:
        if (right->type == PT_CONSTANT) {
            if (right->constant == 0)
                return (mkcon(1.0));
            else if (right->constant == 1)
                return (left);
        }
        break;
    }

    p->type = type;
    p->left = left;
    p->right = right;

    for (i = 0; i < NUM_OPS; i++)
        if (ops[i].number == type)
            break;
    if (i == NUM_OPS) {
        fprintf(stderr, "Internal Error: bad type %d\n", type);
        return (NULL);
    }
    p->function = ops[i].funcptr;
    p->funcname = ops[i].name;

    return (p);
}

static INPparseNode *
mkf(type, arg)
    int type;
    INPparseNode *arg;
{
    INPparseNode *p = (INPparseNode *) MALLOC(sizeof (INPparseNode));
    int i;
    double constval;

    for (i = 0; i < NUM_FUNCS; i++)
        if (funcs[i].number == type)
            break;
    if (i == NUM_FUNCS) {
        fprintf(stderr, "Internal Error: bad type %d\n", type);
        return (NULL);
    }

    if (arg->type == PT_CONSTANT) {
        constval = ((*funcs[i].funcptr) (arg->constant));
        return (mkcon(constval));
    }

    p->type = PT_FUNCTION;
    p->left = arg;

    p->funcnum = i;
    p->function = funcs[i].funcptr;
    p->funcname = funcs[i].name;

    return (p);
}

/* Check for remaining PT_PLACEHOLDERs in the parse tree.  Returns 1 if ok. */

static int
PTcheck(p)
    INPparseNode *p;
{
    switch (p->type) {
    case PT_PLACEHOLDER:
        return (0);

    case PT_CONSTANT:
    case PT_VAR:
        return (1);

    case PT_FUNCTION:
        return (PTcheck(p->left));

    case PT_PLUS:
    case PT_MINUS:
    case PT_TIMES:
    case PT_DIVIDE:
    case PT_POWER:
        return (PTcheck(p->left) && PTcheck(p->right));
    
    default:
        fprintf(stderr, "Internal error: bad node type %d\n", p->type);
        return (0);
    }
}

/* The operator-precedence table for the parser. */

#define G 1 /* Greater than. */
#define L 2 /* Less than. */
#define E 3 /* Equal. */
#define R 4 /* Error. */

static char prectable[11][11] = {
       /* $  +  -  *  /  ^  u- (  )  v  , */
/* $ */ { R, L, L, L, L, L, L, L, R, L, R },
/* + */ { G, G, G, L, L, L, L, L, G, L, G },
/* - */ { G, G, G, L, L, L, L, L, G, L, G },
/* * */ { G, G, G, G, G, L, L, L, G, L, G },
/* / */ { G, G, G, G, G, L, L, L, G, L, G },
/* ^ */ { G, G, G, G, G, L, L, L, G, L, G },
/* u-*/ { G, G, G, G, G, G, G, L, G, L, R },
/* ( */ { R, L, L, L, L, L, L, L, E, L, L },
/* ) */ { G, G, G, G, G, G, G, R, G, R, G },
/* v */ { G, G, G, G, G, G, G, G, G, R, G },
/* , */ { G, L, L, L, L, L, L, L, G, L, G }

} ;

/* Return an expr. */

static INPparseNode *
PTparse(line)
    char **line;
{
    PTelement stack[PT_STACKSIZE];
    int sp = 0, st, i;
    PTelement *top, *next;
    INPparseNode *pn, *lpn, *rpn;