/* terr.f -- translated by f2c (version of 3 February 1990  3:36:42).
   You must link the resulting object file with the libraries:
	-lF77 -lI77 -lm -lc   (in that order)
*/

#include "f2c.h"

/* Common Block Declarations */

struct {
    integer ielmnt, isbckt, nsbckt, iunsat, nunsat, itemps, numtem, isens, 
	    nsens, ifour, nfour, ifield, icode, idelim, icolum, insize, 
	    junode, lsbkpt, numbkp, iorder, jmnode, iur, iuc, ilc, ilr, 
	    numoff, isr, nmoffc, iseq, iseq1, neqn, nodevs, ndiag, iswap, 
	    iequa, macins, lvnim1, lx0, lvn, lynl, lyu, lyl, lx1, lx2, lx3, 
	    lx4, lx5, lx6, lx7, ld0, ld1, ltd, imynl, imvn, lcvn, nsnod, 
	    nsmat, nsval, icnod, icmat, icval, loutpt, lpol, lzer, irswpf, 
	    irswpr, icswpf, icswpr, irpt, jcpt, irowno, jcolno, nttbr, nttar, 
	    lvntmp;
} tabinf_;

#define tabinf_1 tabinf_

struct {
    doublereal omega, time, delta, delold[7], ag[7], vt, xni, egfet, xmu, 
	    sfactr;
    integer mode, modedc, icalc, initf, method, iord, maxord, noncon, iterno, 
	    itemno, nosolv, modac, ipiv, ivmflg, ipostp, iscrch, iofile;
} status_;

#define status_1 status_

struct {
    doublereal twopi, xlog2, xlog10, root2, rad, boltz, charge, ctok, gmin, 
	    reltol, abstol, vntol, trtol, chgtol, eps0, epssil, epsox, pivtol,
	     pivrel;
} knstnt_;

#define knstnt_1 knstnt_

struct {
    doublereal value[200000];
} blank_;

#define blank_1 blank_

