b1dset.c

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

/**********
Copyright 1990 Regents of the University of California.  All rights reserved.
Author: 1985 Hong J. Park, Thomas L. Quarles
**********/

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

int
B1dSetup(GENmodel *inModel, CKTcircuit *ckt)
{
  B1model* model = (B1model*)inModel;
  B1instance *here;
  double DrainSatCurrent;
  double EffectiveLength;
  double GateBulkOverlapCap;
  double GateDrainOverlapCap;
  double GateSourceOverlapCap;
  double SourceSatCurrent;
  double DrainArea;
  double SourceArea;
  double DrainPerimeter;
  double SourcePerimeter;
  double vt0;
  double evbs;
  double lgbs1, lgbs2, lgbs3;
  double czbd, czbs, czbdsw, czbssw;
  double PhiB, MJ, MJSW, PhiBSW;
  double arg, argsw, sarg, sargsw;
  double capbs1, capbs2, capbs3;
  double capbd1, capbd2, capbd3;
  double qg;
  double qb;
  double qd;
  double Vfb;
  double Phi;
  double K1;
  double K2;
  double Vdd;
  double Ugs;
  double Uds;
  double Leff;
  double Eta;
  double Vpb;
  double SqrtVpb;
  double Von;
  double Vth;
  double DrCur;
  double G;
  double A;
  double Arg;
  double Beta;
  double Beta_Vds_0;
  double BVdd;
  double Beta0;
  double VddSquare;
  double C1;
  double C2;
  double VdsSat;
  double Argl1;
  double Argl2;
  double Vc;
  double Term1;
  double K;
  double Args1;
  double Args2;
  double Args3;
  double Warg1;
  double Vcut;
  double N;
  double N0;
  double NB;
  double ND;
  double Warg2;
  double Wds;
  double Wgs;
  double Ilimit;
  double Iexp;
  double Vth0;
  double Arg1;
  double Arg2;
  double Arg3;
  double Arg5;
  double Ent;
  double Vcom;
  double Vgb;
  double Vgb_Vfb;
  double VdsPinchoff;
  double EntSquare;
  double Argl5;
  double Argl6;
  double Argl7;
  double WLCox;
  double Vtsquare;
  int ChargeComputationNeeded;
  double co4v15;
  double VgsVth  = 0.0;
  double lgbd1, lgbd2, lgbd3, evbd; 
  double vbd = 0.0;
  double vgd = 0.0;
  double vgb = 0.0;
  double vds = 0.0;
  double vgs = 0.0;
  double vbs = 0.0;
  double dBVdddVds;
  Dderivs d_Argl6;
  Dderivs d_Vgb, d_Vgb_Vfb, d_EntSquare, d_Argl5, d_Argl7;
  Dderivs d_Arg5, d_Ent, d_Vcom, d_VdsPinchoff;
  Dderivs d_qb, d_qd, d_Vth0, d_Arg1, d_Arg2, d_Arg3;
  Dderivs d_dummy, d_Vth, d_BVdd, d_Warg1, d_qg;
  Dderivs d_N, d_Wds, d_Wgs, d_Iexp;
  Dderivs d_Argl1, d_Argl2, d_Args1, d_Args2, d_Args3;
  Dderivs d_Beta, d_Vc, d_Term1, d_K, d_VdsSat;
  Dderivs d_Beta_Vds_0, d_C1, d_C2, d_Beta0;
  Dderivs d_VgsVth, d_G, d_A, d_Arg, d_DrCur;
  Dderivs d_Ugs, d_Uds, d_Eta, d_Vpb, d_SqrtVpb, d_Von;
  Dderivs d_p, d_q, d_r, d_zero;


  /*  loop through all the B1 device models */
  for( ; model != NULL; model = model->B1nextModel ) {

    /* loop through all the instances of the model */
    for (here = model->B1instances; here != NULL ;
	 here=here->B1nextInstance) {
      
if (here->B1owner != ARCHme) continue;
        
      EffectiveLength=here->B1l - model->B1deltaL * 1.e-6;/* m */
      DrainArea = here->B1m * here->B1drainArea;
      SourceArea = here->B1m * here->B1sourceArea;
      DrainPerimeter = here->B1m * here->B1drainPerimeter;
      SourcePerimeter = here->B1m * here->B1sourcePerimeter;
      if( (DrainSatCurrent=DrainArea*model->B1jctSatCurDensity) 
	  < 1e-15){
	DrainSatCurrent = 1.0e-15;
      }
      if( (SourceSatCurrent=SourceArea*model->B1jctSatCurDensity)
	  <1.0e-15){
	SourceSatCurrent = 1.0e-15;
      }
      GateSourceOverlapCap = model->B1gateSourceOverlapCap * here->B1w * here->B1m;
      GateDrainOverlapCap = model->B1gateDrainOverlapCap * here->B1w * here-> B1m;
      GateBulkOverlapCap = model->B1gateBulkOverlapCap *EffectiveLength;
      vt0 = model->B1type * here->B1vt0;

      vbs = model->B1type * ( 
			     *(ckt->CKTrhsOld+here->B1bNode) -
			     *(ckt->CKTrhsOld+here->B1sNodePrime));
      vgs = model->B1type * ( 
			     *(ckt->CKTrhsOld+here->B1gNode) -
			     *(ckt->CKTrhsOld+here->B1sNodePrime));
      vds = model->B1type * ( 
			     *(ckt->CKTrhsOld+here->B1dNodePrime) -
			     *(ckt->CKTrhsOld+here->B1sNodePrime));

      if(vds >= 0) {
	/* normal mode */
	here->B1mode = 1;
      } else {
	/* inverse mode */
	here->B1mode = -1;
	vds = -vds;
	vgs = vgs + vds; /* these are the local(fake) values now */
	vbs = vbs + vds;
      }

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

      if(vbs <= 0.0 ) {
	lgbs1 = SourceSatCurrent / CONSTvt0 + ckt->CKTgmin;
	lgbs2 = lgbs3 = 0.0;
      } else {
	evbs = exp(vbs/CONSTvt0);
	lgbs1 = SourceSatCurrent*evbs/CONSTvt0 + ckt->CKTgmin;
	lgbs2 = (lgbs1 - ckt->CKTgmin)/(CONSTvt0*2);
	lgbs3 = lgbs2/(CONSTvt0*3);
      }
      if(vbd <= 0.0) {
	lgbd1 = DrainSatCurrent / CONSTvt0 + ckt->CKTgmin;
	lgbd2 = lgbd3 = 0.0;
      } else {
	evbd = exp(vbd/CONSTvt0);
	lgbd1 = DrainSatCurrent*evbd/CONSTvt0 + ckt->CKTgmin;
	lgbd2 = (lgbd1 - ckt->CKTgmin)/(CONSTvt0*2);
	lgbd3 = lgbd2/(CONSTvt0*3);
      }
      /* line 400 */
      /* call B1evaluate to calculate drain current and its 
       * derivatives and charge and capacitances related to gate
       * drain, and bulk
       */
          
      /* check quadratic interpolation for beta0 ; line 360 */
      /*
       * Copyright (c) 1985 Hong J. Park, Thomas L. Quarles
       * modified 1988 Jaijeet Roychowdhury 
       */

      /* This routine evaluates the drain current, its derivatives and the
       * charges associated with the gate,bulk and drain terminal
       * using the B1 (Berkeley Short-Channel IGFET Model) Equations.
       */
      /* as usual p=vgs, q=vbs, r=vds */
      {

        ChargeComputationNeeded  =  1;

	Vfb  =  here->B1vfb;
	Phi  =  here->B1phi;
	K1   =  here->B1K1;
	K2   =  here->B1K2;
	Vdd  =  model->B1vdd;
	d_p.value = 0.0;
	d_p.d1_p = 1.0;
	d_p.d1_q = 0.0;
	d_p.d1_r = 0.0;
	d_p.d2_p2 = 0.0;
	d_p.d2_q2 = 0.0;
	d_p.d2_r2 = 0.0;
	d_p.d2_pq = 0.0;
	d_p.d2_qr = 0.0;
	d_p.d2_pr = 0.0;
	d_p.d3_p3 = 0.0;
	d_p.d3_q3 = 0.0;
	d_p.d3_r3 = 0.0;
	d_p.d3_p2r = 0.0;
	d_p.d3_p2q = 0.0;
	d_p.d3_q2r = 0.0;
	d_p.d3_pq2 = 0.0;
	d_p.d3_pr2 = 0.0;
	d_p.d3_qr2 = 0.0;
	d_p.d3_pqr = 0.0;
	EqualDeriv(&d_q,&d_p);
	EqualDeriv(&d_r,&d_p);
	EqualDeriv(&d_zero,&d_p);
	d_q.d1_p = d_r.d1_p = d_zero.d1_p = 0.0;
	d_q.d1_q = d_r.d1_r = 1.0;
	d_p.value = vgs; d_q.value = vbs; d_r.value = vds;

	if((Ugs  =  here->B1ugs + here->B1ugsB * vbs) <= 0 ) {
	  Ugs = 0;
	  EqualDeriv(&d_Ugs,&d_zero);
	} else {
	  EqualDeriv(&d_Ugs,&d_q);
	  d_Ugs.value = Ugs;
	  d_Ugs.d1_q = here->B1ugsB;
	}
	if((Uds  =  here->B1uds + here->B1udsB * vbs + 
            here->B1udsD*(vds-Vdd)) <= 0 ) {    
	  Uds = 0.0;
	  EqualDeriv(&d_Uds,&d_zero);
	} else {
	  Leff  =  here->B1l * 1.e6 - model->B1deltaL; /* Leff in um */
	  /*const*/
	  Uds  =  Uds / Leff;
	  /* Uds = (here->B1uds + here->B1udsB * vbs here->B1udsD*
	     (vds-Vdd))/Leff */
	  EqualDeriv(&d_Uds,&d_r);
	  d_Uds.value = Uds;
	  d_Uds.d1_r = here->B1udsD/Leff;
	  d_Uds.d1_q = here->B1udsB/Leff;
	}
	Eta  =  here->B1eta + here->B1etaB * vbs + here->B1etaD * 
	  (vds - Vdd);
	EqualDeriv(&d_Eta,&d_r);
	d_Eta.value = Eta;
	d_Eta.d1_r = here->B1etaD;
	d_Eta.d1_q = here->B1etaB;

	if( Eta <= 0 ) {   
	  Eta  = 0; 
	  EqualDeriv(&d_Eta,&d_zero);
	} else if ( Eta > 1 ) {
	  Eta = 1;
	  EqualDeriv(&d_Eta,&d_zero);
	  d_Eta.value = 1.0;
	} 
	if( vbs < 0 ) {
	  Vpb  =  Phi - vbs;
	  EqualDeriv(&d_Vpb,&d_q);
	  d_Vpb.value = Vpb;
	  d_Vpb.d1_q = -1;
	} else {
	  Vpb  =  Phi;
	  EqualDeriv(&d_Vpb,&d_zero);
	  d_Vpb.value = Phi;
	}
	SqrtVpb  =  sqrt( Vpb );
	SqrtDeriv(&d_SqrtVpb,&d_Vpb);
	Von  = Vfb + Phi + K1 * SqrtVpb - K2 * Vpb - Eta * vds;
	EqualDeriv(&d_dummy,&d_r);
	d_dummy.value = -Eta*vds;
	d_dummy.d1_r = -Eta;
	TimesDeriv(&d_Von,&d_Vpb,-K2);
	PlusDeriv(&d_Von,&d_Von,&d_dummy);
	TimesDeriv(&d_dummy,&d_SqrtVpb,K1);
	PlusDeriv(&d_Von,&d_dummy,&d_Von);
	d_Von.value = Von;
	Vth = Von;
	EqualDeriv(&d_Vth,&d_Von);
	VgsVth  =  vgs - Vth;
	TimesDeriv(&d_VgsVth,&d_Vth,-1.0);
	d_VgsVth.value = VgsVth;
	d_VgsVth.d1_p += 1.0;

	G  =   1./(1.744+0.8364 * Vpb);
	TimesDeriv(&d_G,&d_Vpb,0.8364);
	d_G.value += 1.744;
	InvDeriv(&d_G,&d_G);
	G = 1 - G;
	TimesDeriv(&d_G,&d_G,-1.0);
	d_G.value += 1.0;
	A  =  G/SqrtVpb;
	DivDeriv(&d_A,&d_G,&d_SqrtVpb);
	A = 1.0 + 0.5*K1*A;
	TimesDeriv(&d_A,&d_A,0.5*K1);
	d_A.value += 1.0;
	A = MAX( A, 1.0);   /* Modified */
	if (A <= 1.0) {
	  EqualDeriv(&d_A,&d_zero);
	  d_A.value = 1.0;
	}
	Arg  = MAX(( 1 + Ugs * VgsVth), 1.0);
	MultDeriv(&d_dummy,&d_Ugs,&d_VgsVth);
	d_dummy.value += 1.0;
	if (d_dummy.value <= 1.0) {
	  EqualDeriv(&d_Arg,&d_zero);
	  d_Arg.value = 1.0;
	} 
	else
	  EqualDeriv(&d_Arg,&d_dummy);

	if( VgsVth < 0 ) {
	  /* cutoff */
	  DrCur  = 0;
	  EqualDeriv(&d_DrCur,&d_zero);
	  goto SubthresholdComputation;
	}

	/* Quadratic Interpolation for Beta0 (Beta at vgs  =  0, vds=Vds) */

	Beta_Vds_0  =  (here->B1betaZero + here->B1betaZeroB * vbs);
	EqualDeriv(&d_Beta_Vds_0,&d_q);
	d_Beta_Vds_0.value = Beta_Vds_0;
	d_Beta_Vds_0.d1_q = here->B1betaZeroB;
	BVdd  =  (here->B1betaVdd + here->B1betaVddB * vbs);
	EqualDeriv(&d_BVdd,&d_q);
	d_BVdd.value = BVdd;
	d_BVdd.d1_q = here->B1betaVddB;
	dBVdddVds = MAX( here->B1betaVddD, 0.0);


	/* is the above wrong ?!? */



	if( vds > Vdd ) {
	  Beta0  =  BVdd + dBVdddVds * (vds - Vdd);
	  EqualDeriv(&d_Beta0,&d_r); d_Beta0.value = vds - Vdd;
	  TimesDeriv(&d_Beta0,&d_Beta0,dBVdddVds);
	  PlusDeriv(&d_Beta0,&d_Beta0,&d_BVdd);
	  /* dBVdddVds = const */


	  /* and this stuff here? */

	} else {
	  VddSquare  =  Vdd * Vdd; /* const */
	  C1  =  ( -BVdd + Beta_Vds_0 + dBVdddVds * Vdd) / VddSquare;
	  TimesDeriv(&d_C1,&d_BVdd,-1.0);
	  PlusDeriv(&d_C1,&d_C1,&d_Beta_Vds_0);
	  d_C1.value += dBVdddVds * Vdd;
	  TimesDeriv(&d_C1,&d_C1,1/VddSquare);

	  C2  =  2 * (BVdd - Beta_Vds_0) / Vdd - dBVdddVds;
	  TimesDeriv(&d_C2,&d_Beta_Vds_0,-1.0);
	  PlusDeriv(&d_C2,&d_C2,&d_BVdd);
	  TimesDeriv(&d_C2,&d_C2,2/Vdd);
	  d_C2.value -= dBVdddVds;
	  Beta0  =  (C1 * vds + C2) * vds + Beta_Vds_0;
	  MultDeriv(&d_Beta0,&d_r,&d_C1);
	  PlusDeriv(&d_Beta0,&d_Beta0,&d_C2);
	  MultDeriv(&d_Beta0,&d_Beta0,&d_dummy);
	  PlusDeriv(&d_Beta0,&d_Beta0,&d_Beta_Vds_0);

	  /*
	    dBeta0dVds  =  2*C1*vds + C2;
	    dBeta0dVbs  =  dC1dVbs * vds * vds + dC2dVbs * vds + dBeta_Vds_0_dVbs;
	    maybe we'll need these later */
	}

	/*Beta  =  Beta0 / ( 1 + Ugs * VgsVth );*/

	Beta = Beta0 / Arg ;
	DivDeriv(&d_Beta,&d_Beta0,&d_Arg);

	/*VdsSat  = MAX( VgsVth / ( A + Uds * VgsVth ),  0.0);*/

	Vc  =  Uds * VgsVth / A;
	DivDeriv(&d_Vc,&d_VgsVth,&d_A);
	MultDeriv(&d_Vc,&d_Vc,&d_Uds);

	if(Vc < 0.0 ) {
	  EqualDeriv(&d_Vc,&d_zero);
	  Vc=0.0;
	}
	Term1  =  sqrt( 1 + 2 * Vc );
	TimesDeriv(&d_Term1,&d_Vc,2.0);
	d_Term1.value += 1.0;
	SqrtDeriv(&d_Term1,&d_Term1);
	K  =  0.5 * ( 1 + Vc + Term1 );
	PlusDeriv(&d_K,&d_Vc,&d_Term1);
	d_K.value += 1.0;
	TimesDeriv(&d_K,&d_K,0.5);
	VdsSat = VgsVth / ( A * sqrt(K));
	if (VdsSat < 0.0) {
	  EqualDeriv(&d_VdsSat,&d_zero);
	  VdsSat = 0.0;
	}
	else
	  {
	    SqrtDeriv(&d_VdsSat,&d_K);
	    MultDeriv(&d_VdsSat,&d_VdsSat,&d_A);