b3v1ald.c

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

                             + (ckt->CKTgmin/m) * 0.5;


              here->BSIM3v1Acbs += here->BSIM3v1Agbs * (vbs - 0.5); 

          }

          DrainSatCurrent = here->BSIM3v1AdrainArea * model->BSIM3v1AjctSatCurDensity;
	  if (here->BSIM3v1AdrainArea <= 0.0)
              DrainSatCurrent = 1.0e-14;
          if (DrainSatCurrent <= 0.0)
	  {   here->BSIM3v1Agbd = ckt->CKTgmin/m;

              here->BSIM3v1Acbd = here->BSIM3v1Agbd * vbd;
          }

	  else if (vbd <= -3*CONSTvt0)
          {   arg=3*CONSTvt0/(vbd*CONSTe);
              arg = arg * arg * arg;

              here->BSIM3v1Acbd = -DrainSatCurrent*(1+arg)+(ckt->CKTgmin/m)*vbd;
              here->BSIM3v1Agbd = DrainSatCurrent*3*arg/vbd+(ckt->CKTgmin/m);

          }

	  else if (vbd < 0.5)
	  {   evbd = exp(vbd/CONSTvt0);
              here->BSIM3v1Agbd = DrainSatCurrent * evbd / CONSTvt0 
                             + (ckt->CKTgmin/m);

              here->BSIM3v1Acbd = DrainSatCurrent * (evbd - 1.0) 
                             + (ckt->CKTgmin/m) * vbd;

          }
	  else
	  {   evbd = exp(0.5/CONSTvt0);
              here->BSIM3v1Agbd = DrainSatCurrent * evbd / CONSTvt0 
                             + (ckt->CKTgmin/m);

              here->BSIM3v1Acbd = DrainSatCurrent * (evbd - 1.0) 
                             + (ckt->CKTgmin/m) * 0.5;


              here->BSIM3v1Acbd += here->BSIM3v1Agbd * (vbd - 0.5); 

          }

          if (vds >= 0.0)
	  {   /* normal mode */
              here->BSIM3v1Amode = 1;
              Vds = vds;
              Vgs = vgs;
              Vbs = vbs;
          }
	  else
	  {   /* inverse mode */
              here->BSIM3v1Amode = -1;
              Vds = -vds;
              Vgs = vgd;
              Vbs = vbd;
          }

          ChargeComputationNeeded =  
                 ((ckt->CKTmode & (MODEAC | MODETRAN | MODEINITSMSIG)) ||
                 ((ckt->CKTmode & MODETRANOP) && (ckt->CKTmode & MODEUIC)))
                 ? 1 : 0;

	  T0 = Vbs - pParam->BSIM3v1Avbsc - 0.001;
	  T1 = sqrt(T0 * T0 - 0.004 * pParam->BSIM3v1Avbsc);
	  Vbseff = pParam->BSIM3v1Avbsc + 0.5 * (T0 + T1);
	  dVbseff_dVb = 0.5 * (1.0 + T0 / T1);

          if (Vbseff > 0.0)
	  {   T0 = pParam->BSIM3v1Aphi / (pParam->BSIM3v1Aphi + Vbseff);
              Phis = pParam->BSIM3v1Aphi * T0;
              dPhis_dVb = -T0 * T0;
              sqrtPhis = pParam->BSIM3v1Aphis3 / (pParam->BSIM3v1Aphi + 0.5 * Vbseff);
              dsqrtPhis_dVb = -0.5 * sqrtPhis * sqrtPhis / pParam->BSIM3v1Aphis3;
          }
	  else
	  {   Phis = pParam->BSIM3v1Aphi - Vbseff;
              dPhis_dVb = -1.0;
              sqrtPhis = sqrt(Phis);
              dsqrtPhis_dVb = -0.5 / sqrtPhis; 
          }
          Xdep = pParam->BSIM3v1AXdep0 * sqrtPhis / pParam->BSIM3v1AsqrtPhi;
          dXdep_dVb = (pParam->BSIM3v1AXdep0 / pParam->BSIM3v1AsqrtPhi)
		    * dsqrtPhis_dVb;

          Leff = pParam->BSIM3v1Aleff;