/*<       subroutine terr(loct,delnew) >*/
/* Subroutine */ int terr_(loct, delnew)
integer *loct;
doublereal *delnew;
{
    /* Initialized data */

    static doublereal coef[6] = { .5,.2222222222,.1363636364,.096,
	    .07299270073,.0583090379 };
    static doublereal xtwelv = .08333333333;

    /* System generated locals */
    integer i_1;
    doublereal d_1, d_2, d_3, d_4;

    /* Builtin functions */
    double sqrt(), log(), exp();

    /* Local variables */
#define ccap ((doublereal *)&blank_1 + 1)
    static doublereal diff[8];
#define qcap ((doublereal *)&blank_1)
    static doublereal ctol;
    static integer i;
    static doublereal const_;
    static integer istop;
#define nodplc ((integer *)&blank_1)
#define cvalue ((complex *)&blank_1)
    static doublereal deltmp[7], del, tol;

/*<       implicit double precision (a-h,o-z) >*/

/*     this routine estimates the local truncation error for a particular 
*/
/* circuit element.  it then computes the appropriate stepsize which */
/* should be used. */

/* spice version 2g.6  sccsid=tabinf 3/15/83 */
/*<       common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem, >*/
/*<      1   isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize, >*/
/*<      2   junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr, >*/
/*<      3   nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1, >*/
/*<      4   lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd, >*/
/*<      5   imynl,imvn,lcvn,nsnod,nsmat,nsval,icnod,icmat,icval, >*/
/*<      6   loutpt,lpol,lzer,irswpf,irswpr,icswpf,icswpr,irpt,jcpt, >*/
/*<      7   irowno,jcolno,nttbr,nttar,lvntmp >*/
/* spice version 2g.6  sccsid=status 3/15/83 */
/*<       common /status/ omega,time,delta,delold(7),ag(7),vt,xni,egfet, >*/
/*<      1   xmu,sfactr,mode,modedc,icalc,initf,method,iord,maxord,noncon, >*/
/*<      2   iterno,itemno,nosolv,modac,ipiv,ivmflg,ipostp,iscrch,iofile >*/
/* spice version 2g.6  sccsid=knstnt 3/15/83 */
/*<       common /knstnt/ twopi,xlog2,xlog10,root2,rad,boltz,charge,ctok, >*/
/*<      1   gmin,reltol,abstol,vntol,trtol,chgtol,eps0,epssil,epsox, >*/
/*<      2   pivtol,pivrel >*/
/* spice version 2g.6  sccsid=blank 3/15/83 */
/*<       common /blank/ value(200000) >*/
/*<       integer nodplc(64) >*/
/*<       complex cvalue(32) >*/
/*<       equivalence (value(1),nodplc(1),cvalue(1)) >*/


/*<       dimension qcap(1),ccap(1),diff(8),deltmp(7),coef(6) >*/
/*<       equivalence (qcap(1),value(1)),(ccap(1),value(2)) >*/
/*<       data coef / 5.000000000d-1, 2.222222222d-1, 1.363636364d-1, >*/
/*<      1            9.600000000d-2, 7.299270073d-2, 5.830903790d-2 / >*/
/*<       data xtwelv / 8.333333333d-2 / >*/


/*<       tol=reltol*dmax1(dabs(ccap(lx0+loct)),dabs(ccap(lx1+loct)))+abstol >*/
/* Computing MAX */
    d_3 = (d_1 = ccap[tabinf_1.lx0 + *loct - 1], abs(d_1)), d_4 = (d_2 = ccap[
	    tabinf_1.lx1 + *loct - 1], abs(d_2));
    tol = knstnt_1.reltol * max(d_4,d_3) + knstnt_1.abstol;
/*<       ctol=reltol*dmax1(dabs(qcap(lx0+loct)),dabs(qcap(lx1+loct)), >*/
/*<      1   chgtol)/delta >*/
/* Computing MAX */
    d_3 = (d_1 = qcap[tabinf_1.lx0 + *loct - 1], abs(d_1)), d_4 = (d_2 = qcap[
	    tabinf_1.lx1 + *loct - 1], abs(d_2)), d_3 = max(d_4,d_3);
    ctol = knstnt_1.reltol * max(knstnt_1.chgtol,d_3) / status_1.delta;
/*<       tol=dmax1(tol,ctol) >*/
    tol = max(tol,ctol);

/*  determine divided differences */

/*<       go to (6,5,4,3,2,1), iord >*/
    switch (status_1.iord) {
	case 1:  goto L6;
	case 2:  goto L5;
	case 3:  goto L4;
	case 4:  goto L3;
	case 5:  goto L2;
	case 6:  goto L1;
    }
/*<     1 diff(8)=qcap(lx7+loct) >*/
L1:
    diff[7] = qcap[tabinf_1.lx7 + *loct - 1];
/*<     2 diff(7)=qcap(lx6+loct) >*/
L2:
    diff[6] = qcap[tabinf_1.lx6 + *loct - 1];
/*<     3 diff(6)=qcap(lx5+loct) >*/
L3:
    diff[5] = qcap[tabinf_1.lx5 + *loct - 1];
/*<     4 diff(5)=qcap(lx4+loct) >*/
L4:
    diff[4] = qcap[tabinf_1.lx4 + *loct - 1];
/*<     5 diff(4)=qcap(lx3+loct) >*/
L5:
    diff[3] = qcap[tabinf_1.lx3 + *loct - 1];
/*<     6 diff(3)=qcap(lx2+loct) >*/
L6:
    diff[2] = qcap[tabinf_1.lx2 + *loct - 1];
/*<       diff(2)=qcap(lx1+loct) >*/
    diff[1] = qcap[tabinf_1.lx1 + *loct - 1];
/*<       diff(1)=qcap(lx0+loct) >*/
    diff[0] = qcap[tabinf_1.lx0 + *loct - 1];
/*<       istop=iord+1 >*/
    istop = status_1.iord + 1;
/*<       do 10 i=1,istop >*/
    i_1 = istop;
    for (i = 1; i <= i_1; ++i) {
/*<       deltmp(i)=delold(i) >*/
	deltmp[i - 1] = status_1.delold[i - 1];
/*<    10 continue >*/
/* L10: */
    }
/*<    20 do 30 i=1,istop >*/
L20:
    i_1 = istop;
    for (i = 1; i <= i_1; ++i) {
/*<       diff(i)=(diff(i)-diff(i+1))/deltmp(i) >*/
	diff[i - 1] = (diff[i - 1] - diff[i]) / deltmp[i - 1];
/*<    30 continue >*/
/* L30: */
    }
/*<       istop=istop-1 >*/
    --istop;
/*<       if (istop.eq.0) go to 100 >*/
    if (istop == 0) {
	goto L100;
    }
/*<       do 40 i=1,istop >*/
    i_1 = istop;
    for (i = 1; i <= i_1; ++i) {
/*<       deltmp(i)=deltmp(i+1)+delold(i) >*/
	deltmp[i - 1] = deltmp[i] + status_1.delold[i - 1];
/*<    40 continue >*/
/* L40: */
    }
/*<       go to 20 >*/
    goto L20;

/*  diff(1) contains divided difference */

/*<   100 const=coef(iord) >*/
L100:
    const_ = coef[status_1.iord - 1];
/*<       if ((method.eq.1).and.(iord.eq.2)) const=xtwelv >*/
    if (status_1.method == 1 && status_1.iord == 2) {
	const_ = xtwelv;
    }
/*<       del=trtol*tol/dmax1(abstol,const*dabs(diff(1))) >*/
/* Computing MAX */
    d_1 = knstnt_1.abstol, d_2 = const_ * abs(diff[0]);
    del = knstnt_1.trtol * tol / max(d_2,d_1);
/*<       if (iord.eq.1) go to 200 >*/
    if (status_1.iord == 1) {
	goto L200;
    }
/*<       if (iord.ge.3) go to 150 >*/
    if (status_1.iord >= 3) {
	goto L150;
    }
/*<       del=dsqrt(del) >*/
    del = sqrt(del);
/*<       go to 200 >*/
    goto L200;
/*<   150 del=dexp(dlog(del)/dble(iord)) >*/
L150:
    del = exp(log(del) / (doublereal) status_1.iord);
/*<   200 delnew=dmin1(delnew,del) >*/
L200:
    *delnew = min(*delnew,del);
/*<       return >*/
    return 0;
/*<       end >*/
} /* terr_ */

#undef cvalue
#undef nodplc
#undef qcap
#undef ccap


