diff --git a/src/spicelib/analysis/dctran.c b/src/spicelib/analysis/dctran.c index fe208f516..731a04631 100644 --- a/src/spicelib/analysis/dctran.c +++ b/src/spicelib/analysis/dctran.c @@ -5,8 +5,7 @@ Modified: 2000 AlansFixes Modified: 2023 XSPICE breakpoint fix for shared ngspice by Vyacheslav Shevchuk **********/ -/* subroutine to do DC TRANSIENT analysis - --- ONLY, unlike spice2 routine with the same name! */ +/* subroutine to do DC TRANSIENT analysis */ #include "ngspice/ngspice.h" #include "ngspice/cktdefs.h" @@ -21,15 +20,12 @@ extern struct dbcomm *dbs; #include "ngspice/ftedebug.h" #ifdef XSPICE -/* gtri - add - wbk - Add headers */ #include "ngspice/miftypes.h" - #include "ngspice/evt.h" #include "ngspice/enh.h" #include "ngspice/mif.h" #include "ngspice/evtproto.h" #include "ngspice/ipctiein.h" -/* gtri - end - wbk - Add headers */ #endif #ifdef CLUSTER @@ -100,13 +96,12 @@ DCtran(CKTcircuit *ckt, int ltra_num; CKTnode *node; #ifdef XSPICE -/* gtri - add - wbk - 12/19/90 - Add IPC stuff */ + /* 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 */ // Fix for sharedsync olddelta: When DCTran processes // either analog or XSPICE breakpoint, then it subtracts delta from @@ -154,6 +149,8 @@ DCtran(CKTcircuit *ckt, ckt->CKTtimePoints = TMALLOC(double, ckt->CKTtimeListSize); /* end LTRA code addition */ + /* Reset the breakpoint list. + Add breakpoints at 0 time and TSTOP (final) time */ if(ckt->CKTbreaks) FREE(ckt->CKTbreaks); ckt->CKTbreaks = TMALLOC(double, 2); if(ckt->CKTbreaks == NULL) return(E_NOMEM); @@ -166,16 +163,18 @@ DCtran(CKTcircuit *ckt, #endif #ifdef XSPICE -/* gtri - begin - wbk - 12/19/90 - Modify setting of CKTminBreak */ - /* 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 */ + /* Modify setting of CKTminBreak + Set to 10 times delmin (minimum step time). */ + if(ckt->CKTminBreak == 0) + ckt->CKTminBreak = 10.0 * ckt->CKTdelmin; + #else - if(ckt->CKTminBreak==0) ckt->CKTminBreak=ckt->CKTmaxStep*5e-5; + if(ckt->CKTminBreak == 0) + ckt->CKTminBreak = ckt->CKTmaxStep * 5e-5; #endif #ifdef XSPICE -/* gtri - add - wbk - 12/19/90 - Add IPC stuff and set anal_init and anal_type */ + /* 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; @@ -183,7 +182,6 @@ DCtran(CKTcircuit *ckt, g_mif_info.circuit.anal_type = MIF_DC; g_mif_info.circuit.anal_init = MIF_TRUE; -/* gtri - end - wbk */ #endif error = CKTnames(ckt,&numNames,&nameList); if(error) return(error); @@ -196,7 +194,7 @@ DCtran(CKTcircuit *ckt, tfree(nameList); if(error) return(error); - /* initialize CKTsoaCheck `warn' counters */ + /* initialize CKTsoaCheck 'warn' counters */ if (ckt->CKTsoaCheck) error = CKTsoaInit(); @@ -225,13 +223,12 @@ DCtran(CKTcircuit *ckt, } #ifdef XSPICE -/* gtri - begin - wbk - set a breakpoint at end of supply ramping time */ - /* must do this after CKTtime set to 0 above */ + /* set a breakpoint at end of supply ramping time. + Must do this after CKTtime is 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 */ + /* Call EVTop if event-driven instances exist */ if(ckt->evt->counts.num_insts != 0) { /* use new DCOP algorithm */ converged = EVTop(ckt, @@ -242,10 +239,9 @@ DCtran(CKTcircuit *ckt, EVTdump(ckt, IPC_ANAL_DCOP, 0.0); EVTop_save(ckt, MIF_FALSE, 0.0); - -/* gtri - end - wbk - Call EVTop if event-driven instances exist */ } else #endif + /* Get operating point for analogue circuit */ converged = CKTop(ckt, (ckt->CKTmode & MODEUIC) | MODETRANOP | MODEINITJCT, (ckt->CKTmode & MODEUIC) | MODETRANOP | MODEINITFLOAT, @@ -280,40 +276,31 @@ DCtran(CKTcircuit *ckt, return(converged); } #ifdef XSPICE -/* gtri - add - wbk - 12/19/90 - Add IPC stuff */ - - /* Send the operating point results for Mspice compatibility */ + /* IPC stuff: Send the operating point results */ if(g_ipc.enabled) { ipc_send_dcop_prefix(); CKTdump(ckt, 0.0, job->TRANplot); ipc_send_dcop_suffix(); } -/* gtri - end - wbk */ - -/* gtri - add - wbk - 12/19/90 - set anal_init and anal_type */ - + /* 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 */ - + /* 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 */ #endif ckt->CKTstat->STATtimePts ++; ckt->CKTorder = 1; + /* Initialze CKTdeltaOld with maximum value */ for(i=0;i<7;i++) { ckt->CKTdeltaOld[i]=ckt->CKTmaxStep; } - ckt->CKTdelta = delta; /* delta set in line 135 */ + ckt->CKTdelta = delta; /* delta set in line 130 */ #ifdef STEPDEBUG (void)printf("delta initialized to %g\n",ckt->CKTdelta); #endif @@ -339,9 +326,9 @@ DCtran(CKTcircuit *ckt, ckt->CKTorder = save2; } #endif - - ckt->CKTmode = (ckt->CKTmode&MODEUIC) | MODETRAN | MODEINITTRAN; /* modeinittran set here */ + ckt->CKTmode = (ckt->CKTmode&MODEUIC) | MODETRAN | MODEINITTRAN; + /* Reset the gear variable coefficient matrix. */ ckt->CKTag[0]=ckt->CKTag[1]=0; if (ckt->CKTstate1 && ckt->CKTstate0) { memcpy(ckt->CKTstate1, ckt->CKTstate0, @@ -360,14 +347,16 @@ DCtran(CKTcircuit *ckt, #ifdef CLUSTER CLUsetup(ckt); #endif + /* End of (restart || ckt->CKTtime == 0) */ } else { - /* saj As traninit resets CKTmode */ + /* traninit resets CKTmode */ ckt->CKTmode = (ckt->CKTmode&MODEUIC) | MODETRAN | MODEINITPRED; - /* saj */ + INIT_STATS(); - if(ckt->CKTminBreak==0) ckt->CKTminBreak=ckt->CKTmaxStep*5e-5; + if(ckt->CKTminBreak == 0) + ckt->CKTminBreak = ckt->CKTmaxStep * 5e-5; firsttime=0; - /* To get rawfile working saj*/ + /* To get rawfile working */ error = SPfrontEnd->OUTpBeginPlot (NULL, NULL, NULL, NULL, 0, @@ -377,8 +366,7 @@ DCtran(CKTcircuit *ckt, fprintf(stderr, "Couldn't relink rawfile\n"); return error; } - /* end saj*/ - goto resume; + goto resume; /* line 517 */ } /* 650 */ @@ -403,7 +391,8 @@ DCtran(CKTcircuit *ckt, error = CKTaccept(ckt); /* check if current breakpoint is outdated; if so, clear */ - if (ckt->CKTtime > ckt->CKTbreaks[0]) CKTclrBreak(ckt); + if (ckt->CKTtime > ckt->CKTbreaks[0]) + CKTclrBreak(ckt); if (ckt->CKTsoaCheck) error = CKTsoaCheck(ckt); @@ -432,8 +421,8 @@ DCtran(CKTcircuit *ckt, return(error); } #ifdef XSPICE -/* gtri - modify - wbk - 12/19/90 - Send IPC stuff */ + /* Send IPC stuff */ if ((g_ipc.enabled) || wantevtdata) { /* Send event-driven results */ @@ -477,8 +466,8 @@ DCtran(CKTcircuit *ckt, g_ipc.last_time = ckt->CKTtime; } + /* End of Send IPC stuff*/ } else -/* gtri - modify - wbk - 12/19/90 - Send IPC stuff */ #endif #ifdef CLUSTER CLUoutput(ckt); @@ -487,15 +476,14 @@ DCtran(CKTcircuit *ckt, || (!(ckt->CKTmode&MODEUIC) && ckt->CKTtime >= ckt->CKTinitTime)) CKTdump(ckt, ckt->CKTtime, job->TRANplot); #ifdef XSPICE -/* 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 */ + /* 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 */ #endif ckt->CKTstat->STAToldIter = ckt->CKTstat->STATnumIter; - /* check for the end of the tran simulation, either by< stop time given, + /* Check for the end of the tran simulation, either by stop time given, or final time has been reached. */ if (have_autostop) /* time consuming autostop check only, when variable 'autostop' has been set @@ -518,7 +506,7 @@ DCtran(CKTcircuit *ckt, if (flag_autostop) fprintf(stdout, "\nNote: Autostop after %e s, all measurement conditions are fulfilled.\n", ckt->CKTtime); - /* Final return from tran*/ + /* Final return from tran upon success */ return(OK); } if(SPfrontEnd->IFpauseTest()) { @@ -556,14 +544,12 @@ resume: #endif - /* are we at a breakpoint, or indistinguishably close? */ - /* if ((ckt->CKTtime == ckt->CKTbreaks[0]) || (ckt->CKTbreaks[0] - */ + /* Are we at a breakpoint, or indistinguishably close? */ if ( AlmostEqualUlps( ckt->CKTtime, ckt->CKTbreaks[0], 100 ) || ckt->CKTbreaks[0] - ckt->CKTtime <= ckt->CKTdelmin) { - /* first timepoint after a breakpoint - cut integration order */ - /* and limit timestep to .1 times minimum of time to next breakpoint, - * and previous timestep - */ + /* First timepoint after a breakpoint: cut integration order + and limit timestep to .1 times minimum of time to next breakpoint, + and previous timestep. */ ckt->CKTorder = 1; #ifdef STEPDEBUG if( (ckt->CKTdelta > .1*ckt->CKTsaveDelta) || @@ -612,11 +598,10 @@ resume: #ifdef XSPICE -/* gtri - begin - wbk - Add Breakpoint stuff */ - + /* More Breakpoint stuff */ if(ckt->CKTtime + ckt->CKTdelta >= g_mif_info.breakpoint.current) { - /* If next time > temporary breakpoint, force it to the breakpoint */ - /* And mark that timestep was set by temporary breakpoint */ + /* If next time > temporary breakpoint, force it to the breakpoint, + and mark that timestep was set by temporary breakpoint */ ckt->CKTsaveDelta = ckt->CKTdelta; ckt->CKTdelta = g_mif_info.breakpoint.current - ckt->CKTtime; g_mif_info.breakpoint.last = ckt->CKTtime + ckt->CKTdelta; @@ -625,12 +610,8 @@ resume: g_mif_info.breakpoint.last = 1.0e30; } -/* gtri - end - wbk - Add Breakpoint stuff */ - -/* gtri - begin - wbk - Modify Breakpoint stuff */ /* Throw out any permanent breakpoint with time <= current time or in the - * very near future, unless it the final stop break. - */ + very near future, unless it is the final stop break. */ #ifdef STEPDEBUG printf(" brk_pt: %g ckt_time: %g ckt_min_break: %g\n", ckt->CKTbreaks[0], ckt->CKTtime, ckt->CKTminBreak); @@ -654,8 +635,6 @@ resume: ckt->CKTdelta = ckt->CKTbreaks[0] - ckt->CKTtime; } -/* gtri - end - wbk - Modify Breakpoint stuff */ - #ifdef SHARED_MODULE /* Either directly go to next time step, or modify ckt->CKTdelta depending on synchronization requirements. sharedsync() returns 0. */ @@ -663,17 +642,16 @@ resume: ckt->CKTdelmin, 0, &ckt->CKTstat->STATrejected, 0); #endif -/* gtri - begin - wbk - Do event solution */ - + /* Do event solution */ if(ckt->evt->counts.num_insts > 0) { - /* if time = 0 and op_alternate was specified as false during */ - /* dcop analysis, call any changed instances to let them */ - /* post their outputs with their associated delays */ + /* If time = 0 and op_alternate was specified as false during + dcop analysis, call any changed instances to let them + post their outputs with their associated delays */ if((ckt->CKTtime == 0.0) && (! ckt->evt->options.op_alternate)) EVTiter(ckt); - /* while there are events on the queue with event time <= next */ + /* While there are events on the queue with event time <= next */ /* projected analog time, process them */ while((g_mif_info.circuit.evt_step = EVTnext_time(ckt)) <= (ckt->CKTtime + ckt->CKTdelta)) { @@ -685,14 +663,14 @@ resume: EVTdequeue(ckt, g_mif_info.circuit.evt_step); EVTiter(ckt); - /* If any instances have forced an earlier */ - /* next analog time, cut the delta */ + /* If any instances have forced an earlier + next analog time, cut the delta */ if(ckt->CKTbreaks[0] < g_mif_info.breakpoint.current) if(ckt->CKTbreaks[0] > ckt->CKTtime + ckt->CKTminBreak) g_mif_info.breakpoint.current = ckt->CKTbreaks[0]; if(g_mif_info.breakpoint.current < ckt->CKTtime + ckt->CKTdelta) { - /* Breakpoint must be > last accepted timepoint */ - /* and >= current event time */ + /* Breakpoint must be > last accepted timepoint + and >= current event time */ if(g_mif_info.breakpoint.current > ckt->CKTtime + ckt->CKTminBreak && g_mif_info.breakpoint.current >= g_mif_info.circuit.evt_step) { ckt->CKTsaveDelta = ckt->CKTdelta; @@ -704,7 +682,6 @@ resume: } /* end while next event time <= next analog time */ } /* end if there are event instances */ -/* gtri - end - wbk - Do event solution */ #else /* no XSPICE */ #ifdef CLUSTER @@ -739,14 +716,10 @@ resume: redostep = 1; #endif #ifdef XSPICE -/* gtri - add - wbk - 4/17/91 - Fix Berkeley bug */ -/* This is needed here to allow CAPask to output currents */ -/* during Transient analysis. A grep for CKTcurrentAnalysis */ -/* indicates that it should not hurt anything else ... */ - +/* This is needed here to allow CAPask to output currents + during Transient analysis. */ ckt->CKTcurrentAnalysis = DOING_TRAN; -/* gtri - end - wbk - 4/17/91 - Fix Berkeley bug */ xspice_breakpoints_processed = 0; #endif olddelta=ckt->CKTdelta; @@ -763,33 +736,25 @@ resume: save_mode = ckt->CKTmode; save_order = ckt->CKTorder; #ifdef XSPICE -/* gtri - begin - wbk - Add Breakpoint stuff */ /* Initialize temporary breakpoint to infinity */ g_mif_info.breakpoint.current = 1.0e30; -/* gtri - end - wbk - Add Breakpoint stuff */ - - -/* gtri - begin - wbk - add convergence problem reporting flags */ + /* Add convergence problem reporting flags */ /* delta is forced to equal delmin on last attempt near line 650 */ if(ckt->CKTdelta <= ckt->CKTdelmin) ckt->enh->conv_debug.last_NIiter_call = MIF_TRUE; else ckt->enh->conv_debug.last_NIiter_call = MIF_FALSE; -/* gtri - begin - wbk - add convergence problem reporting flags */ +/* Call all hybrids */ -/* gtri - begin - wbk - Call all hybrids */ - -/* gtri - begin - wbk - Set evt_step */ - + /* Set evt_step */ if(ckt->evt->counts.num_insts > 0) { g_mif_info.circuit.evt_step = ckt->CKTtime; } -/* gtri - end - wbk - Set evt_step */ #endif - + /* Central solver step */ converged = NIiter(ckt,ckt->CKTtranMaxIter); ckt->CKTstat->STATtimePts ++; @@ -805,6 +770,7 @@ resume: return(converged); } + /* If no convergence in Central solver step */ if(converged != 0) { #ifndef CLUSTER #ifndef SHARED_MODULE @@ -825,10 +791,12 @@ resume: ckt->CKTorder = 1; #ifdef XSPICE -/* gtri - begin - wbk - Add Breakpoint stuff */ - + /* Add Breakpoint stuff */ } else if(g_mif_info.breakpoint.current < ckt->CKTtime) { - /* Force backup if temporary breakpoint is < current time */ + /* Force backup if temporary breakpoint is < current time: + - save old delta + - retract time by old delta + - new delta by difference between by retracted time and breakpoint */ past_breakpoint: ckt->CKTsaveDelta = ckt->CKTdelta; @@ -841,11 +809,13 @@ resume: ckt->CKTmode = (ckt->CKTmode&MODEUIC)|MODETRAN | MODEINITTRAN; } ckt->CKTorder = 1; - -/* gtri - end - wbk - Add Breakpoint stuff */ #endif } else { + /* If converged: + - go to next time step if this was the first time. + - If not the first time step, don not accept, but check the truncation errors, + and reduce delta accordingly, thenm redo the step, to bound the error. */ if (firsttime) { #ifdef WANT_SENSE2 if(ckt->CKTsenInfo && (ckt->CKTsenInfo->SENmode & TRANSEN)){ @@ -863,15 +833,16 @@ resume: #endif firsttime = 0; #if !defined CLUSTER && !defined SHARED_MODULE - goto nextTime; /* no check on - * first time point - */ + /* no check on first time point */ + goto nextTime; /* line 373 */ #else redostep = 0; goto chkStep; #endif } newdelta = ckt->CKTdelta; + /* Scan through all devices, estimate the truncation error, + and reduce the time step, if necessary.*/ error = CKTtrunc(ckt,&newdelta); if(error) { UPDATE_STATS(DOING_TRAN); @@ -893,13 +864,12 @@ resume: EVTcall_hybrids(ckt); if (g_mif_info.breakpoint.current < ckt->CKTtime) { /* A hybrid requested a breakpoint in the past. */ - goto past_breakpoint; } } #endif - - if ((ckt->CKTorder == 1) && (ckt->CKTmaxOrder > 1)) { /* don't rise the order for backward Euler */ + /* don't raise the order for backward Euler */ + if ((ckt->CKTorder == 1) && (ckt->CKTmaxOrder > 1)) { newdelta = ckt->CKTdelta; ckt->CKTorder = 2; error = CKTtrunc(ckt, &newdelta); @@ -950,13 +920,16 @@ resume: #endif #if !defined CLUSTER && !defined SHARED_MODULE - /* go to 650 - trapezoidal */ - goto nextTime; + /* trapezoidal */ + goto nextTime; /* line 373 */ #else redostep = 0; goto chkStep; #endif } else { + /* not (newdelta > .9 * ckt->CKTdelta): reject the step + - redo the time + - apply the new (reduced) delta */ #ifndef CLUSTER #ifndef SHARED_MODULE ckt->CKTtime = ckt->CKTtime -ckt->CKTdelta; @@ -968,11 +941,13 @@ resume: ckt->CKTdelta = newdelta; #ifdef STEPDEBUG (void)printf( - "delta set to truncation error result:point rejected\n"); + "delta set to truncation error result: point rejected\n"); #endif } } - + /* Set the new delta to delmin (minimum delta allowed). However: + If the new delta has been less than the minimum delta + for the second time, bail out with 'Timestep too small'. */ if (ckt->CKTdelta <= ckt->CKTdelmin) { if (olddelta > ckt->CKTdelmin) { ckt->CKTdelta = ckt->CKTdelmin; @@ -990,8 +965,7 @@ resume: } } #ifdef XSPICE -/* gtri - begin - wbk - Do event backup */ - + /* Do event backup */ if(ckt->evt->counts.num_insts > 0) { #ifdef SHARED_MODULE double discard_start_time = ckt->CKTtime + ckt->CKTdelta; @@ -1011,7 +985,6 @@ resume: #endif } -/* gtri - end - wbk - Do event backup */ #endif #ifdef CLUSTER chkStep: