dctran.sav07nov91

支持数字元件仿真的SPICE插件

/*
 * Copyright (c) 1985 Thomas L. Quarles
 */

/* subroutine to do DC TRANSIENT analysis    
        --- ONLY, unlike spice2 routine with the same name! */

#include "prefix.h"
#include <stdio.h>
#include <math.h>
/* gtri - add - wbk - Add headers */
#include "MIFtypes.h"  /* Don't know why, but this is required here! */
/* gtri - end - wbk - Add headers */
#include "TRANdefs.h"
#include "CKTdefs.h"
#include "util.h"
#include "SPerror.h"

/* gtri - add - wbk - Add headers */
#include "EVT.h"
#include "MIF.h"
#include "EVTproto.h"
#include "sim_ipc_tiein.h"
/* gtri - end - wbk - Add headers */

#include "suffix.h"
void bcopy();   /* shouldn't need this, but not defined in a sys. header file */

RCSID("DCtran.c $Revision: 1.8 $ on $Date: 91/08/27 15:01:03 $")

int
DCtran(ckt,restart) 
    CKTcircuit *ckt;
    int restart;    /* forced restart flag */
{
    int i;
    double olddelta;
    double delta;
    double new;
    double *temp;
    double startdTime;
    double startsTime;
    double startTime;
    int startIters;
    int converged;
    int firsttime;
    int error;
#ifdef SENSDEBUG
    FILE *outsen;
#endif /* SENSDEBUG */
    int save_order;
    int save;
    int save2;
    int size;
    long save_mode;
    long save1;
    static char *msg = "Timestep too small";
    IFuid timeUid;
    IFuid *nameList;
    int numNames;
#ifdef SHORTMACRO
    double mt;  /* temporary so macro call won't cross line boundry */
#endif

/* gtri - add - wbk - 12/19/90 - Add IPC stuff */
    Ipc_Boolean_t  ipc_firsttime = IPC_TRUE;
    Ipc_Boolean_t  ipc_secondtime = IPC_FALSE;
    Ipc_Boolean_t  ipc_delta_cut = IPC_FALSE;
    double         ipc_last_time = 0.0;
    double         ipc_last_delta = 0.0;
/* gtri - end - wbk - 12/19/90 - Add IPC stuff */



    if(restart || ckt->CKTtime == 0) {
        delta=MIN(ckt->CKTfinalTime/50,ckt->CKTstep)/10;
        if(ckt->CKTbreaks) FREE(ckt->CKTbreaks);
        ckt->CKTbreaks=(double *)MALLOC(2*sizeof(double));
        if(ckt->CKTbreaks == (double *)NULL) return(E_NOMEM);
        *(ckt->CKTbreaks)=0;
        *(ckt->CKTbreaks+1)=ckt->CKTfinalTime;
        ckt->CKTbreakSize=2;

/* gtri - begin - wbk - 12/19/90 - Modify setting of CKTminBreak */
/*      if(ckt->CKTminBreak==0) ckt->CKTminBreak=ckt->CKTmaxStep*5e-5; */
        /* Set to 10 times delmin for ATESSE 1 compatibity */
        if(ckt->CKTminBreak==0) ckt->CKTminBreak = 10.0 * ckt->CKTdelmin;
/* gtri - end - wbk - 12/19/90 - Modify setting of CKTminBreak */

/* gtri - add - wbk - 12/19/90 - Add IPC stuff and set anal_init and anal_type */

        /* Tell the beginPlot routine what mode we're in */
        g_ipc.anal_type = IPC_ANAL_TRAN;

        /* Tell the code models what mode we're in */
        g_mif_info.circuit.anal_type = MIF_DC;

        g_mif_info.circuit.anal_init = MIF_TRUE;

/* gtri - end - wbk */

        error = CKTnames(ckt,&numNames,&nameList);
        if(error) return(error);
        (*(SPfrontEnd->IFnewUid))((GENERIC *)ckt,&timeUid,(IFuid)NULL,
                "time", UID_OTHER, (GENERIC **)NULL);
        (*(SPfrontEnd->OUTpBeginPlot))((GENERIC *)ckt,(GENERIC*)ckt->CKTcurJob,
                ckt->CKTcurJob->JOBname,timeUid,IF_REAL,numNames,nameList,
                IF_REAL,&(((TRANan*)ckt->CKTcurJob)->TRANplot));

        ckt->CKTtime = 0;
        ckt->CKTdelta = 0;
        ckt->CKTbreak=1;
        firsttime = 1;
        save_mode = (ckt->CKTmode&MODEUIC)|MODETRANOP | MODEINITJCT;
        save_order = ckt->CKTorder;

/* gtri - begin - wbk - set a breakpoint at end of supply ramping time */
        /* must do this after CKTtime set to 0 above */
        if(ckt->enh->ramp.ramptime > 0.0)
            CKTsetBreak(ckt, ckt->enh->ramp.ramptime);
/* gtri - end - wbk - set a breakpoint at end of supply ramping time */

/* gtri - begin - wbk - Call EVTop if event-driven instances exist */
    if(ckt->evt->counts.num_insts == 0) {
        /* If no event-driven instances, do what SPICE normally does */
        converged = CKTop(ckt,
                (ckt->CKTmode & MODEUIC)|MODETRANOP| MODEINITJCT,
                (ckt->CKTmode & MODEUIC)|MODETRANOP| MODEINITFLOAT,
                ckt->CKTdcMaxIter);
        if(converged != 0) return(converged);
    }
    else {
        /* Else, use new DCOP algorithm */
        converged = EVTop(ckt,
                (ckt->CKTmode & MODEUIC) | MODETRANOP | MODEINITJCT,
                (ckt->CKTmode & MODEUIC)|MODETRANOP| MODEINITFLOAT,
                ckt->CKTdcMaxIter,
                MIF_TRUE);
        EVTop_save(ckt, MIF_FALSE, 0.0);
        if(converged != 0)
            return(converged);
    }
/* gtri - end - wbk - Call EVTop if event-driven instances exist */

/* gtri - add - wbk - 12/19/90 - Add IPC stuff */

        /* Send the operating point results for Mspice compatibility */
        if(g_ipc.enabled) {
            ipc_send_dcopb();
            CKTdump(ckt,(double)0,(((TRANan*)ckt->CKTcurJob)->TRANplot));
            ipc_send_enddcop();
        }

/* gtri - end - wbk */


/* gtri - add - wbk - 12/19/90 - set anal_init and anal_type */

        g_mif_info.circuit.anal_init = MIF_TRUE;

        /* Tell the code models what mode we're in */
        g_mif_info.circuit.anal_type = MIF_TRAN;

/* gtri - end - wbk */

/* gtri - begin - wbk - Add Breakpoint stuff */

        /* Initialize the temporary breakpoint variables to infinity */
        g_mif_info.breakpoint.current = 1.0e30;
        g_mif_info.breakpoint.last    = 1.0e30;

/* gtri - end - wbk - Add Breakpoint stuff */

        ckt->CKTstat->STATtimePts ++;
        ckt->CKTorder=1;
        for(i=0;i<7;i++) {
            ckt->CKTdeltaOld[i]=ckt->CKTmaxStep;
        }
        ckt->CKTdelta = delta;
#ifdef STEPDEBUG
        (void)printf("delta initialized to %g\n",ckt->CKTdelta);
#endif
        ckt->CKTsaveDelta = ckt->CKTfinalTime/50;

        if(ckt->CKTsenInfo && (ckt->CKTsenInfo->SENmode & TRANSEN)){
#ifdef SENSDEBUG
            printf("\nTransient Sensitivity Results\n\n");
            CKTsenPrint(ckt);
#endif /* SENSDEBUG */
            save = ckt->CKTsenInfo->SENmode;
            ckt->CKTsenInfo->SENmode = TRANSEN;
            save1 = ckt->CKTmode;
            save2 = ckt->CKTorder;
            ckt->CKTmode = save_mode;
            ckt->CKTorder = save_order;
            if(error = CKTsenDCtran(ckt)) return(error);
            ckt->CKTmode = save1;
            ckt->CKTorder = save2;
        }

        ckt->CKTmode = (ckt->CKTmode&MODEUIC)|MODETRAN | MODEINITTRAN;
        ckt->CKTag[0]=ckt->CKTag[1]=0;
        bcopy((char *)ckt->CKTstate0,(char *)ckt->CKTstate1,
                ckt->CKTnumStates*sizeof(double));

        if(ckt->CKTsenInfo && (ckt->CKTsenInfo->SENmode & TRANSEN)){
            size = SMPmatSize(ckt->CKTmatrix);
            for(i = 1; i<=size ; i++)
                *(ckt->CKTrhsOp + i) = *(ckt->CKTrhsOld + i);
        }

        startTime=(*(SPfrontEnd->IFseconds))();
        startIters = ckt->CKTstat->STATnumIter;
        startdTime = ckt->CKTstat->STATdecompTime;
        startsTime = ckt->CKTstat->STATsolveTime;
    } else {
        startTime=(*(SPfrontEnd->IFseconds))();
        startIters = ckt->CKTstat->STATnumIter;
        startdTime = ckt->CKTstat->STATdecompTime;
        startsTime = ckt->CKTstat->STATsolveTime;
        firsttime=0;
        goto resume;
    }

/* 650 */
nextTime:
    error = CKTaccept(ckt);
#ifdef STEPDEBUG
    printf("accepted at %g\n",ckt->CKTtime);
#endif /* STEPDEBUG */
    ckt->CKTstat->STATaccepted ++;
    ckt->CKTbreak=0;
    if(error)  {
        ckt->CKTcurrentAnalysis = DOING_TRAN;
        ckt->CKTstat->STATtranTime += (*(SPfrontEnd->IFseconds))()-startTime;
        ckt->CKTstat->STATtranIter += ckt->CKTstat->STATnumIter - startIters;
        ckt->CKTstat->STATtranDecompTime += ckt->CKTstat->STATdecompTime -
                startdTime;
        ckt->CKTstat->STATtranSolveTime += ckt->CKTstat->STATsolveTime -
                startsTime;
        return(error);
    }

/* gtri - modify - wbk - 12/19/90 - Send IPC stuff */

    if(g_ipc.enabled) {

        /* Test to see if delta was cut by a breakpoint, */
        /* a non-convergence, or a too large truncation error */
        if(ipc_firsttime)
            ipc_delta_cut = IPC_FALSE;
        else if(ckt->CKTtime < (ipc_last_time + (0.999 * ipc_last_delta)))
            ipc_delta_cut = IPC_TRUE;
        else
            ipc_delta_cut = IPC_FALSE;

        /* Record the data required to check for delta cuts */
        ipc_last_time = ckt->CKTtime;
        ipc_last_delta = MIN(ckt->CKTdelta, ckt->CKTmaxStep);

        /* Send results data if time since last dump is greater */
        /* than 'mintime', or if first or second timepoints, */
        /* or if delta was cut */
        if( (ckt->CKTtime >= (g_ipc.mintime + g_ipc.last_time)) ||
            ipc_firsttime || ipc_secondtime || ipc_delta_cut ) {

            ipc_send_datab(ckt->CKTtime);
            CKTdump(ckt,ckt->CKTtime,
                  (((TRANan*)ckt->CKTcurJob)->TRANplot));
            ipc_send_enddata();

            if(ipc_firsttime) {
                ipc_firsttime = IPC_FALSE;
                ipc_secondtime = IPC_TRUE;
            }
            else if(ipc_secondtime)
                ipc_secondtime = IPC_FALSE;

            g_ipc.last_time = ckt->CKTtime;
        }
    }
    else {
        /* Original code */
        if(ckt->CKTtime >= ckt->CKTinitTime) {
            CKTdump(ckt,ckt->CKTtime,
                 (((TRANan*)ckt->CKTcurJob)->TRANplot));
        }
    }

/* gtri - modify - wbk - 12/19/90 - Send IPC stuff */

/* gtri - begin - wbk - Update event queues/data for accepted timepoint */

    /* Note: this must be done AFTER sending results to SI so it can't */
    /* go next to CKTaccept() above */
    if(ckt->evt->counts.num_insts > 0)
        EVTaccept(ckt, ckt->CKTtime);

/* gtri - end - wbk - Update event queues/data for accepted timepoint */

    ckt->CKTstat->STAToldIter = ckt->CKTstat->STATnumIter;
    if(FABS(ckt->CKTtime - ckt->CKTfinalTime) < ckt->CKTminBreak) {
        /*printf(" done:  time is %g, final time is %g, and tol is %g\n",*/
        /*ckt->CKTtime,ckt->CKTfinalTime,ckt->CKTminBreak);*/
        (*(SPfrontEnd->OUTendPlot))( (((TRANan*)ckt->CKTcurJob)->TRANplot));
        ckt->CKTcurrentAnalysis = 0;
        ckt->CKTstat->STATtranTime += (*(SPfrontEnd->IFseconds))()-startTime;
        ckt->CKTstat->STATtranIter += ckt->CKTstat->STATnumIter - startIters;
        ckt->CKTstat->STATtranDecompTime += ckt->CKTstat->STATdecompTime -
                startdTime;
        ckt->CKTstat->STATtranSolveTime += ckt->CKTstat->STATsolveTime -
                startsTime;
        if(ckt->CKTsenInfo && (ckt->CKTsenInfo->SENmode & TRANSEN)){
            ckt->CKTsenInfo->SENmode = save;
#ifdef SENSDEBUG
            fclose(outsen);
#endif /* SENSDEBUG */
        }
        return(OK);
    }
    if( (*(SPfrontEnd->IFpauseTest))() ) {
        /* user requested pause... */
        ckt->CKTcurrentAnalysis = DOING_TRAN;
        ckt->CKTstat->STATtranTime += (*(SPfrontEnd->IFseconds))()-startTime;
        ckt->CKTstat->STATtranIter += ckt->CKTstat->STATnumIter - startIters;
        ckt->CKTstat->STATtranDecompTime += ckt->CKTstat->STATdecompTime -
                startdTime;
        ckt->CKTstat->STATtranSolveTime += ckt->CKTstat->STATsolveTime -
                startsTime;
        return(E_PAUSE);
    }
resume:
#ifdef STEPDEBUG
    if( (ckt->CKTdelta <= ckt->CKTfinalTime/50) && 
            (ckt->CKTdelta <= ckt->CKTmaxStep)) {
            ;
    } else {
        if(ckt->CKTfinalTime/50<ckt->CKTmaxStep) {
            (void)printf("limited by Tstop/50\n");
        } else {
            (void)printf("limited by Tmax\n");
        }
    }
#endif
    ckt->CKTdelta = 
            MIN(ckt->CKTdelta,ckt->CKTmaxStep);

/* gtri - begin - wbk - Modify Breakpoint stuff */
/* Allow delta to expand back to previous value after a temporary breakpoint */
/*  if(ckt->CKTtime == *(ckt->CKTbreaks)) {  */
    if( (ckt->CKTtime == *(ckt->CKTbreaks)) ||
        (ckt->CKTtime == g_mif_info.breakpoint.last) ) {
/* gtri - end - wbk - Modify Breakpoint stuff */

        /* first timepoint after a breakpoint - cut integration order */
        /* and limit timestep to .1 times minimum of time to next breakpoint,