/* Vth Calculation */
          if ((T1 = 1.0 + pParam->BSIM3v1Advt2 * Vbseff) < 1.0e-10)
	      T1 = 1.0E-10;
          if ((T2 = 1.0 + pParam->BSIM3v1Advt2w * Vbseff) < 1.0E-10)
	      T2 = 1.0E-10;

          T3 = sqrt(Xdep);
          lt1 = model->BSIM3v1Afactor1 * T3 * T1;
          dlt1_dVb = model->BSIM3v1Afactor1 * (0.5 / T3 * T1 * dXdep_dVb
                   + T3 * pParam->BSIM3v1Advt2);

          ltw = model->BSIM3v1Afactor1 * T3 * T2;
          dltw_dVb = model->BSIM3v1Afactor1 * (0.5 / T3 * T2 * dXdep_dVb
                   + T3 * pParam->BSIM3v1Advt2w);

          T0 = -0.5 * pParam->BSIM3v1Advt1 * Leff / lt1;
          if (T0 > -EXP_THRESHOLD)
          {   T1 = exp(T0);
              dT1_dVb = -T0 / lt1 * T1 * dlt1_dVb;
          }
          else
          {   T1 = MIN_EXP;
              dT1_dVb = 0.0;
          }

          Theta0 = T1 * (1.0 + 2.0 * T1);
          dTheta0_dVb = (1.0 + 4.0 * T1) * dT1_dVb;

          here->BSIM3v1Athetavth = pParam->BSIM3v1Advt0 * Theta0;
          T0 = pParam->BSIM3v1Avbi - pParam->BSIM3v1Aphi;
          Delt_vth = here->BSIM3v1Athetavth * T0;
          dDelt_vth_dVb = pParam->BSIM3v1Advt0 * dTheta0_dVb * T0;

          T0 = -0.5 * pParam->BSIM3v1Advt1w * pParam->BSIM3v1Aweff * Leff / ltw;
          if (T0 > -EXP_THRESHOLD)
          {   T1 = exp(T0);
              dT1_dVb = -T0 / ltw * T1 * dltw_dVb;
          }
          else
          {   T1 = MIN_EXP;
              dT1_dVb = 0.0;
          }

          T2 = T1 * (1.0 + 2.0 * T1);
          dT2_dVb = (1.0 + 4.0 * T1) * dT1_dVb;

          T0 = pParam->BSIM3v1Advt0w * T2;
          T1 = pParam->BSIM3v1Avbi - pParam->BSIM3v1Aphi;
          T2 = T0 * T1;
          dT2_dVb = pParam->BSIM3v1Advt0w * dT2_dVb * T1;

          TempRatio =  ckt->CKTtemp / model->BSIM3v1Atnom - 1.0;
          T0 = sqrt(1.0 + pParam->BSIM3v1Anlx / Leff);
          T1 = pParam->BSIM3v1Ak1 * (T0 - 1.0) * pParam->BSIM3v1AsqrtPhi
             + (pParam->BSIM3v1Akt1 + pParam->BSIM3v1Akt1l / Leff
             + pParam->BSIM3v1Akt2 * Vbseff) * TempRatio;
          tmp2 = model->BSIM3v1Atox / (pParam->BSIM3v1Aweff
	       + pParam->BSIM3v1Aw0) * pParam->BSIM3v1Aphi;

          dDIBL_Sft_dVd = (pParam->BSIM3v1Aeta0 + pParam->BSIM3v1Aetab
                        * Vbseff) * pParam->BSIM3v1Atheta0vb0;
          DIBL_Sft = dDIBL_Sft_dVd * Vds;

          Vth = model->BSIM3v1Atype * pParam->BSIM3v1Avth0 + pParam->BSIM3v1Ak1 
              * (sqrtPhis - pParam->BSIM3v1AsqrtPhi) - pParam->BSIM3v1Ak2 
              * Vbseff - Delt_vth - T2 + (pParam->BSIM3v1Ak3 + pParam->BSIM3v1Ak3b
              * Vbseff) * tmp2 + T1 - DIBL_Sft;

          here->BSIM3v1Avon = Vth; 

          dVth_dVb = pParam->BSIM3v1Ak1 * dsqrtPhis_dVb - pParam->BSIM3v1Ak2
                   - dDelt_vth_dVb - dT2_dVb + pParam->BSIM3v1Ak3b * tmp2
                   - pParam->BSIM3v1Aetab * Vds * pParam->BSIM3v1Atheta0vb0
                   + pParam->BSIM3v1Akt2 * TempRatio;
          dVth_dVd = -dDIBL_Sft_dVd; 

