1 /***************************************************************************
2 JSPICE3 adaptation of Spice3f2 - Copyright (c) Stephen R. Whiteley 1992
3 Copyright 1990 Regents of the University of California. All rights reserved.
4 Authors: 1987 Thomas L. Quarles
5 1993 Stephen R. Whiteley
6 ****************************************************************************/
7
8 /* subroutine to do TRANSIENT analysis */
9
10 #include "spice.h"
11 #include <stdio.h>
12 #include <math.h>
13 #include "srcdefs.h"
14 #include "sperror.h"
15 #include "outdata.h"
16 #include "trandefs.h"
17 #include "optdefs.h"
18 #include "tskdefs.h"
19 #include "cpdefs.h"
20 #include "util.h"
21 #include "cktext.h"
22 #include "niext.h"
23
24 #ifdef __STDC__
25 static int tran_dcoperation(CKTcircuit*,struct sOUTdata*,int);
26 static void tran_setbp(CKTcircuit*,struct sTRANint*);
27 static int tran_accept
28 (CKTcircuit*,struct sTRANint*,STATistics*,struct sOUTdata*,bool*,bool*);
29 static int tran_step(CKTcircuit*,struct sTRANint*,STATistics*);
30 static int interpolate(CKTcircuit*,struct sTRANint*,int*);
31 static int is_trunc_dev(CKTcircuit*);
32 #else
33 static int tran_dcoperation();
34 static void tran_setbp();
35 static int tran_accept();
36 static int tran_step();
37 static int interpolate();
38 static int is_trunc_dev();
39 #endif
40
41
42 int
TRANan(cktp,restart)43 TRANan(cktp,restart)
44
45 GENERIC *cktp;
46 int restart; /* forced restart flag */
47 {
48 CKTcircuit *ckt = (CKTcircuit *)cktp;
49 TRANAN *job = (TRANAN*)ckt->CKTcurJob;
50 struct sTRANint *tran = &job->TS;
51 static struct sOUTdata outd;
52 int i, error, code;
53 double delmax, tt0;
54 STATistics *stat = (STATistics*)ckt->CKTstat;
55 extern char *kw_nousertp, *kw_hitusertp, *kw_polydegree, *kw_nojjtp;
56
57 tran->startTime = (*(SPfrontEnd->IFseconds))();
58 if (restart) {
59
60 ckt->CKTmode = job->TRANmode;
61
62 if (ckt->CKTjjPresent && !(ckt->CKTmode & MODEUIC)) {
63 WARNmsg("Using initial conditions");
64 ckt->CKTmode |= MODEUIC;
65 }
66 ckt->CKTfinalTime = job->TRANfinalTime;
67 ckt->CKTstep = job->TRANstep;
68 ckt->CKTinitTime = job->TRANinitTime;
69 ckt->CKTmaxStep = job->TRANmaxStep;
70
71 if (ckt->CKTstep <= 0 ||
72 ckt->CKTinitTime < 0 ||
73 ckt->CKTfinalTime <= 0 ||
74 ckt->CKTmaxStep < 0 ||
75 ckt->CKTfinalTime <= ckt->CKTinitTime) return (E_PARMVAL);
76
77 if (ckt->CKTmaxStep == 0) {
78 if (!is_trunc_dev(ckt))
79 /* no timestep predicting devices in circuit,
80 * set to user step
81 */
82 ckt->CKTmaxStep = ckt->CKTstep;
83 else
84 /* 50 really is too small */
85 ckt->CKTmaxStep =
86 (ckt->CKTfinalTime-ckt->CKTinitTime)/50;
87 }
88 tran->delmax = .02 * (ckt->CKTfinalTime - ckt->CKTinitTime);
89 tran->delmax = (1.01)*MIN(tran->delmax,ckt->CKTmaxStep);
90
91 tran->nointerp = false;
92 tran->hitusertp = false;
93 tran->polydegree = 1;
94 cp_getvar(kw_nousertp,VT_BOOL,(char*)&tran->nointerp);
95 cp_getvar(kw_hitusertp,VT_BOOL,(char*)&tran->hitusertp);
96 cp_getvar(kw_polydegree,VT_NUM,(char*)&tran->polydegree);
97 if (job->DC.eltName[0] && tran->nointerp) {
98 WARNmsg(
99 "Note: DCsource given, \"set nousertp\" ignored.");
100 tran->nointerp = false;
101 }
102 if (tran->nointerp && tran->hitusertp) {
103 WARNmsg(
104 "Note: nousertp set, \"set hitusertp\" ignored.");
105 tran->hitusertp = false;
106 }
107
108 if (ckt->CKTjjPresent) {
109 tran->nojjtp = false;
110 cp_getvar(kw_nojjtp,VT_BOOL,(char*)&tran->nojjtp);
111 if (tran->nojjtp)
112 ckt->CKTjjPresent = false;
113 }
114 else
115 tran->nojjtp = true;
116
117 if (job->DC.eltName[0] != NULL) {
118 /* DC source was given */
119 code = CKTtypelook("Source");
120
121 for (i = 0; i <= job->DC.nestLevel; i++) {
122 SRCinstance *here;
123 here = (SRCinstance*)NULL;
124 error = CKTfndDev((GENERIC*)ckt,&code,(GENERIC**)&here,
125 job->DC.eltName[i], (GENERIC *)NULL,(IFuid)NULL);
126 if (error) {
127 (*(SPfrontEnd->IFerror))(ERR_FATAL,
128 "DCtrCurv: source %s not in circuit",
129 &(job->DC.eltName[i]));
130 return (E_NOTFOUND);
131 }
132 job->DC.elt[i] = (GENinstance*)here;
133 job->DC.vsave[i] = here->SRCdcValue;
134 }
135 }
136 error = CKTnames(ckt,&outd.numNames,&outd.dataNames);
137 if (error)
138 return (error);
139
140 (*(SPfrontEnd->IFnewUid))((GENERIC *)ckt,&outd.refName,(IFuid)NULL,
141 "time", UID_OTHER, (GENERIC **)NULL);
142
143 if (tran->nointerp)
144 outd.numPts = 1;
145 else
146 outd.numPts =
147 (ckt->CKTfinalTime - ckt->CKTinitTime)/ckt->CKTstep + 2;
148
149 outd.circuitPtr = (GENERIC *)ckt;
150 outd.analysisPtr = (GENERIC*)ckt->CKTcurJob;
151 outd.analName = ckt->CKTcurJob->JOBname;
152 outd.refType = IF_REAL;
153 outd.dataType = IF_REAL;
154 outd.plotPtr = &job->TRANplot;
155 outd.initValue = ckt->CKTinitTime;
156 outd.finalValue = ckt->CKTfinalTime;
157 outd.step = ckt->CKTstep;
158 outd.count = 0;
159 (*(SPfrontEnd->OUTbeginPlot))((GENERIC*)&outd);
160
161 FREE(outd.dataNames);
162 }
163 tran->startIters = stat->STATnumIter;
164 tran->startdTime = stat->STATdecompTime;
165 tran->startsTime = stat->STATsolveTime;
166
167 if (ckt->CKTminBreak == 0)
168 ckt->CKTminBreak = ckt->CKTmaxStep * 5e-5;
169 if (ckt->CKTdelmin == 0)
170 ckt->CKTdelmin = ckt->CKTmaxStep * 1e-9;
171 error = DCTloop(tran_dcoperation,ckt,restart,&job->DC,&outd);
172 if (error == E_PAUSE)
173 return (error);
174
175 tt0 = (*(SPfrontEnd->IFseconds))();
176 (*(SPfrontEnd->OUTendPlot)) (job->TRANplot);
177 tran->startTime += (*(SPfrontEnd->IFseconds))() - tt0;
178 return (error);
179 }
180
181
182 /* ARGSUSED */
183 static int
tran_dcoperation(ckt,outd,restart)184 tran_dcoperation(ckt,outd,restart)
185
186 CKTcircuit *ckt;
187 struct sOUTdata *outd;
188 int restart;
189 {
190 double *temp, tt;
191 STATistics *stat = (STATistics*)ckt->CKTstat;
192 struct sTRANint *tran = &((TRANAN*)ckt->CKTcurJob)->TS;
193 SRCmodel *model;
194 SRCinstance *here;
195 bool done = false;
196 bool afterpause = false;
197 int i, error = 0, type;
198
199 afterpause = !restart;
200 if (restart) {
201
202 ckt->CKTtime = 0;
203 ckt->CKTdelta = 0;
204 ckt->CKTbreak = 1;
205
206 /* set up breakpoint table */
207 if (ckt->CKTbreaks)
208 FREE(ckt->CKTbreaks);
209 ckt->CKTbreaks = (double *)MALLOC(2*sizeof(double));
210 if (!ckt->CKTbreaks)
211 return(E_NOMEM);
212 *(ckt->CKTbreaks) = 0;
213 *(ckt->CKTbreaks + 1) = ckt->CKTfinalTime;
214 ckt->CKTbreakSize = 2;
215
216 type = CKTtypelook("Source");
217 if (type >= 0) {
218 model = (SRCmodel*)ckt->CKThead[type];
219 for (; model; model = model->SRCnextModel) {
220 for (here = model->SRCinstances; here;
221 here = here->SRCnextInstance) {
222 if (here->SRCtree)
223 (*here->SRCtree->IFinit)
224 (here->SRCtree,ckt->CKTstep,ckt->CKTfinalTime);
225 }
226 }
227 }
228
229 tran->dumpit = false;
230
231 tran->tcheck = (ckt->CKTinitTime || !tran->nointerp) ?
232 ckt->CKTinitTime : ckt->CKTfinalTime;
233
234 tran->firsttime = true;
235
236 /* set initial conditions */
237 error = CKTic(ckt);
238 if (error)
239 return (error);
240
241 if (!(ckt->CKTmode & MODEUIC)) {
242 error = CKTop(ckt,
243 MODETRANOP|MODEINITJCT,
244 MODETRANOP|MODEINITFLOAT,
245 ckt->CKTdcMaxIter);
246 }
247 else {
248 error = CKTsetic(ckt);
249 if (error)
250 return (error);
251 /* niiter() used to just do this and return */
252 temp = ckt->CKTrhs;
253 ckt->CKTrhs = ckt->CKTrhsOld;
254 ckt->CKTrhsOld = temp;
255 ckt->CKTmode = MODEUIC|MODETRANOP|MODEINITJCT;
256 error = CKTload(ckt);
257 }
258 if (error)
259 return (error);
260
261 stat->STATtimePts++;
262 ckt->CKTorder = 1;
263
264 if (ckt->CKTinitDelta > 0)
265 ckt->CKTdelta = ckt->CKTinitDelta;
266 else
267 ckt->CKTdelta = .1 * MIN(tran->delmax,ckt->CKTstep);
268
269 for (i = 0; i < 7; i++) {
270 ckt->CKTdeltaOld[i] = tran->delmax;
271 }
272 ckt->CKTdeltaOld[0] = ckt->CKTdelta;
273 ckt->CKTsaveDelta = ckt->CKTdelta;
274
275 ckt->CKTmode = (ckt->CKTmode&MODEUIC)|MODETRAN | MODEINITTRAN;
276 /* modeinittran set here */
277 ckt->CKTag[0] = ckt->CKTag[1] = 0;
278 bcopy((char *)ckt->CKTstate0,(char *)ckt->CKTstate1,
279 ckt->CKTnumStates*sizeof(double));
280
281 if (ckt->CKTtime >= ckt->CKTinitTime)
282 tran->dumpit = true;
283
284 tran->tend = ckt->CKTfinalTime - ckt->CKTminBreak;
285 }
286
287 for (;;) {
288 error = tran_accept(ckt,tran,stat,outd,&done,&afterpause);
289 if (error || done) break;
290
291 stat->STAToldIter = stat->STATnumIter;
292 if (ckt->CKTdelta > tran->delmax) ckt->CKTdelta = tran->delmax;
293 tran_setbp(ckt,tran);
294
295 /* rotate delta vector */
296 for (i = 5; i >= 0; i--) {
297 ckt->CKTdeltaOld[i+1] = ckt->CKTdeltaOld[i];
298 }
299 ckt->CKTdeltaOld[0] = ckt->CKTdelta;
300
301 /* rotate states */
302 temp = ckt->CKTstates[ckt->CKTmaxOrder+1];
303 for (i = ckt->CKTmaxOrder; i >= 0; i--) {
304 ckt->CKTstates[i+1] = ckt->CKTstates[i];
305 }
306 ckt->CKTstates[0] = temp;
307
308 error = tran_step(ckt,tran,stat);
309 if (error) break;
310 }
311
312 tt = (*(SPfrontEnd->IFseconds))();
313 stat->STATtranTime += tt - tran->startTime;
314 stat->STATtranIter += stat->STATnumIter - tran->startIters;
315 stat->STATtranDecompTime += stat->STATdecompTime - tran->startdTime;
316 stat->STATtranSolveTime += stat->STATsolveTime - tran->startsTime;
317 return (error);
318 }
319
320
321 static void
tran_setbp(ckt,tran)322 tran_setbp(ckt,tran)
323
324 CKTcircuit *ckt;
325 struct sTRANint *tran;
326 {
327 double ttmp;
328
329 /*
330 * Breakpoint handling scheme:
331 * When a timepoint t is accepted (by CKTaccept), clear all previous
332 * breakpoints, because they will never be needed again.
333 *
334 * t may itself be a breakpoint, or indistinguishably close. DON'T
335 * clear t itself; recognise it as a breakpoint and act accordingly
336 *
337 * if t is not a breakpoint, limit the timestep so that the next
338 * breakpoint is not crossed
339 */
340
341 if (tran->nojjtp) {
342
343 /* are we at a breakpoint, or indistinguishably close? */
344 if ((ckt->CKTtime == *(ckt->CKTbreaks)) || (*(ckt->CKTbreaks) -
345 (ckt->CKTtime) <= ckt->CKTdelmin)) {
346 /* first timepoint after a breakpoint - cut integration order
347 * and limit timestep to .1 times minimum of time to next
348 * breakpoint, and previous timestep
349 */
350 ckt->CKTorder = 1;
351
352 ttmp = *(ckt->CKTbreaks+1) - *(ckt->CKTbreaks);
353 if (ckt->CKTsaveDelta < ttmp)
354 ttmp = ckt->CKTsaveDelta;
355 ttmp *= 0.1;
356 if (ttmp < ckt->CKTdelta)
357 ckt->CKTdelta = ttmp;
358
359 if (tran->firsttime) {
360 ckt->CKTdelta *= 0.1;
361 }
362
363 /* don't want to get below delmin for no reason */
364 ttmp = 2.0*ckt->CKTdelmin;
365 if (ckt->CKTdelta < ttmp)
366 ckt->CKTdelta = ttmp;
367 }
368 else if (ckt->CKTtime + ckt->CKTdelta >= *(ckt->CKTbreaks)) {
369 ckt->CKTsaveDelta = ckt->CKTdelta;
370 ckt->CKTdelta = *(ckt->CKTbreaks) - ckt->CKTtime;
371 ckt->CKTbreak = 1; /* why? the current pt. is not a bkpt. */
372 }
373 }
374
375 if (tran->hitusertp || tran->nointerp) {
376 if (ckt->CKTtime + ckt->CKTdelta >= tran->tcheck) {
377 if (ckt->CKTtime < tran->tcheck) {
378 ckt->CKTsaveDelta = ckt->CKTdelta;
379 ckt->CKTdelta = tran->tcheck - ckt->CKTtime;
380 tran->dumpit = true;
381 }
382 tran->tcheck += ckt->CKTstep;
383 if (tran->tcheck > ckt->CKTfinalTime || tran->nointerp)
384 tran->tcheck = ckt->CKTfinalTime;
385 }
386 }
387 }
388
389
390 static int
tran_accept(ckt,tran,stat,outd,done,afterpause)391 tran_accept(ckt,tran,stat,outd,done,afterpause)
392
393 CKTcircuit *ckt;
394 struct sTRANint *tran;
395 STATistics *stat;
396 struct sOUTdata *outd;
397 bool *done, *afterpause;
398 {
399 TRANAN *job = (TRANAN*)ckt->CKTcurJob;
400 int error;
401 double tt0;
402
403 /* afterpause if true if resuming. If we paused during
404 * interpolation, go back and finish up, otherwise just
405 * return.
406 */
407 if (*afterpause) {
408 *afterpause = false;
409 if (tran->hitusertp || tran->nointerp)
410 return (OK);
411 }
412 else {
413 error = CKTaccept(ckt);
414 /* check if current breakpoint is outdated; if so, clear */
415 if (ckt->CKTtime > *(ckt->CKTbreaks)) CKTclrBreak(ckt);
416
417 stat->STATaccepted++;
418 ckt->CKTbreak = 0;
419 if (error) {
420 ckt->CKTcurrentAnalysis = DOING_TRAN;
421 return (error);
422 }
423 }
424
425 tt0 = (*(SPfrontEnd->IFseconds))();
426 if (tran->hitusertp || tran->nointerp) {
427 if (tran->dumpit) {
428 CKTdump(ckt,ckt->CKTtime,job->TRANplot);
429 outd->count++;
430 if (!tran->nointerp)
431 tran->dumpit = false;
432 }
433 if (tran->nointerp) {
434 if (ckt->CKTtime >= tran->tend)
435 *SPfrontEnd->OUTendit = true;
436 }
437 else {
438 if (tran->tcheck >= tran->tend)
439 *SPfrontEnd->OUTendit = true;
440 }
441 }
442 else {
443 error = interpolate(ckt,tran,&outd->count);
444 if (error == E_NOCHANGE)
445 *SPfrontEnd->OUTendit = true;
446 }
447 tran->startTime += (*(SPfrontEnd->IFseconds))() - tt0;
448
449 if (*SPfrontEnd->OUTendit) {
450 *SPfrontEnd->OUTendit = false;
451 ckt->CKTcurrentAnalysis = 0;
452 *done = true;
453 return (OK);
454 }
455
456 if ( error == E_PAUSE || (*(SPfrontEnd->IFpauseTest))() ) {
457 /* user requested pause... */
458 ckt->CKTcurrentAnalysis = DOING_TRAN;
459 return (E_PAUSE);
460 }
461 return (OK);
462 }
463
464
465 static int
tran_step(ckt,tran,stat)466 tran_step(ckt,tran,stat)
467
468 CKTcircuit *ckt;
469 struct sTRANint *tran;
470 STATistics *stat;
471 {
472 double olddelta, new;
473 int i, error, converged;
474
475 for (;;) {
476
477 olddelta = ckt->CKTdelta;
478 ckt->CKTtime += ckt->CKTdelta;
479 ckt->CKTdeltaOld[0] = ckt->CKTdelta;
480 NIcomCof(ckt);
481
482 (void)DEVpredNew(ckt);
483
484 converged = NIiter(ckt,ckt->CKTtranMaxIter);
485 stat->STATtimePts++;
486
487 ckt->CKTmode = (ckt->CKTmode&MODEUIC)|MODETRAN | MODEINITPRED;
488
489 if (tran->firsttime) {
490 for (i = 0; i < ckt->CKTnumStates; i++) {
491 *(ckt->CKTstate2+i) = *(ckt->CKTstate1+i);
492 *(ckt->CKTstate3+i) = *(ckt->CKTstate1+i);
493 }
494 }
495 if (converged != 0) {
496 ckt->CKTtime -= ckt->CKTdelta;
497 stat->STATrejected++;
498 ckt->CKTdelta *= 0.125;
499
500 if (tran->firsttime) {
501 ckt->CKTmode =
502 (ckt->CKTmode&MODEUIC)|MODETRAN | MODEINITTRAN;
503 }
504 ckt->CKTorder = 1;
505 }
506 else {
507 if (tran->firsttime) {
508 tran->firsttime = false;
509 /* no check on first time point */
510 break;
511 }
512
513 if (!tran->nojjtp) {
514 ckt->CKTorder++;
515 if (ckt->CKTorder > ckt->CKTmaxOrder)
516 ckt->CKTorder = ckt->CKTmaxOrder;
517 break;
518 }
519
520 new = ckt->CKTdelta;
521 error = CKTtrunc(ckt,&new);
522 if (error) {
523 ckt->CKTcurrentAnalysis = DOING_TRAN;
524 return (error);
525 }
526 /* spice2g had .5 here, why change?
527 if (new > .9*ckt->CKTdelta) {
528 */
529 if (new > .5*ckt->CKTdelta) {
530 if (ckt->CKTorder == 1) {
531 new = ckt->CKTdelta;
532 ckt->CKTorder = 2;
533 error = CKTtrunc(ckt,&new);
534 if (error) {
535 ckt->CKTcurrentAnalysis = DOING_TRAN;
536 return (error);
537 }
538 if (new <= 1.05*ckt->CKTdelta) {
539 ckt->CKTorder = 1;
540 }
541 }
542 /* time point OK */
543 ckt->CKTdelta = new;
544 break;
545 }
546 else {
547 /*
548 printf("rejected: time=%g delta=%g new=%g\n",ckt->CKTtime,ckt->CKTdelta,new);
549 */
550 ckt->CKTtime -= ckt->CKTdelta;
551 stat->STATrejected++;
552 ckt->CKTdelta = new;
553 }
554 }
555 if (ckt->CKTdelta <= ckt->CKTdelmin) {
556 if (olddelta > ckt->CKTdelmin) {
557 ckt->CKTdelta = ckt->CKTdelmin;
558 }
559 else {
560 ckt->CKTcurrentAnalysis = DOING_TRAN;
561 errMsg = CKTtrouble(ckt, "Timestep too small");
562 return (E_TIMESTEP);
563
564 }
565 }
566 }
567 return (OK);
568 }
569
570
571 static int
interpolate(ckt,tran,count)572 interpolate(ckt,tran,count)
573
574 /* Previously, all raw timepoint data was saved in the output data
575 * structures, and interpolation was performed when the analysis was
576 * complete. This can take huge amounts of storage if Josephson junctions
577 * are present. Here, we (by default) save only data at user time points,
578 * performing interpolation "on the fly".
579 */
580 CKTcircuit *ckt;
581 struct sTRANint *tran;
582 int *count;
583 {
584 TRANAN *job = (TRANAN*)ckt->CKTcurJob;
585 double d1,d2,d3,d4,diff[4];
586 double c2,c3,c4;
587 double finalmax,finalmin;
588 double time, delta;
589 double *temp;
590 int i,size;
591
592 /* use CKTrhs to store the interpolated values, since it is not
593 * used for anything - it will be cleared before the next time
594 * point.
595 */
596
597 if (tran->tcheck > ckt->CKTtime)
598 return (OK);
599 size = spGetSize(ckt->CKTmatrix,1);
600 delta = .5*ckt->CKTdelmin;
601 finalmax = ckt->CKTfinalTime + delta;
602 finalmin = ckt->CKTfinalTime - delta;
603 time = ckt->CKTtime + delta;
604
605 if (tran->polydegree <= 1) {
606
607 c2 = ckt->CKTdeltaOld[0];
608
609 for (; tran->tcheck <= time; tran->tcheck += ckt->CKTstep) {
610 if ( (*(SPfrontEnd->IFpauseTest))() )
611 return(E_PAUSE);
612
613 d1 = ckt->CKTtime - tran->tcheck;
614 d2 = d1 - c2;
615
616 /* first order extrapolation factors */
617 diff[0] = d2/(-c2);
618 diff[1] = d1/(c2);
619
620 for (i = 1; i <= size; i++) {
621 ckt->CKTrhs[i] =
622 diff[0]*ckt->CKTsols[0][i] +
623 diff[1]*ckt->CKTsols[1][i];
624 }
625 if (tran->tcheck < finalmax) {
626 temp = ckt->CKTrhsOld;
627 ckt->CKTrhsOld = ckt->CKTrhs;
628 CKTdump(ckt,tran->tcheck,job->TRANplot);
629 ckt->CKTrhsOld = temp;
630 (*count)++;
631 if (finalmin < tran->tcheck)
632 return (E_NOCHANGE);
633 }
634 else
635 return (E_NOCHANGE);
636 }
637 }
638 else if (tran->polydegree == 2) {
639
640 c2 = ckt->CKTdeltaOld[0];
641 c3 = c2 + ckt->CKTdeltaOld[1];
642
643 for (; tran->tcheck <= time; tran->tcheck += ckt->CKTstep) {
644 if ( (*(SPfrontEnd->IFpauseTest))() )
645 return(E_PAUSE);
646
647 d1 = ckt->CKTtime - tran->tcheck;
648 d2 = d1 - c2;
649 d3 = d1 - c3;
650
651 /* second order extrapolation factors */
652 diff[0] = d2*d3/((c2)*(c3));
653 diff[1] = d1*d3/((c2)*(c2-c3));
654 diff[2] = d1*d2/((c3)*(c3-c2));
655
656 for (i = 1; i <= size; i++) {
657 ckt->CKTrhs[i] =
658 diff[0]*ckt->CKTsols[0][i] +
659 diff[1]*ckt->CKTsols[1][i] +
660 diff[2]*ckt->CKTsols[2][i];
661 }
662 if (tran->tcheck < finalmax) {
663 temp = ckt->CKTrhsOld;
664 ckt->CKTrhsOld = ckt->CKTrhs;
665 CKTdump(ckt,tran->tcheck,job->TRANplot);
666 ckt->CKTrhsOld = temp;
667 (*count)++;
668 if (finalmin < tran->tcheck)
669 return (E_NOCHANGE);
670 }
671 else
672 return (E_NOCHANGE);
673 }
674 }
675 else if (tran->polydegree >= 3) {
676
677 c2 = ckt->CKTdeltaOld[0];
678 c3 = c2 + ckt->CKTdeltaOld[1];
679 c4 = c3 + ckt->CKTdeltaOld[2];
680
681 for (; tran->tcheck <= time; tran->tcheck += ckt->CKTstep) {
682 if ( (*(SPfrontEnd->IFpauseTest))() )
683 return(E_PAUSE);
684
685 d1 = ckt->CKTtime - tran->tcheck;
686 d2 = d1 - c2;
687 d3 = d1 - c3;
688 d4 = d1 - c4;
689
690 /* third order extrapolation factors */
691 diff[0] = d2*d3*d4/((-c2)*(c3)*(c4));
692 diff[1] = d1*d3*d4/((c2)*(c2-c3)*(c2-c4));
693 diff[2] = d1*d2*d4/((c3)*(c3-c2)*(c3-c4));
694 diff[3] = d1*d2*d3/((c4)*(c4-c2)*(c4-c3));
695
696 for (i = 1; i <= size; i++) {
697 ckt->CKTrhs[i] =
698 diff[0]*ckt->CKTsols[0][i] +
699 diff[1]*ckt->CKTsols[1][i] +
700 diff[2]*ckt->CKTsols[2][i] +
701 diff[3]*ckt->CKTsols[3][i];
702 }
703 if (tran->tcheck < finalmax) {
704 temp = ckt->CKTrhsOld;
705 ckt->CKTrhsOld = ckt->CKTrhs;
706 CKTdump(ckt,tran->tcheck,job->TRANplot);
707 ckt->CKTrhsOld = temp;
708 (*count)++;
709 if (finalmin < tran->tcheck)
710 return (E_NOCHANGE);
711 }
712 else
713 return (E_NOCHANGE);
714 }
715 }
716 if (tran->tcheck > finalmin)
717 return (E_NOCHANGE);
718 return (OK);
719 }
720
721 static int
is_trunc_dev(ckt)722 is_trunc_dev(ckt)
723 CKTcircuit *ckt;
724 {
725 extern SPICEdev *DEVices[];
726 struct sCKTmodHead *mh;
727
728 for (mh = ckt->CKTheadList; mh != NULL; mh = mh->next) {
729 if (DEVices[mh->type]->DEVtrunc != NULL)
730 return (TRUE);
731 }
732 return (FALSE);
733 }
734