home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 13 / 13.iso / p / p064 / 3.ddi / MOS3.C < prev    next >
Encoding:
C/C++ Source or Header  |  1991-07-01  |  15.0 KB  |  566 lines

  1. /****************************************************************************/
  2. /*        Copyright (C) 1986, 1987, 1988, 1989, 1990, 1991                  */
  3. /*          By MicroSim Corporation, All Rights Reserved                    */
  4. /****************************************************************************/
  5. /* mos3.c
  6.  *   $Revision:   1.8  $
  7.  *   $Author:   whb  $
  8.  *   $Date:   04 Mar 1991 13:55:22  $ */
  9.  
  10. #include "option1.h"
  11. #include "mos.fd"
  12.  
  13. /* these are set in MOS3eval and used there, and later in MOS3cap */
  14.  
  15. static double
  16.     Cox,
  17.     Eta,
  18.     Short,
  19.     Narrow,
  20.     Qspof;
  21.  
  22. #define PARAM(a)     ((double)Model->a)
  23. #define DELTA_MOS    PARAM(M_delta)
  24. #define ETA        PARAM(M_eta)
  25. #define GAMMA        PARAM(M_gamma)
  26. #define KAPPA        PARAM(M_kappa)
  27. #define KP        PARAM(M_kp)
  28. #define LAMBDA        PARAM(M_lambda)
  29. #define LD        PARAM(M_ld)
  30. #define NFS        PARAM(M_nfs)
  31. #define PHI        PARAM(M_phi)
  32. #define THETA        PARAM(M_theta)
  33. #define TOX        PARAM(M_tox)
  34. #define TYPE        (Model->M_type)
  35. #define UO        PARAM(M_uo)
  36. #define VBI        PARAM(M_vbi)
  37. #define VMAX        PARAM(M_vmax)
  38. #define XD        PARAM(M_xd)
  39. #define XJ        PARAM(M_xj)
  40. /* */
  41. #ifndef cc_noline
  42. #line 4 "MOS3eval"
  43. #endif
  44.  
  45. void MOS3eval(    /* evaluate MOS level 3 (semi-empirical, small geometry) model */
  46.   Model,            /* (pointer to) model */
  47.   Len, Wid, Mdev,        /* effective device size */
  48.   Vgs, Vds, Vbs,        /* terminal voltages */
  49.   Von, Vdsat,            /* (pointer to) device voltages (returned) */
  50.   Id,                /* (pointer to) drain current (returned) */
  51.   Gds, Gm, Gmbs,        /* (pointer to) device conductances (returned) */
  52.   Charge,            /* do charge calculations (YES/NO) */
  53.   Cggb, Cgdb, Cgsb,        /* (pointer to) device capacitances (returned) */
  54.   Cbgb, Cbdb, Cbsb,
  55.   Qgate, Qchan, Qbulk        /* (pointer to) device charges (returned) */
  56.   )
  57.   struct M_ *Model;
  58.   double Len, Wid, Mdev, Vgs, Vds, Vbs, *Von, *Vdsat, *Id, *Gds, *Gm, *Gmbs;
  59.   int Charge;
  60.   double *Cggb, *Cgdb, *Cgsb, *Cbgb, *Cbdb, *Cbsb, *Qgate, *Qchan, *Qbulk;
  61. /*
  62. 86-08-28 pwt    - creation
  63. 87-10-10 pwt    - re-work as part of BSIM addition
  64. 88-01-05 pwt    - parameter M_tox now oxide thickness, not Cox
  65. 89-04-15 pwt    - add device multiplier, to fix narrow-channel calculation
  66. 89-06-12 whjb   - Made use of EXP() macro
  67. 90-02-12 sv    - add local temperatures to models.
  68. *  whjb - 04 Mar 91 - Corrected bug in weak inversion code
  69. */
  70. { double
  71.     vgsx, vdsx,            /* voltages */
  72.     beta, gds0, gdsx, vdsat1,
  73.     vth,            /* threshold voltage & derivatives */
  74.     dvtdvd, dvtdvb,
  75.     vpof,            /* pinch-off voltage */
  76.     cdo, cd1, cdnorm,        /* drain currents factors */
  77.     dcodvb,            /* & derivatives */
  78.     onxl,            /* = 1/Len */
  79.     fgate, onfg, us,        /* mobility vs. gate voltage factors */
  80.     dfgdvg, dfgdvd, dfgdvb,    /* & derivatives */
  81.     vdsc, onvdsc,        /* saturation voltage factors */
  82.     dvsdvg, dvsdvd, dvsdvb,    /* & derivatives */
  83.     fdrain,            /* Vmax related currents */
  84.     dfddvg, dfddvd, dfddvb,    /* & derivatives */
  85.     xn, dxndvb,            /* weak inversion: Von = Vth + xn*VT */
  86.     dvodvb, dvodvd,
  87.     dfsdvb,            /* short-channel effect derivatives */
  88.     fbody, onfbdy, dfbdvb,    /* body effect & derivatives */
  89.     qbonco, dqbdvb,        /* bulk charge & derivatives */
  90.     phibs, sqphbs, dsqdvb;    /* factors & derivatives */
  91. #define    vt    VT
  92.  
  93.   static double        /* empirical constants */
  94.     coeff0 =  0.0631353,
  95.     coeff1 =  0.8013292,
  96.     coeff2 = -0.01110777;
  97.  
  98.  
  99. /* note: charge flag == NO ... charge calculations are never done!
  100.  *
  101.  * after tests are running, (maybe) use code similar to MOS2 flag setting
  102.  * NOTE: I have noticed problems with the values of VPOF and VDSAT1
  103.  */
  104.  
  105.   beta   = KP*Wid/Len;
  106.   Cox    = Wid*Len*EPSOX/TOX;
  107.   Narrow = DELTA_MOS/(Wid/Mdev);
  108.  
  109. /* reference Id equations to source and charge equations to bulk */
  110.  
  111.   *Vdsat = *Qgate = *Qchan = *Qbulk = *Cgdb = *Cbdb = 0.;
  112.  
  113.   onxl = 1/Len;
  114.   Eta  = ETA*onxl*onxl*onxl;
  115.  
  116. /* square root term */
  117.  
  118.   if ( Vbs > 0. ) {
  119.     double sqphis = sqrt(PHI);
  120.  
  121.     sqphbs = sqphis/(1 + Vbs/(PHI+PHI) );
  122.     phibs  = sqphbs*sqphbs;
  123.     dsqdvb = -.5*phibs/(PHI*sqphis);
  124.     }
  125.   else {
  126.     phibs  = PHI-Vbs;
  127.     sqphbs = sqrt(phibs);
  128.     dsqdvb = -.5/sqphbs;
  129.     }
  130.  
  131. /* short-channel effect factor */
  132.  
  133.   if ( XJ == 0. || XD == 0. ) {
  134.     dfsdvb = 0.;
  135.     Short  = 1.;
  136.     }
  137.   else {
  138.     double
  139.       wps    = XD*sqphbs,
  140.       onxj   = 1/XJ,
  141.       xjonxl = XJ*onxl,
  142.       djonxj = LD*onxj,
  143.       wponxj = wps*onxj,
  144.       wconxj = (coeff2*wponxj + coeff1)*wponxj + coeff0,
  145.       arga   = wconxj + djonxj,
  146.       argc   = wponxj/( 1 + wponxj ),
  147.       argb   = sqrt( 1 - argc*argc ),
  148.       dwpdvb = XD*dsqdvb,
  149.       dadvb  = ( coeff1 + coeff2*( wponxj + wponxj ))*dwpdvb*onxj,
  150.       dbdvb  = -argc*argc*( 1 - argc )*dwpdvb/(argb*wps);
  151.  
  152.     dfsdvb = -xjonxl*( dadvb*argb + arga*dbdvb );
  153.     Short  = 1 - xjonxl*( arga*argb - djonxj );
  154.     }
  155.  
  156.   { double                /* threshold voltage effects */
  157.       gammas, fbodys, vbix;
  158.  
  159. /* body effect*/
  160.  
  161.     gammas = GAMMA*Short;
  162.     fbodys = .25*gammas/sqphbs;
  163.     fbody  = fbodys + Narrow;
  164.     onfbdy = 1/( 1 + fbody );
  165.     dfbdvb = -fbodys*dsqdvb/sqphbs + fbodys*dfsdvb/Short;
  166.     qbonco = gammas*sqphbs + Narrow*phibs;
  167.     dqbdvb = gammas*dsqdvb + GAMMA*dfsdvb*sqphbs - Narrow;
  168.  
  169. /* static feedback effect */
  170.  
  171.     vbix = TYPE*VBI - Eta*Vds;
  172.  
  173. /* threshold voltage */
  174.  
  175.     vth    = vbix + qbonco;
  176.     dvtdvd = -Eta;
  177.     dvtdvb = dqbdvb;
  178.     }
  179.  
  180. /* joint weak inversion and strong inversion */
  181.  
  182.   *Von = vth;
  183.   if ( NFS == 0.) {
  184.  
  185.     if ( Vgs <= *Von ) {        /* cut-off region */
  186.       *Id = *Gm = *Gds = *Gmbs = 0.;
  187.  
  188.       if ( Charge == YES ) {
  189.         MOS3cap(Model,Vds,Vbs,vpof,*Von,vdsat1,Cggb,Cgdb,Cgsb,Cbgb,Cbdb,Cbsb,Qgate,Qchan,Qbulk);
  190. /***                          ^^^^      ^^^^^^ undefined here ***/
  191.         }
  192.       Qspof = 0.;
  193.  
  194.       return;        /*** done ***/
  195.     } }
  196.   else {
  197.     double
  198. /*      csonco = CHARGE*NFS*Len*Wid/Cox,*/
  199. /*      csonco = CHARGE*NFS/COX,*/
  200.       csonco = CHARGE*NFS*TOX/EPSOX,
  201.       cdonco = .5*qbonco/phibs;
  202.  
  203.     xn     = 1 + csonco + cdonco;
  204.     *Von   = vth + vt*xn;
  205.     dxndvb = .5*dqbdvb/phibs - qbonco*dsqdvb/(phibs*sqphbs);
  206.     dvodvd = dvtdvd;
  207.     dvodvb = dvtdvb + vt*dxndvb;
  208.     }
  209.  
  210. /* device is on */
  211.  
  212.   vgsx = MAX(Vgs,*Von);
  213.  
  214. /* mobility modulation by gate voltage */
  215.  
  216.   onfg   = 1 + THETA*(vgsx-vth);
  217.   fgate  = 1/onfg;
  218.   us     =  UO*fgate;
  219.   dfgdvg = -THETA*fgate*fgate;
  220.   dfgdvd = -dfgdvg*dvtdvd;
  221.   dfgdvb = -dfgdvg*dvtdvb;
  222.  
  223. /* saturation voltage */
  224.  
  225.   *Vdsat =
  226.   vpof   = (vgsx-vth)*onfbdy;
  227.  
  228.   if ( VMAX > 0.) {
  229.     double arga, argb, dvsdga;
  230.  
  231.     vdsc   = Len*VMAX/us;
  232.     onvdsc = 1/vdsc;
  233.     arga   = ( vgsx - vth )*onfbdy;
  234.     argb   = sqrt( arga*arga + vdsc*vdsc );
  235.     *Vdsat = arga + vdsc - argb;
  236.     dvsdga = ( 1 - arga/argb )*onfbdy;
  237.     dvsdvg = dvsdga - ( 1 - vdsc/argb )*vdsc*dfgdvg*onfg;
  238.     dvsdvd = -dvsdvg*dvtdvd;
  239.     dvsdvb = -dvsdvg*dvtdvb - arga*dvsdga*dfbdvb;
  240.     }
  241.   else {
  242.     dvsdvg = onfbdy;
  243.     dvsdvd = -dvsdvg*dvtdvd;
  244.     dvsdvb = -dvsdvg*dvtdvb - (*Vdsat)*dfbdvb*onfbdy;
  245.     }
  246.  
  247. /* current factors in linear region */
  248.  
  249.   vdsx = MIN(Vds,*Vdsat);
  250.   if ( vdsx == 0. ) {            /* special case of Vds = 0 */
  251.     *Id = *Gm = *Gmbs = 0.;
  252.     beta *= fgate;
  253.     *Gds  = beta*(vgsx-vth);
  254.     if ( NFS != 0. && Vgs < *Von )
  255.       *Gds *= EXP((Vgs - *Von)/(vt*xn));
  256.  
  257.     if ( Charge == YES )
  258.       MOS3cap(Model,Vds,Vbs,vpof,*Von,vdsat1,Cggb,Cgdb,Cgsb,Cbgb,Cbdb,Cbsb,Qgate,Qchan,Qbulk);
  259. /***                                  ^^^^^^ undefined here ***/
  260.     Qspof = 0.;
  261.  
  262.     return;    /* done */
  263.     }
  264.  
  265.   cdo    = vgsx - vth - .5*( 1 + fbody )*vdsx;
  266.   dcodvb = -dvtdvb - .5*dfbdvb*vdsx;
  267.  
  268. /* normalized drain current */
  269.  
  270.   cdnorm = cdo*vdsx;
  271.   *Gm    = vdsx;
  272.   *Gds   = -dvtdvd*vdsx;
  273.   gdsx   = vgsx - vth - ( 1 + fbody )*vdsx;
  274.   *Gmbs  = dcodvb*vdsx;
  275.  
  276. /* drain current without velocity saturation effect */
  277.  
  278.   cd1   = beta*cdnorm;
  279.   beta *= fgate;
  280.   *Id   = beta*cdnorm;
  281.   *Gm   = beta*(*Gm)  + dfgdvg*cd1;
  282.   *Gds  = beta*(*Gds);
  283.   gdsx  = beta*(gdsx) + dfgdvd*cd1;
  284.   *Gmbs = beta*(*Gmbs);
  285.  
  286. /* velocity saturation factor */
  287.  
  288.   if ( VMAX != 0.) {
  289.     double arga, fd2;
  290.  
  291.     fdrain = 1/( 1 + vdsx*onvdsc );
  292.     fd2    = fdrain*fdrain;
  293.     arga   = fd2*vdsx*onvdsc*onfg;
  294.     dfddvg = -dfgdvg*arga;
  295.     dfddvd = -dfgdvd*arga - fd2*onvdsc;
  296.     dfddvb = -dfgdvb*arga;
  297.  
  298. /* drain current with velocity saturation effect */
  299.  
  300.     *Gm   = fdrain*(*Gm)   + dfddvg**Id;
  301.     *Gds  = fdrain*(*Gds);
  302.     gdsx  = fdrain*(gdsx)  + dfddvd**Id;
  303.     *Gmbs = fdrain*(*Gmbs) + dfddvb**Id;
  304.     *Id   = fdrain*(*Id);
  305.     beta *= fdrain;
  306.     }
  307.  
  308. /* channel length modulation */
  309.  
  310.   if ( Vds > *Vdsat && (VMAX == 0. || LAMBDA != 0.) ) {
  311.     double delxl, dldvd, ddldvg, ddldvd, ddldvb;
  312.  
  313.     if ( VMAX == 0.) {
  314.       delxl  = sqrt(KAPPA*(Vds - *Vdsat)*LAMBDA);
  315.       dldvd  = .5*delxl/(Vds - *Vdsat);
  316.       ddldvg = 0.;
  317.       ddldvd = -dldvd;
  318.       ddldvb = 0.;
  319.       }
  320.     else {
  321.       double
  322.     idsat,
  323.         gdsat, gdoncd, gdonfd, gdonfg,
  324.         emax,  emoncd, emongd,
  325.         dgdvg,  dgdvd,  dgdvb,
  326.         demdvg, demdvd, demdvb, dldem,
  327.         arga, argb, argc;
  328.  
  329.       idsat  = *Id;
  330.       gdsat  = idsat*( 1 - fdrain )*onvdsc;
  331.       gdsat  = MAX(1E-12,gdsat);
  332.       gdoncd = gdsat/idsat;
  333.       gdonfd = gdsat/( 1 - fdrain );
  334.       gdonfg = gdsat*onfg;
  335.       dgdvg  = gdoncd*(*Gm)   - gdonfd*dfddvg + gdonfg*dfgdvg;
  336.       dgdvd  = gdoncd*(*Gds)  - gdonfd*dfddvd + gdonfg*dfgdvd;
  337.       dgdvb  = gdoncd*(*Gmbs) - gdonfd*dfddvb + gdonfg*dfgdvb;
  338.  
  339.       emax   = idsat*onxl/gdsat;
  340.       emoncd = emax/idsat;
  341.       emongd = emax/gdsat;
  342.       demdvg = emoncd*(*Gm)   - emongd*dgdvg;
  343.       demdvd = emoncd*(*Gds)  - emongd*dgdvd;
  344.       demdvb = emoncd*(*Gmbs) - emongd*dgdvb;
  345.  
  346.       arga   = .5*emax*LAMBDA;
  347.       argc   = KAPPA*LAMBDA;
  348.       argb   = sqrt( arga*arga + argc*(Vds - *Vdsat) );
  349.       delxl  = argb - arga;
  350.       dldvd  = argc/(argb+argb);
  351.       dldem  = .5*( arga/argb - 1 )*LAMBDA;
  352.       ddldvg = dldem*demdvg;
  353.       ddldvd = dldem*demdvd - dldvd;
  354.       ddldvb = dldem*demdvb;
  355.       }
  356.  
  357. /* punch through approximation */
  358.  
  359.     if ( delxl > .5*Len ) {
  360.       double arga;
  361.  
  362.       delxl  = Len - .25*Len*Len/delxl;
  363.       arga   = 4*(Len-delxl)*(Len-delxl)/(Len*Len);
  364.       ddldvg *= arga;
  365.       ddldvd *= arga;
  366.       ddldvb *= arga;
  367.       dldvd  *= arga;
  368.       }
  369.  
  370. /* saturation region */
  371.  
  372.     { double diddl,
  373.         dlonxl = delxl*onxl,
  374.         xlfact = 1/( 1 - dlonxl );
  375.  
  376.       *Id   *= xlfact;
  377.       diddl  = *Id/(Len-delxl);
  378.       *Gm    = xlfact*(*Gm)   + diddl*ddldvg;
  379.       *Gds   = xlfact*(*Gds);
  380.       gds0   = xlfact*(gdsx)  + diddl*ddldvd;
  381.       *Gmbs  = xlfact*(*Gmbs) + diddl*ddldvb;
  382.       *Gm   += gds0*dvsdvg;
  383.       *Gmbs += gds0*dvsdvb;
  384.       gdsx   = gds0*dvsdvd + diddl*dldvd;
  385.       }
  386.  
  387.     }
  388.   *Gds += gdsx;
  389.  
  390. /* finish strong inversion case */
  391.  
  392.   if ( Vgs < *Von ) {            /* weak inversion */
  393.     double gms, gmw,
  394.     onxn  = 1/xn,
  395.     argg  = onxn/vt,
  396.     argt  = (Vgs - *Von)*argg,
  397.     expg  = EXP(argt),
  398.     cdson = *Id;
  399.  
  400.     *Id = cdson*expg;
  401.     gmw = *Id*argg;
  402.     gms = *Gm*expg;
  403.     *Gm = gmw;
  404.     if ( Vds > *Vdsat )
  405.       *Gm += gds0*dvsdvg*(*Id/cdson);
  406.  
  407.     *Gds  = (*Id/cdson)*(*Gds)  + (gms-gmw)*dvodvd;
  408.     *Gmbs = (*Id/cdson)*(*Gmbs) + (gms-gmw)*dvodvb - gmw*(Vgs - *Von)*onxn*dxndvb;
  409.     }
  410.  
  411. /* charge computation */
  412.  
  413.   if ( Charge == NO ) Qspof = 0.;
  414.   else {
  415.     if ( Vgs > vth ) {
  416.       Mqspof(Model,Vds,Vbs,Vgs,vpof,*Von,*Vdsat,vdsat1,Cggb,Cgdb,Cgsb,Cbgb,Cbdb,Cbsb,Qgate,Qchan,Qbulk);
  417. /***                                            ^^^^^^ undefined here ***/
  418.       }
  419.     else {
  420.       MOS3cap(Model,Vds,Vbs,vpof,*Von,vdsat1,Cggb,Cgdb,Cgsb,Cbgb,Cbdb,Cbsb,Qgate,Qchan,Qbulk);
  421. /***                                  ^^^^^^ undefined here ***/
  422.       Qspof = 0.;
  423.     } }
  424.  
  425.   } /* done */
  426.  
  427. /* */
  428. #ifndef cc_noline
  429. #line 4 "MOS3cap"
  430. #endif
  431.  
  432. void
  433.   MOS3cap(Model, Vds, Vbs, Vgs, Von, Vdsat,
  434.     Cggb, Cgdb, Cgsb,
  435.     Cbgb, Cbdb, Cbsb,
  436.     Qgate, Qchan, Qbulk
  437.     )
  438.  
  439. /* calculate MOS level 3 capacitance and charge */
  440.  
  441.   struct M_    *Model;
  442.  
  443.   double
  444.     Vds, Vbs, Vgs, Von, Vdsat;
  445.  
  446.   double
  447.     *Qgate, *Qchan, *Qbulk,
  448.     *Cggb, *Cgdb, *Cgsb,
  449.     *Cbgb, *Cbdb, *Cbsb;
  450.  
  451. /* Author
  452.   PWT    86-08-29    creation
  453.   PWT    87-10-10    re-work as part of BSIM addition
  454. */
  455. { double vgb, vfb, phibs;
  456.   double gammas, qbonco;
  457.   double vbix, vth;
  458.   double vbi = TYPE*VBI;
  459.  
  460. /*** these are never initialized ***/
  461.   double sqphbs, dsqdvb, dfsdvb;
  462.  
  463. /* charge equations are referenced to bulk */
  464.  
  465.   vgb = Vgs - Vbs;
  466.   vfb = vbi - PHI;
  467. /*  onxl   = 1/Len,    not used */
  468.   phibs  = sqphbs*sqphbs;
  469.  
  470. /* body effect */
  471.  
  472.   gammas = GAMMA*Short;
  473.   qbonco = gammas*sqphbs + Narrow*phibs;
  474.  
  475. /* static feedback effect */
  476.  
  477.   vbix = vbi - Eta*Vds;
  478.  
  479. /* threshold voltage */
  480.  
  481.   vth = vbix + qbonco;
  482.  
  483. /* branch according to region of operation */
  484.  
  485.   if ( Vgs > vth ) {
  486.     double
  487.       fbodys = .5*gammas/( sqphbs + sqphbs ),
  488.       fbody  = fbodys + Narrow,
  489.       onfbdy = 1/( 1 + fbody ),
  490.       dfbdvb = -fbodys*dsqdvb/sqphbs + fbodys*dfsdvb/Short,
  491.       dqbdvb = gammas*dsqdvb + GAMMA*dfsdvb*sqphbs - Narrow,
  492.       dvtdvd = -Eta,
  493.       dvtdvb = dqbdvb,
  494.       vgsx   = MAX(Vgs,Von),
  495.       vdsx   = MIN(Vds,Vdsat);
  496.  
  497.     if ( vdsx == 0. ) {            /* special case */
  498.       *Qgate = Cox*(Vgs-vbi);
  499.       *Qbulk = -Cox*qbonco;
  500.       *Cggb  = Cox;
  501.       *Cgdb  = -Cox*( .5 + dvtdvd );
  502.       *Cgsb  = -Cox*( .5 - dvtdvb );
  503.       *Cbgb  = 0.;
  504.       *Cbdb  = -.5*Cox*fbody;
  505.       *Cbsb  = Cox*( dqbdvb + .5*fbody );
  506.       }
  507.     else {
  508.       double cdo, dcodvg, dcodvd, dcodvb;
  509.       double arga, dadco, dadvd, dadfb;
  510.  
  511.       cdo    = vgsx - vth - .5*(1+fbody)*vdsx;
  512.       dcodvg = 1.;
  513.       if ( Vds < Vdsat ) dcodvd = -dvtdvd - .5*(1+fbody);
  514.       dcodvb = -dvtdvb - .5*dfbdvb*vdsx;
  515.  
  516. /* charge terms */
  517.  
  518.       arga  = (1+fbody)*vdsx*vdsx/(12*cdo);
  519.       dadco = -arga/cdo;
  520.       if ( Vds < Vdsat ) dadvd = arga/vdsx;
  521.       dadfb = arga*onfbdy;
  522.  
  523. /* gate charge */
  524.  
  525.       *Qgate = Cox*( Vgs - vbix - .5*vdsx + arga );
  526.       *Cggb  = Cox*( 1 + dadco*dcodvg );
  527.       if ( Vds < Vdsat ) *Cgdb = Cox*( -dvtdvd - .5 + dadvd + dadco*dcodvd);
  528.       *Cgsb  = -*Cggb - *Cgdb - Cox*( dadco*dcodvb + dadfb*dfbdvb );
  529.  
  530. /* bulk charge */
  531.  
  532.       arga  = arga*fbody;
  533.       dadco = dadco*fbody;
  534.       if ( Vds < Vdsat ) dadvd = dadvd*fbody;
  535.       dadfb = dadfb*( 1 + fbody + fbody );
  536.  
  537.       *Qbulk = -Cox*( qbonco + .5*fbody*vdsx - arga );
  538.       *Cbgb  = Cox*dadco*dcodvg;
  539.       if ( Vds < Vdsat ) *Cbdb = -Cox*( .5*fbody - dadvd - dadco*dcodvd );
  540.       *Cbsb  = -*Cbgb - *Cbdb + Cox*( dqbdvb + ( .5*vdsx - dadfb )*dfbdvb - dadco*dcodvb );
  541.       }
  542.     }
  543.   else {                /* charge terms of Vgs < vth */
  544.     if ( vgb > vfb ) {
  545.       double
  546.         gamma2 = .5*gammas,
  547.         arga   = sqrt( gamma2*gamma2 + (vgb-vfb) );
  548.  
  549.       *Qgate = gammas*Cox*( arga - gamma2 );
  550.       *Cggb  = .5*Cox*gammas/arga;
  551.       }
  552.     else {
  553.       *Qgate = Cox*(vgb-vfb);
  554.       *Cggb  = Cox;
  555.       }
  556.     *Qbulk = -*Qgate;
  557.     *Cbgb  = -*Cggb;
  558.     *Cgdb  =  *Cgsb = *Cbdb = *Cbsb = 0.;
  559.     }
  560.  
  561.   *Qchan = -( *Qgate + *Qbulk );
  562.  
  563. /* done */
  564.  
  565.   }
  566.