/* Calculate n */
          tmp2 = pParam->BSIM3v1Anfactor * EPSSI / Xdep;
          tmp3 = pParam->BSIM3v1Acdsc + pParam->BSIM3v1Acdscb * Vbseff
               + pParam->BSIM3v1Acdscd * Vds;
	  n = 1.0 + (tmp2 + tmp3 * Theta0 + pParam->BSIM3v1Acit) / model->BSIM3v1Acox;
	  if (n > 1.0)
	  {   dn_dVb = (-tmp2 / Xdep * dXdep_dVb + tmp3 * dTheta0_dVb
                     + pParam->BSIM3v1Acdscb * Theta0) / model->BSIM3v1Acox;  
              dn_dVd = pParam->BSIM3v1Acdscd * Theta0 / model->BSIM3v1Acox;            
	  }
	  else
	  {   n = 1.0;
	      dn_dVb = dn_dVd = 0.0;  
	  }

/* Poly Gate Si Depletion Effect */
	  T0 = pParam->BSIM3v1Avfb + pParam->BSIM3v1Aphi;
          if (model->BSIM3v1AngateGiven && (Vgs > T0))
	  {   T1 = 1.0e6 * Charge_q * EPSSI * pParam->BSIM3v1Angate
                 / (model->BSIM3v1Acox * model->BSIM3v1Acox);
              T4 = sqrt(1.0 + 2.0 * (Vgs - T0) / T1);
              T2 = T1 * (T4 - 1.0);
              T3 = 0.5 * T2 * T2 / T1;

              if (T3 < 1.12)
	      {   Vgs_eff = T0 + T2;
                  dVgs_eff_dVg = 1.0 / T4;
              }
	      else
	      {   Vgs_eff = Vgs - 1.12;
                  dVgs_eff_dVg = 1.0;
              }
          }
	  else
	  {   Vgs_eff = Vgs;
              dVgs_eff_dVg = 1.0;
          }
          Vgst = Vgs_eff - Vth;

/* Effective Vgst (Vgsteff) Calculation */
          Vtm = model->BSIM3v1Avtm;

          T10 = 2.0 * n * Vtm;
          VgstNVt = Vgst / T10;
          if (VgstNVt < -EXP_THRESHOLD) 
	  {   T1 = T10 * MIN_EXP; 
	      dT1_dVg = dT1_dVd = dT1_dVb = 0.0;
	  }
          else if (VgstNVt > EXP_THRESHOLD)
	  {   T1 = Vgst;
              dT1_dVg = dVgs_eff_dVg;
              dT1_dVd = -dVth_dVd;
              dT1_dVb = -dVth_dVb;
          }
	  else
	  {   ExpVgst = exp(VgstNVt);
              T1 = T10 * log(1.0 + ExpVgst);
              dT1_dVg = ExpVgst / (1.0 + ExpVgst);
              dT1_dVb = -dT1_dVg * (dVth_dVb + Vgst / n * dn_dVb)
		      + T1 / n * dn_dVb; 
	      dT1_dVd = -dT1_dVg * (dVth_dVd + Vgst / n * dn_dVd)
		      + T1 / n * dn_dVd;
	      dT1_dVg *= dVgs_eff_dVg;
          }

	  T2 = model->BSIM3v1Atox / (pParam->BSIM3v1Aweff + pParam->BSIM3v1Aw0);
          ExpArg = (2.0 * pParam->BSIM3v1Avoff - Vgst) / T10;
          if (ExpArg < -EXP_THRESHOLD)
	  {   T2 = 1.0;
              dT2_dVg = dT2_dVd = dT2_dVb = 0.0;

          }
	  else if (ExpArg > EXP_THRESHOLD)
	  {   T2 = 1.0 + 2.0 * n * model->BSIM3v1Acox / pParam->BSIM3v1Acdep0
		 * MAX_EXP;
              dT2_dVg = dT2_dVd = dT2_dVb = 0.0;
          }
	  else
	  {   dT2_dVg = -model->BSIM3v1Acox / Vtm / pParam->BSIM3v1Acdep0
		      * exp(ExpArg);
              T2 = 1.0 - T10 * dT2_dVg;
              dT2_dVd = -dT2_dVg * (dVth_dVd - 2.0 * Vtm
		      * ExpArg * dn_dVd) + (T2 - 1.0) / n * dn_dVd;
              dT2_dVb = -dT2_dVg * (dVth_dVb - 2.0 * Vtm
		      * ExpArg * dn_dVb) + (T2 - 1.0) / n * dn_dVb;
	      dT2_dVg *= dVgs_eff_dVg;
          }

          Vgsteff = T1 / T2;
          dVgsteff_dVg = (T2 * dT1_dVg - T1 * dT2_dVg) / (T2 * T2);
          dVgsteff_dVd = (T2 * dT1_dVd - T1 * dT2_dVd) / (T2 * T2);
          dVgsteff_dVb = (T2 * dT1_dVb - T1 * dT2_dVb) / (T2 * T2);

