b3v1ald.c

ngspice又一个电子CAD仿真软件代码.功能更全

/**********
Copyright 1990 Regents of the University of California.  All rights reserved.
Author: 1991 JianHui Huang and Min-Chie Jeng.
File: b3v1ald.c          1/3/92
Modified by Mansun Chan  (1995)
Modified by Paolo Nenzi 2002
**********/

#include "ngspice.h"
#include "cktdefs.h"
#include "bsim3v1adef.h"
#include "trandefs.h"
#include "const.h"
#include "sperror.h"
#include "devdefs.h"
#include "suffix.h"

#define MAX_EXP 5.834617425e14
#define MIN_EXP 1.713908431e-15
#define EXP_THRESHOLD 34.0
#define EPSOX 3.453133e-11
#define EPSSI 1.03594e-10
#define Charge_q 1.60219e-19
#define DELTA_1 0.02
#define DELTA_2 0.02
#define DELTA_3 0.02
#define DELTA_4 0.02


int
BSIM3v1Aload(GENmodel *inModel, CKTcircuit *ckt)
{
BSIM3v1Amodel *model = (BSIM3v1Amodel*)inModel;
BSIM3v1Ainstance *here;
double SourceSatCurrent, DrainSatCurrent;
double ag0, qgd, qgs, qgb, von, cbhat, VgstNVt, ExpVgst = 0.0;
double cdrain, cdhat, cdreq, ceqbd, ceqbs, ceqqb, ceqqd, ceqqg, ceq, geq;
double czbd, czbdsw, czbs, czbssw, evbd, evbs, arg, sarg;
double delvbd, delvbs, delvds, delvgd, delvgs;
double Vfbeff, dVfbeff_dVg, dVfbeff_dVd, dVfbeff_dVb, V3, V4;
double gcbdb, gcbgb, gcbsb, gcddb, gcdgb, gcdsb, gcgdb, gcggb, gcgsb, gcsdb;
double gcsgb, gcssb, PhiB, PhiBSW, MJ, MJSW;
double vbd, vbs, vds, vgb, vgd, vgs, vgdo, xfact;
double qgate = 0.0, qbulk = 0.0, qdrn = 0.0, qsrc, cqgate, cqbulk, cqdrn;
double Vds, Vgs, Vbs, Gmbs, FwdSum, RevSum;
double Vgs_eff, Vfb, dVfb_dVb, dVfb_dVd;
double Phis, dPhis_dVb, sqrtPhis, dsqrtPhis_dVb, Vth, dVth_dVb, dVth_dVd;
double Vgst, dVgs_eff_dVg;
double n, dn_dVb, Vtm;
double ExpArg;
double Denomi, dDenomi_dVg, dDenomi_dVd, dDenomi_dVb;
double ueff, dueff_dVg, dueff_dVd, dueff_dVb; 
double Esat, Vdsat;
double EsatL, dEsatL_dVg, dEsatL_dVd, dEsatL_dVb;
double dVdsat_dVg, dVdsat_dVb, dVdsat_dVd, Vasat; 
double dVasat_dVg, dVasat_dVb, dVasat_dVd, Va, dVa_dVd, dVa_dVg, dVa_dVb; 
double Vbseff, dVbseff_dVb, VbseffCV, dVbseffCV_dVb; 
double Arg1, One_Third_CoxWL, Two_Third_CoxWL, CoxWL; 
double T0, dT0_dVg, dT0_dVd, dT0_dVb;
double T1, dT1_dVg, dT1_dVd, dT1_dVb;
double T2, dT2_dVg, dT2_dVd, dT2_dVb;
double T3, dT3_dVg, dT3_dVd, dT3_dVb;
double T4, dT4_dVg, dT4_dVd, dT4_dVb;
double T5;
double T6;
double T7;
double T8;
double T10;
double Abulk, dAbulk_dVb, Abulk0, dAbulk0_dVb;
double VACLM, dVACLM_dVg, dVACLM_dVd, dVACLM_dVb;
double VADIBL, dVADIBL_dVg, dVADIBL_dVd, dVADIBL_dVb;
double Xdep, dXdep_dVb, lt1, dlt1_dVb, ltw, dltw_dVb, Delt_vth, dDelt_vth_dVb;
double Theta0, dTheta0_dVb;
double TempRatio, tmp1, tmp2, tmp3, tmp4;
double DIBL_Sft, dDIBL_Sft_dVd, Pmos_factor;
double tempv, a1;

double Vgsteff, dVgsteff_dVg, dVgsteff_dVd, dVgsteff_dVb; 
double Vdseff, dVdseff_dVg, dVdseff_dVd, dVdseff_dVb; 
double VdseffCV, dVdseffCV_dVg, dVdseffCV_dVd, dVdseffCV_dVb; 
double diffVds;
double dAbulk_dVg, dn_dVd ;
double beta, dbeta_dVg, dbeta_dVd, dbeta_dVb;
double gche, dgche_dVg, dgche_dVd, dgche_dVb;
double fgche1, dfgche1_dVg, dfgche1_dVd, dfgche1_dVb;
double fgche2, dfgche2_dVg, dfgche2_dVd, dfgche2_dVb;
double Idl, dIdl_dVg, dIdl_dVd, dIdl_dVb;
double Idsa, dIdsa_dVg, dIdsa_dVd, dIdsa_dVb;
double Ids, Gm, Gds, Gmb;
double Isub, Gbd, Gbg, Gbb;
double VASCBE, dVASCBE_dVg, dVASCBE_dVd, dVASCBE_dVb;
double CoxWovL;
double Rds, dRds_dVg, dRds_dVb, WVCox, WVCoxRds;
double Vgst2Vtm, VdsatCV, dVdsatCV_dVg, dVdsatCV_dVb;
double Leff = 0.0, Weff, dWeff_dVg, dWeff_dVb;
double AbulkCV, dAbulkCV_dVb;
double qgdo, qgso, cgdo, cgso;

double qcheq, qdef, gqdef, cqdef, cqcheq, gtau_diff, gtau_drift;
double gcqdb,gcqsb,gcqgb,gcqbb;
double dxpart, sxpart;

double gbspsp, gbbdp, gbbsp, gbspg, gbspb, gbspdp; 
double gbdpdp, gbdpg, gbdpb, gbdpsp; 
double Cgg, Cgd, Cgb;
double Csg, Csd, Csb, Cbg, Cbd, Cbb;
double Cgg1, Cgb1, Cgd1, Cbg1, Cbb1, Cbd1, Qac0, Qsub0;
double dQac0_dVg, dQac0_dVd, dQac0_dVb, dQsub0_dVg, dQsub0_dVd, dQsub0_dVb;
   
double m = 0.0;

struct bsim3v1aSizeDependParam *pParam;
int ByPass, Check, ChargeComputationNeeded = 0, error;

for (; model != NULL; model = model->BSIM3v1AnextModel)
{    for (here = model->BSIM3v1Ainstances; here != NULL; 
          here = here->BSIM3v1AnextInstance)
     {    
     
          if (here->BSIM3v1Aowner != ARCHme)
	          continue;

     
          Check = 1;
          ByPass = 0;
	  pParam = here->pParam;
          if ((ckt->CKTmode & MODEINITSMSIG))
	  {   vbs = *(ckt->CKTstate0 + here->BSIM3v1Avbs);
              vgs = *(ckt->CKTstate0 + here->BSIM3v1Avgs);
              vds = *(ckt->CKTstate0 + here->BSIM3v1Avds);
              qdef = *(ckt->CKTstate0 + here->BSIM3v1Aqcdump);
          }
	  else if ((ckt->CKTmode & MODEINITTRAN))
	  {   vbs = *(ckt->CKTstate1 + here->BSIM3v1Avbs);
              vgs = *(ckt->CKTstate1 + here->BSIM3v1Avgs);
              vds = *(ckt->CKTstate1 + here->BSIM3v1Avds);
              qdef = *(ckt->CKTstate1 + here->BSIM3v1Aqcdump);
          }
	  else if ((ckt->CKTmode & MODEINITJCT) && !here->BSIM3v1Aoff)
	  {   vds = model->BSIM3v1Atype * here->BSIM3v1AicVDS;
              vgs = model->BSIM3v1Atype * here->BSIM3v1AicVGS;
              vbs = model->BSIM3v1Atype * here->BSIM3v1AicVBS;
              qdef = 0.0;

              if ((vds == 0.0) && (vgs == 0.0) && (vbs == 0.0) && 
                  ((ckt->CKTmode & (MODETRAN | MODEAC|MODEDCOP |
                   MODEDCTRANCURVE)) || (!(ckt->CKTmode & MODEUIC))))
	      {   vbs = 0.0;
                  vgs = pParam->BSIM3v1Avth0 + 0.1;
                  vds = 0.1;
              }
          }
	  else if ((ckt->CKTmode & (MODEINITJCT | MODEINITFIX)) && 
                  (here->BSIM3v1Aoff)) 
          {    qdef = vbs = vgs = vds = 0.0;
	  }
          else
	  {
#ifndef PREDICTOR
               if ((ckt->CKTmode & MODEINITPRED))
	       {   xfact = ckt->CKTdelta / ckt->CKTdeltaOld[1];
                   *(ckt->CKTstate0 + here->BSIM3v1Avbs) = 
                         *(ckt->CKTstate1 + here->BSIM3v1Avbs);
                   vbs = (1.0 + xfact)* (*(ckt->CKTstate1 + here->BSIM3v1Avbs))
                         - (xfact * (*(ckt->CKTstate2 + here->BSIM3v1Avbs)));
                   *(ckt->CKTstate0 + here->BSIM3v1Avgs) = 
                         *(ckt->CKTstate1 + here->BSIM3v1Avgs);
                   vgs = (1.0 + xfact)* (*(ckt->CKTstate1 + here->BSIM3v1Avgs))
                         - (xfact * (*(ckt->CKTstate2 + here->BSIM3v1Avgs)));
                   *(ckt->CKTstate0 + here->BSIM3v1Avds) = 
                         *(ckt->CKTstate1 + here->BSIM3v1Avds);
                   vds = (1.0 + xfact)* (*(ckt->CKTstate1 + here->BSIM3v1Avds))
                         - (xfact * (*(ckt->CKTstate2 + here->BSIM3v1Avds)));
                   *(ckt->CKTstate0 + here->BSIM3v1Avbd) = 
                         *(ckt->CKTstate0 + here->BSIM3v1Avbs)
                         - *(ckt->CKTstate0 + here->BSIM3v1Avds);
                   qdef = (1.0 + xfact)* (*(ckt->CKTstate1 + here->BSIM3v1Aqcdump))
                           -(xfact * (*(ckt->CKTstate2 + here->BSIM3v1Aqcdump)));
               }
	       else
	       {
#endif /* PREDICTOR */
                   vbs = model->BSIM3v1Atype
                       * (*(ckt->CKTrhsOld + here->BSIM3v1AbNode)
                       - *(ckt->CKTrhsOld + here->BSIM3v1AsNodePrime));
                   vgs = model->BSIM3v1Atype
                       * (*(ckt->CKTrhsOld + here->BSIM3v1AgNode) 
                       - *(ckt->CKTrhsOld + here->BSIM3v1AsNodePrime));
                   vds = model->BSIM3v1Atype
                       * (*(ckt->CKTrhsOld + here->BSIM3v1AdNodePrime)
                       - *(ckt->CKTrhsOld + here->BSIM3v1AsNodePrime));
                   qdef = *(ckt->CKTrhsOld + here->BSIM3v1AqNode);
#ifndef PREDICTOR
               }
#endif /* PREDICTOR */

               vbd = vbs - vds;
               vgd = vgs - vds;
               vgdo = *(ckt->CKTstate0 + here->BSIM3v1Avgs)
                    - *(ckt->CKTstate0 + here->BSIM3v1Avds);
               delvbs = vbs - *(ckt->CKTstate0 + here->BSIM3v1Avbs);
               delvbd = vbd - *(ckt->CKTstate0 + here->BSIM3v1Avbd);
               delvgs = vgs - *(ckt->CKTstate0 + here->BSIM3v1Avgs);
               delvds = vds - *(ckt->CKTstate0 + here->BSIM3v1Avds);
               delvgd = vgd - vgdo;

               if (here->BSIM3v1Amode >= 0) 
	       {   cdhat = here->BSIM3v1Acd - here->BSIM3v1Agbd * delvbd 
                         + here->BSIM3v1Agmbs * delvbs + here->BSIM3v1Agm * delvgs
                         + here->BSIM3v1Agds * delvds;
	       }
               else
	       {   cdhat = here->BSIM3v1Acd - (here->BSIM3v1Agbd - here->BSIM3v1Agmbs)
                         * delvbd - here->BSIM3v1Agm * delvgd 
                         + here->BSIM3v1Agds * delvds;
               
	       }
               cbhat = here->BSIM3v1Acbs + here->BSIM3v1Acbd + here->BSIM3v1Agbd
                     * delvbd + here->BSIM3v1Agbs * delvbs;

#ifndef NOBYPASS
           /* following should be one big if connected by && all over
            * the place, but some C compilers can't handle that, so
            * we split it up here to let them digest it in stages
            */

               if ((!(ckt->CKTmode & MODEINITPRED)) && (ckt->CKTbypass))
               if ((fabs(delvbs) < (ckt->CKTreltol * MAX(fabs(vbs),
                   fabs(*(ckt->CKTstate0+here->BSIM3v1Avbs))) + ckt->CKTvoltTol)))
               if ((fabs(delvbd) < (ckt->CKTreltol * MAX(fabs(vbd),
                   fabs(*(ckt->CKTstate0+here->BSIM3v1Avbd))) + ckt->CKTvoltTol)))
               if ((fabs(delvgs) < (ckt->CKTreltol * MAX(fabs(vgs),
                   fabs(*(ckt->CKTstate0+here->BSIM3v1Avgs))) + ckt->CKTvoltTol)))
               if ((fabs(delvds) < (ckt->CKTreltol * MAX(fabs(vds),
                   fabs(*(ckt->CKTstate0+here->BSIM3v1Avds))) + ckt->CKTvoltTol)))
               if ((fabs(cdhat - here->BSIM3v1Acd) < ckt->CKTreltol 
                   * MAX(fabs(cdhat),fabs(here->BSIM3v1Acd)) + ckt->CKTabstol))
	       {   tempv = MAX(fabs(cbhat),fabs(here->BSIM3v1Acbs 
                         + here->BSIM3v1Acbd)) + ckt->CKTabstol;
                   if ((fabs(cbhat - (here->BSIM3v1Acbs + here->BSIM3v1Acbd))
                       < ckt->CKTreltol * tempv))
		   {   /* bypass code */
                       vbs = *(ckt->CKTstate0 + here->BSIM3v1Avbs);
                       vbd = *(ckt->CKTstate0 + here->BSIM3v1Avbd);
                       vgs = *(ckt->CKTstate0 + here->BSIM3v1Avgs);
                       vds = *(ckt->CKTstate0 + here->BSIM3v1Avds);
                       qcheq = *(ckt->CKTstate0 + here->BSIM3v1Aqcheq);

                       vgd = vgs - vds;
                       vgb = vgs - vbs;

                       cdrain = here->BSIM3v1Amode * (here->BSIM3v1Acd
                              + here->BSIM3v1Acbd); 
                       if ((ckt->CKTmode & (MODETRAN | MODEAC)) || 
                           ((ckt->CKTmode & MODETRANOP) && 
                           (ckt->CKTmode & MODEUIC)))
		       {   ByPass = 1;
                           goto line755;
                       }
		       else
		       {   goto line850;
		       }
                   }
               }

#endif /*NOBYPASS*/
               von = here->BSIM3v1Avon;
               if (*(ckt->CKTstate0 + here->BSIM3v1Avds) >= 0.0)
	       {   vgs = DEVfetlim(vgs, *(ckt->CKTstate0+here->BSIM3v1Avgs), von);
                   vds = vgs - vgd;
                   vds = DEVlimvds(vds, *(ckt->CKTstate0 + here->BSIM3v1Avds));
                   vgd = vgs - vds;

               }
	       else
	       {   vgd = DEVfetlim(vgd, vgdo, von);
                   vds = vgs - vgd;
                   vds = -DEVlimvds(-vds, -(*(ckt->CKTstate0+here->BSIM3v1Avds)));
                   vgs = vgd + vds;
               }

               if (vds >= 0.0)
	       {   vbs = DEVpnjlim(vbs, *(ckt->CKTstate0 + here->BSIM3v1Avbs),
                                   CONSTvt0, model->BSIM3v1Avcrit, &Check);
                   vbd = vbs - vds;

               }
	       else
	       {   vbd = DEVpnjlim(vbd, *(ckt->CKTstate0 + here->BSIM3v1Avbd),
                                   CONSTvt0, model->BSIM3v1Avcrit, &Check); 
                   vbs = vbd + vds;
               }
          }

          /* determine DC current and derivatives */
          vbd = vbs - vds;
          vgd = vgs - vds;
          vgb = vgs - vbs;

          m=here->BSIM3v1Am;

          SourceSatCurrent = here->BSIM3v1AsourceArea*model->BSIM3v1AjctSatCurDensity;
	  if (here->BSIM3v1AsourceArea <= 0.0)
              SourceSatCurrent = 1.0e-14;
          if (SourceSatCurrent <= 0.0)
	  {   here->BSIM3v1Agbs = ckt->CKTgmin/m;
              here->BSIM3v1Acbs = here->BSIM3v1Agbs * vbs;
          }

	  else if (vbs <= -3*CONSTvt0)
          {   arg=3*CONSTvt0/(vbs*CONSTe);
              arg = arg * arg * arg;
              here->BSIM3v1Acbs = -SourceSatCurrent*(1+arg)+(ckt->CKTgmin/m)*vbs;
              here->BSIM3v1Agbs = SourceSatCurrent*3*arg/vbs+(ckt->CKTgmin/m);
          }

	  else if (vbs < 0.5)
	  {   evbs = exp(vbs / CONSTvt0);
              here->BSIM3v1Agbs = SourceSatCurrent * evbs / CONSTvt0 
                             + (ckt->CKTgmin/m);

              here->BSIM3v1Acbs = SourceSatCurrent * (evbs - 1.0) 
                             + (ckt->CKTgmin/m) * vbs;

          }
	  else
	  {   evbs = exp(0.5 / CONSTvt0);
              here->BSIM3v1Agbs = SourceSatCurrent * evbs / CONSTvt0 
                             + (ckt->CKTgmin/m);

              here->BSIM3v1Acbs = SourceSatCurrent * (evbs - 1.0)