/* Calculate Effective Channel Geometry */
          Weff = pParam->BSIM3v1Aweff - 2.0 * (pParam->BSIM3v1Adwg * Vgsteff 
               + pParam->BSIM3v1Adwb * (sqrtPhis - pParam->BSIM3v1AsqrtPhi)); 
          dWeff_dVg = -2.0 * pParam->BSIM3v1Adwg;
          dWeff_dVb = -2.0 * pParam->BSIM3v1Adwb * dsqrtPhis_dVb;

          if (Weff < 1.0e-8)
	  {   Weff = 1.0e-8;
              dWeff_dVg = dWeff_dVb = 0;
          }

          Rds = pParam->BSIM3v1Ards0 * (1.0 + pParam->BSIM3v1Aprwg * Vgsteff 
              + pParam->BSIM3v1Aprwb * (sqrtPhis - pParam->BSIM3v1AsqrtPhi));
	  if (Rds > 0.0)
	  {   dRds_dVg = pParam->BSIM3v1Ards0 * pParam->BSIM3v1Aprwg;
              dRds_dVb = pParam->BSIM3v1Ards0 * pParam->BSIM3v1Aprwb * dsqrtPhis_dVb;
	  }
	  else
	  {   Rds = dRds_dVg = dRds_dVb = 0.0;
	  }
	  
          WVCox = Weff * pParam->BSIM3v1Avsattemp * model->BSIM3v1Acox;
          WVCoxRds = WVCox * Rds; 

/* Calculate Abulk */
          T0 = 1.0 / (1.0 + pParam->BSIM3v1Aketa * Vbseff);
          dT0_dVb = -pParam->BSIM3v1Aketa * T0 * T0;

          T1 = 0.5 * pParam->BSIM3v1Ak1 / sqrtPhis;
          dT1_dVb = -T1 / sqrtPhis * dsqrtPhis_dVb;

          tmp1 = Leff + 2.0 * sqrt(pParam->BSIM3v1Axj * Xdep);
          T5 = Leff / tmp1; 
          tmp2 = pParam->BSIM3v1Aa0 *T5;
          tmp3 = pParam->BSIM3v1Aweff + pParam->BSIM3v1Ab1; 
          tmp4 = pParam->BSIM3v1Ab0 / tmp3;
          T2 = tmp2 + tmp4;
          dT2_dVb = -tmp2 / tmp1 * sqrt(pParam->BSIM3v1Axj / Xdep) * dXdep_dVb;
          T6 = T5 * T5;
          T7 = T5 * T6;
          Abulk0 = T0 * (1.0 + T1 * T2); 
                   if (Abulk0 < 0.01)
                      Abulk0= 0.01; 
          T8 = pParam->BSIM3v1Aags * pParam->BSIM3v1Aa0 * T7;
          dAbulk_dVg = -T1 * T0 * T8;
          Abulk = Abulk0 + dAbulk_dVg * Vgsteff; 
                   if (Abulk < 0.01)
                        Abulk= 0.01;
                      
          dAbulk0_dVb = T0 * T1 * dT2_dVb + T0 * T2 * dT1_dVb   
               	      + (1.0 + T1 * T2) * dT0_dVb;
          dAbulk_dVb = dAbulk0_dVb - T8 * Vgsteff * (T1 * (3.0 * T0 * dT2_dVb
		     / tmp2 + dT0_dVb) + T0 * dT1_dVb); 
/* Mobility calculation */

          if (model->BSIM3v1AmobMod == 1)
	  {   T0 = Vgsteff + Vth + Vth;
              T2 = pParam->BSIM3v1Aua + pParam->BSIM3v1Auc * Vbseff;
              T3 = T0 / model->BSIM3v1Atox;
              Denomi = 1.0 + T3 * (T2
                     + pParam->BSIM3v1Aub * T3);
              T1 = T2 / model->BSIM3v1Atox + 2.0 * pParam->BSIM3v1Aub * T3
                 / model->BSIM3v1Atox;
              dDenomi_dVg = T1;
              dDenomi_dVd = T1 * 2.0 * dVth_dVd;
              dDenomi_dVb = T1 * 2.0 * dVth_dVb + pParam->BSIM3v1Auc * T3;
          }
	  else if (model->BSIM3v1AmobMod == 2)
	  {   Denomi = 1.0 + Vgsteff / model->BSIM3v1Atox * (pParam->BSIM3v1Aua
		     + pParam->BSIM3v1Auc * Vbseff + pParam->BSIM3v1Aub * Vgsteff
                     / model->BSIM3v1Atox);
              T1 = (pParam->BSIM3v1Aua + pParam->BSIM3v1Auc * Vbseff) / model->BSIM3v1Atox