apps/d_mos4.model

   1 /* $Id: d_mos4.model,v 26.92 2008/08/23 05:40:00 al Exp $ -*- C++ -*-
   2  * Copyright (C) 2001 Albert Davis
   3  * Author: Albert Davis <aldavis@gnu.org>
   4  *
   5  * This file is part of "Gnucap", the Gnu Circuit Analysis Package
   6  *
   7  * This program is free software; you can redistribute it and/or modify
   8  * it under the terms of the GNU General Public License as published by
   9  * the Free Software Foundation; either version 3, or (at your option)
  10  * any later version.
  11  *
  12  * This program is distributed in the hope that it will be useful,
  13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  15  * GNU General Public License for more details.
  16  *
  17  * You should have received a copy of the GNU General Public License
  18  * along with this program; if not, write to the Free Software
  19  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
  20  * 02110-1301, USA.
  21  *------------------------------------------------------------------
  22  * MOS BSIM1 model.
  23  * derived from Spice3f4,Copyright 1990 Regents of the University of California
  24  * 1985 Hong J. Park, Thomas L. Quarles
  25  * Recoded for Gnucap model compiler, Al Davis, 2000
  26  */
  27 h_headers {
  28 #include "d_mos_base.h"
  29 }
  30 cc_headers {
  31 }
  32 /*--------------------------------------------------------------------------*/
  33 model BUILT_IN_MOS4 {
  34   level 4;
  35   public_keys {
  36     nmos4 polarity=pN;
  37     pmos4 polarity=pP;
  38   }
  39   dev_type BUILT_IN_MOS;
  40   inherit BUILT_IN_MOS_BASE;
  41   independent {
  42     override {
  43       double mjsw "" final_default=.33;
  44       double pb "" final_default=0.1 quiet_min=0.1;
  45       double pbsw "" final_default=pb quiet_min=0.1;
  46       double cjo "" default=0.0;
  47       int cmodel "CMODEL" print_test="cmodel!=1"
  48         calculate="((!cmodel)?1:cmodel)";
  49       int mos_level "back-annotate for diode" name=DIODElevel
  50         print_test="mos_level != LEVEL" default=LEVEL;
  51     }
  52     raw_parameters {
  53       double dl_u "Channel length reduction"
  54         name=DL default=0.;
  55       double dw_u "Channel width reduction"
  56         name=DW default=0.;
  57       double tox_u "Gate oxide thickness"
  58         name=TOX default=0. quiet_min=1e-20;
  59       double vdd "Supply voltage to specify mus"
  60         name=VDD default=0.;
  61       double wdf "Default width of source drain diffusion (ignored)"
  62         name=WDF default=0.;
  63       double dell "Length reduction of source drain diff (ignored)"
  64         name=DELL default=0.;
  65       double temp "temperature (ignored)"
  66         name=TEMP default=300.15;
  67       double xpart "Flag for channel charge partitioning"
  68         name=XPART default=0.;
  69     }
  70     calculated_parameters {
  71       double dl "" calculate="dl_u*MICRON2METER";
  72       double dw "" calculate="dw_u*MICRON2METER";
  73       double tox "" calculate="tox_u*MICRON2METER";
  74       double cox "" calculate="3.453e-11 /*E_OX*/ / tox";
  75     }
  76   }
  77   size_dependent {
  78     raw_parameters {
  79       double phi "Strong inversion surface potential"
  80         name=PHI default=0. positive quiet_min=.1;
  81       double vfb "Flat band voltage"
  82         name=VFB default=0.;
  83       double k1 "Bulk effect coefficient 1"
  84         name=K1 default=0. positive;
  85       double k2 "Bulk effect coefficient 2"
  86         name=K2 default=0. positive;
  87       double eta "VDS dependence of threshold voltage"
  88         name=ETA default=0.;
  89       double etaB "VBS dependence of eta (x2e)"
  90         name=X2E default=0.;
  91       double etaD "VDS dependence of eta (x3e)"
  92         name=X3E default=0.;
  93       double mobZero "Zero field mobility at VDS=0 VGS=VTH (muz)"
  94         name=MUZ default=0.;
  95       double mobZeroB "VBS dependence of muz (x2mz)"
  96         name=X2MZ default=0.;
  97       double mobVdd "Mobility @ VDS=VDD VGS=VTH, chan len modulation"
  98         name=MUS default=0.;
  99       double mobVddB "VBS dependence of mus (x2ms)"
 100         name=X2MS default=0.;
 101       double mobVddD "VDS dependence of mus (x3ms)"
 102         name=X3MS default=0.;
 103       double ugs "VGS dependence of mobility (u0)"
 104         name=U0 default=0.;
 105       double ugsB "VBS dependence of u0 (x2u0)"
 106         name=X2U0 default=0.;
 107       double uds "VDS dependence of mobility, velocity saturation"
 108         name=U1 default=0.;
 109       double udsB "VBS dependence of u1 (x2u1)"
 110         name=X2U1 default=0.;
 111       double udsD "VDS dependence of u1 (x3u1)"
 112         name=X3U1 default=0.;
 113       double n0 "Subthreshold slope (n0)"
 114         name=N0 default=0.;
 115       double nB "VBS dependence of subthreshold slope (nb)"
 116         name=NB default=0.;
 117       double nD "VDS dependence of subthreshold slope (nd)"
 118         name=ND default=0.;
 119     }
 120     calculated_parameters {
 121       double betaZero  "Beta at vds = 0 and vgs = Vth"
 122         calculate="mobZero  * CoxWoverL";
 123       double betaZeroB "Vbs dependence of BetaZero"
 124         calculate="mobZeroB * CoxWoverL";
 125       double betaVdd   "Beta at vds=Vdd and vgs=Vth"
 126         calculate="mobVdd   * CoxWoverL";
 127       double betaVddB  "Vbs dependence of BVdd"
 128         calculate="mobVddB  * CoxWoverL";
 129       double betaVddD  "Vds dependence of BVdd"
 130         calculate="mobVddD  * CoxWoverL";
 131       double vt0 "" calculate="vfb + phi + k1 * sqrt(phi) - k2 * phi";
 132     }
 133     code_pre {
 134       l_eff -= m->dl;
 135       w_eff -= m->dw;
 136       double L = l_eff/MICRON2METER;
 137       double W = w_eff/MICRON2METER;
 138       double CoxWoverL = 1e-4 * m->cox * w_eff / l_eff;
 139     }
 140     override {
 141       double cgate "" calculate="m->cox * w_eff * l_eff";
 142     }
 143   }
 144   /*-----------------------------------------------------------------------*/
 145   tr_eval {
 146     trace3("", d->vds, d->vgs, d->vbs);
 147     /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 148     d->reverse_if_needed();
 149     /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 150     trace2("", m->dl, m->dw);
 151     trace3("", s->ugs, s->ugsB, d->vbs);
 152     double Ugs = s->ugs + s->ugsB * d->vbs;
 153     double dUgsdVbs;
 154     if(Ugs <= 0) {
 155       untested();
 156       Ugs = 0;
 157       dUgsdVbs = 0.0;
 158     }else{
 159       dUgsdVbs = s->ugsB;
 160     }
 161     trace2("", Ugs, dUgsdVbs);
 162
 163     double Uds = s->uds + s->udsB * d->vbs + s->udsD * (d->vds - m->vdd);
 164     double dUdsdVbs;
 165     double dUdsdVds;
 166     if(Uds <= 0) {
 167       untested();
 168       Uds = 0.0;
 169       dUdsdVbs = dUdsdVds = 0.0;
 170     }else{
 171       double Leff = s->l_eff * 1e6; /* Leff in um */
 172       Uds  =  Uds / Leff;
 173       dUdsdVbs = s->udsB / Leff;
 174       dUdsdVds = s->udsD / Leff;
 175     }
 176     trace3("", Uds, dUdsdVbs, dUdsdVds);
 177
 178     double Vpb;
 179     if(d->vbs <= 0) {
 180       Vpb = s->phi - d->vbs;
 181       d->sbfwd = false;
 182     }else{
 183       Vpb = s->phi;
 184       d->sbfwd = true;
 185     }
 186     double SqrtVpb = sqrt(Vpb);
 187     trace2("", Vpb, SqrtVpb);
 188     /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 189     /* threshold voltage */
 190     double dVthdVbs;
 191     double dVthdVds;
 192     {
 193       double Eta = s->eta + s->etaB * d->vbs + s->etaD * (d->vds - m->vdd);
 194       double dEtadVds;
 195       double dEtadVbs;
 196       if(Eta <= 0) {
 197         Eta  = 0;
 198         dEtadVds = dEtadVbs = 0.0 ;
 199       }else if (Eta > 1) {
 200         untested();
 201         Eta = 1;
 202         dEtadVds = dEtadVbs = 0;
 203       }else{
 204         untested();
 205         dEtadVds = s->etaD;
 206         dEtadVbs = s->etaB;
 207       }
 208       trace3("", Eta, dEtadVds, dEtadVbs);
 209       d->von = s->vfb + s->phi + s->k1 * SqrtVpb - s->k2 * Vpb - Eta * d->vds;
 210       dVthdVds = -Eta - dEtadVds * d->vds;
 211       dVthdVbs = s->k2 - 0.5 * s->k1 / SqrtVpb - dEtadVbs * d->vds;
 212       d->vgst  = d->vgs - d->von;
 213       trace4("", d->von, dVthdVds, dVthdVbs, d->vgst);
 214     }
 215     /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 216     double G = 1. - 1./(1.744+0.8364 * Vpb);
 217     double A = 1. + 0.5*G*s->k1/SqrtVpb;
 218     A = std::max(A, 1.0);   /* Modified */
 219     double Arg = std::max((1 + Ugs * d->vgst), 1.0);
 220     double dGdVbs = -0.8364 * (1-G) * (1-G);
 221     double dAdVbs = 0.25 * s->k1 / SqrtVpb *(2*dGdVbs + G/Vpb);
 222     trace3("", G, A, Arg);
 223     trace2("", dGdVbs, dAdVbs);
 224     /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 225     /* ids and derivatives calculation */
 226     if (d->vgst < 0) {
 227       d->cutoff = true;
 228       /* cutoff */
 229       d->ids = 0;
 230       d->gmf = 0;
 231       d->gds = 0;
 232       d->gmbf = 0;
 233       d->vdsat = 0;
 234       trace4("cutoff", d->ids, d->gmf, d->gds, d->gmbf);
 235     }else{
 236       d->cutoff = false;
 237       /* Quadratic Interpolation for Beta0 (Beta at d->vgs  =  0, vds=Vds) */
 238       double Beta0;
 239       double dBeta0dVds;
 240       double dBeta0dVbs;
 241       {
 242         trace2("", m->tox, m->cox*1e-4);
 243         trace3("", s->betaVdd, s->betaVddB, d->vbs);
 244         double BetaVdd = (s->betaVdd + s->betaVddB * d->vbs);
 245         double dBetaVdd_dVds = std::max(s->betaVddD, 0.0); /* Modified */
 246         trace2("", BetaVdd, dBetaVdd_dVds);
 247         if(d->vds > m->vdd) {
 248           Beta0      = BetaVdd + dBetaVdd_dVds * (d->vds - m->vdd);
 249           dBeta0dVds = dBetaVdd_dVds;
 250           dBeta0dVbs = s->betaVddB;
 251           trace3("vds>vdd", Beta0, dBeta0dVds, dBeta0dVbs);
 252         }else{
 253           double Beta_Vds_0 = (s->betaZero + s->betaZeroB * d->vbs);
 254           double VddSquare = m->vdd * m->vdd;
 255           double C1 = (-BetaVdd + Beta_Vds_0+dBetaVdd_dVds*m->vdd) / VddSquare;
 256           double C2 = 2 * (BetaVdd - Beta_Vds_0) / m->vdd - dBetaVdd_dVds;
 257           trace4("", Beta_Vds_0, VddSquare, C1, C2);
 258           double dBeta_Vds_0_dVbs = s->betaZeroB;
 259           double dBetaVdd_dVbs = s->betaVddB;
 260           double dC1dVbs = (dBeta_Vds_0_dVbs - dBetaVdd_dVbs) / VddSquare;
 261           double dC2dVbs = dC1dVbs * (-2) * m->vdd;
 262           trace4("", dBeta_Vds_0_dVbs, dBetaVdd_dVbs, dC1dVbs, dC2dVbs);
 263           Beta0    = (C1 * d->vds + C2) * d->vds + Beta_Vds_0;
 264           dBeta0dVds = 2*C1*d->vds + C2;
 265           dBeta0dVbs = dC1dVbs * d->vds * d->vds
 266             + dC2dVbs * d->vds + dBeta_Vds_0_dVbs;
 267           trace3("vds<vdd", Beta0, dBeta0dVds, dBeta0dVbs);
 268         }
 269       }
 270
 271       /*Beta  =  Beta0 / ( 1 + Ugs * d->vgst );*/
 272       double Beta            = Beta0 / Arg ;
 273       double dBetadVgs = -Beta * Ugs / Arg;
 274       double dBetadVds = dBeta0dVds / Arg - dBetadVgs * dVthdVds ;
 275       double dBetadVbs = dBeta0dVbs / Arg
 276         + Beta * Ugs * dVthdVbs / Arg - Beta * d->vgst * dUgsdVbs / Arg;
 277       trace4("", Beta, dBetadVgs, dBetadVds, dBetadVbs);
 278
 279       /*d->vdsat  = std::max(d->vgst / ( A + Uds * d->vgst ),  0.0);*/
 280       double Vc = Uds * d->vgst / A;
 281       if(Vc < 0.0) {
 282         untested();
 283         Vc=0.0;
 284       }
 285
 286       double Term1 = sqrt(1 + 2 * Vc);
 287       double K = 0.5 * (1 + Vc + Term1);
 288       d->vdsat = std::max(d->vgst / (A * sqrt(K)) , 0.0);
 289       trace4("", Vc, Term1, K, d->vdsat);
 290
 291       if(d->vds < d->vdsat) {           /* Triode Region */
 292         d->saturated = false;
 293         /*Argl1  =  1 + Uds * d->vds;*/
 294         double Argl1 = std::max((1 + Uds * d->vds), 1.);
 295         double Argl2 = d->vgst - 0.5 * A * d->vds;
 296         trace2("", Argl1, Argl2);
 297         d->ids = Beta * Argl2 * d->vds / Argl1;
 298         d->gmf   = (dBetadVgs * Argl2 * d->vds + Beta * d->vds) / Argl1;
 299         d->gds  = (dBetadVds * Argl2 * d->vds
 300                    + Beta * (d->vgst - d->vds * dVthdVds - A * d->vds)
 301                    - d->ids * (d->vds * dUdsdVds + Uds)) /  Argl1;
 302         d->gmbf = (dBetadVbs * Argl2 * d->vds
 303                   + Beta * d->vds * (-dVthdVbs - 0.5 * d->vds * dAdVbs)
 304                   - d->ids * d->vds * dUdsdVbs) / Argl1;
 305         trace4("triode", d->ids, d->gmf, d->gds, d->gmbf);
 306       }else{                    /* Pinchoff (Saturation) Region */
 307         d->saturated = true;
 308         double Args1   =  1. + 1. / Term1;
 309         double dVcdVgs =  Uds / A;
 310         double dVcdVds =  d->vgst * dUdsdVds / A - dVcdVgs * dVthdVds;
 311         double dVcdVbs =  (d->vgst * dUdsdVbs
 312                            - Uds * (dVthdVbs + d->vgst * dAdVbs / A )) / A;
 313         double dKdVc   =  0.5 * Args1;
 314         double dKdVgs  =  dKdVc * dVcdVgs;
 315         double dKdVds  =  dKdVc * dVcdVds;
 316         double dKdVbs  =  dKdVc * dVcdVbs;
 317         double Args2   =  d->vgst / A / K;
 318         double Args3   =  Args2 * d->vgst;
 319         trace3("", Args1, Args2, Args3);
 320         trace3("", dVcdVgs, dVcdVds, dVcdVbs);
 321         trace4("", dKdVc, dKdVgs, dKdVds, dKdVbs);
 322         d->ids =  0.5 * Beta * Args3;
 323         d->gmf = 0.5 * Args3 * dBetadVgs + Beta * Args2 - d->ids * dKdVgs / K;
 324         d->gds = 0.5*Args3*dBetadVds - Beta*Args2*dVthdVds - d->ids*dKdVds/K;
 325         d->gmbf = 0.5 * dBetadVbs * Args3 - Beta * Args2 *dVthdVbs
 326           - d->ids * (dAdVbs / A + dKdVbs / K);
 327         trace4("sat", d->ids, d->gmf, d->gds, d->gmbf);
 328       }
 329     }
 330     /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 331     /* SubthresholdComputation */
 332
 333     /* The following 'if' statement has been modified so that subthreshold  *
 334      * current computation is always executed unless N0 >= 200. This should *
 335      * get rid of the Ids kink seen on Ids-Vgs plots at low Vds.            *
 336      *                                                Peter M. Lee          *
 337      *                                                6/8/90                *
 338      *  Old 'if' statement:  (reversed)                                     *
 339      *  if( (N0 >=  200) || (d->vgst < Vcut ) || (d->vgst > (-0.5*Vcut)))   */
 340
 341     //double Vcut  =  - 40. * s->n0 * t->vt ;
 342     if (s->n0 < 200) {
 343       double N = s->n0 + s->nB*d->vbs + s->nD*d->vds; /* subthreshold slope */
 344       trace4("", s->n0, s->nB, s->nD, N);
 345       if (N < 0.5) {
 346         untested();
 347         N = 0.5;
 348       }
 349       const double ref_temp = 300.15; // ignore real temp for spice compatibility
 350       const double vt0 = ref_temp * P_K_Q;
 351       const double Vtsquare = vt0 * vt0 ;
 352       const double nvt0 = N * vt0;
 353       double Warg1 = exp(-d->vds / vt0);
 354       double Wds   = 1 - Warg1;
 355       double Wgs   = exp( d->vgst / nvt0);
 356       double Warg2  = 6.04965 * Vtsquare * s->betaZero;
 357       double Ilimit = 4.5 * Vtsquare * s->betaZero;
 358       double Iexp   = Warg2 * Wgs * Wds;
 359       d->ids += Ilimit * Iexp / (Ilimit + Iexp);
 360       double Temp1  = Ilimit / (Ilimit + Iexp);
 361       Temp1  =  Temp1 * Temp1;
 362       double Temp3  = Ilimit / (Ilimit + Wgs * Warg2);
 363       Temp3 = Temp3 * Temp3 * Warg2 * Wgs;
 364       /*    if ( Temp3 > Ilimit ) Temp3=Ilimit;*/
 365       d->gmf  += Temp1 * Iexp / nvt0;
 366       /* gds term has been modified to prevent blow up at Vds=0 */
 367       d->gds += Temp3
 368         * (Wds / nvt0 * (dVthdVds + d->vgst * s->nD / N)
 369            + Warg1 / vt0);
 370       d->gmbf -= Temp1 * Iexp * (dVthdVbs + d->vgst * s->nB / N) / nvt0;
 371       trace3("", vt0, Vtsquare, nvt0);
 372       trace4("", Warg1, Wds, Wgs, Warg2);
 373       trace4("", Ilimit, Iexp, Temp1, Temp3);
 374       trace4("sub", d->ids, d->gmf, d->gds, d->gmbf);
 375     }else{
 376       untested();
 377     }
 378     /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
 379     /* Some Limiting of DC Parameters */
 380     if(d->ids < 0.0) d->ids = 0.0;
 381     if(d->gmf  < 0.0) d->gmf  = 0.0;
 382     if(d->gds < 0.0) d->gds = 0.0;
 383     if(d->gmbf < 0.0) d->gmbf = 0.0;
 384     trace4("final", d->ids, d->gmf, d->gds, d->gmbf);
 385
 386     trace3("", G, A, s->phi);
 387     trace1("", m->xpart);
 388
 389     double Vth0 = s->vfb + s->phi + s->k1 * SqrtVpb; // almost same as d->von
 390     double Vgs_Vth = d->vgs - Vth0; // almost same as d->vgst
 391     trace2("", Vth0, Vgs_Vth);
 392     double Arg1 = A * d->vds;
 393     double Arg2 = Vgs_Vth - 0.5 * Arg1;
 394     double Arg3 = d->vds - Arg1;
 395     trace3("", Arg1, Arg2, Arg3);
 396     /*double*/ dVthdVbs = -0.5 * s->k1 / SqrtVpb;
 397     /*dbl*/ dAdVbs = 0.5 * s->k1 * (0.5*G/Vpb - 0.8364*(1-G)*(1-G)) / SqrtVpb;
 398     trace2("", dVthdVbs, dAdVbs);
 399     double Ent = std::max(Arg2,1.0e-8);
 400     double dEntdVds = -0.5 * A;
 401     double dEntdVbs = -dVthdVbs - 0.5 * d->vds * dAdVbs;
 402     trace3("", Ent, dEntdVds, dEntdVbs);
 403     double VdsPinchoff = std::max(Vgs_Vth / A, 0.0);
 404     double Vgb  = d->vgs - d->vbs ;
 405     double Vgb_Vfb  = Vgb - s->vfb;
 406     trace3("", VdsPinchoff, Vgb, Vgb_Vfb);
 407
 408     if(Vgb_Vfb < 0) {                   /* Accumulation Region */
 409       untested();
 410       d->qgate = s->cgate * Vgb_Vfb;
 411       d->qbulk = -d->qgate;
 412       d->qdrn  = 0. ;
 413       d->cggb = s->cgate;
 414       d->cgdb = 0.;
 415       d->cgsb = 0.;
 416       d->cbgb = -s->cgate;
 417       d->cbdb = 0.;
 418       d->cbsb = 0.;
 419       d->cdgb = 0.;
 420       d->cddb = 0.;
 421       d->cdsb = 0.;
 422       trace4("", d->qgate, d->cggb, d->cgdb, d->cgsb);
 423       trace4("", d->qbulk, d->cbgb, d->cbdb, d->cbsb);
 424       trace4("acc", d->qdrn, d->cdgb, d->cddb, d->cdsb);
 425     }else if (d->vgs < Vth0) {          /* Subthreshold Region */
 426       d->qgate = 0.5*s->cgate*s->k1*s->k1*(-1+sqrt(1+4*Vgb_Vfb/(s->k1*s->k1)));
 427       d->cggb = s->cgate / sqrt(1 + 4 * Vgb_Vfb / (s->k1 * s->k1));
 428       d->cgdb = d->cgsb = 0.;
 429       d->qbulk = -d->qgate;
 430       d->cbgb = -d->cggb;
 431       d->cbdb = d->cbsb = 0.0;
 432       d->qdrn = 0.;
 433       d->cdgb = d->cddb = d->cdsb = 0.0;
 434       trace4("", d->qgate, d->cggb, d->cgdb, d->cgsb);
 435       trace4("", d->qbulk, d->cbgb, d->cbdb, d->cbsb);
 436       trace4("sub", d->qdrn, d->cdgb, d->cddb, d->cdsb);
 437     }else if (d->vds < VdsPinchoff) {   /* triode region */
 438       double EntSquare = Ent * Ent;
 439       trace1("tri", EntSquare);
 440       double Argl1 = 1.2e1 * EntSquare;
 441       double Argl2 = 1.0 - A;
 442       double Argl3 = Arg1 * d->vds;
 443       trace3("", Argl1, Argl2, Argl3);
 444       double Argl5;
 445       if (Ent > 1.0e-8) {
 446         Argl5 = Arg1 / Ent;
 447       }else{
 448         untested();
 449         Argl5 = 2.0;
 450       }
 451       double Argl7 = Argl5 / 1.2e1;
 452       double Argl8 = 6.0 * Ent;
 453       trace3("", Argl5, Argl7, Argl8);
 454
 455       d->qgate = s->cgate
 456         * (d->vgs - s->vfb - s->phi - 0.5 * d->vds + d->vds * Argl7);
 457       d->cggb = s->cgate * (1.0 - Argl3 / Argl1);
 458       d->cgdb = s->cgate * (-0.5 + Arg1 / Argl8 - Argl3 * dEntdVds / Argl1);
 459       double cgbb = s->cgate * (d->vds*d->vds*dAdVbs*Ent-Argl3*dEntdVbs)/Argl1;
 460       d->cgsb = -(d->cggb + d->cgdb + cgbb);
 461       trace4("", d->qgate, d->cggb, d->cgdb, d->cgsb);
 462
 463       d->qbulk = s->cgate * (-Vth0 + s->vfb + s->phi + 0.5*Arg3 - Arg3*Argl7);
 464       d->cbgb = s->cgate * Argl3 * Argl2 / Argl1;
 465       d->cbdb = s->cgate * Argl2 * (0.5 - Arg1/Argl8 + Argl3 * dEntdVds/Argl1);
 466       double cbbb = -s->cgate * (dVthdVbs + 0.5 * d->vds * dAdVbs
 467                +d->vds*d->vds*((1.0-2.0*A)*dAdVbs*Ent-Argl2*A*dEntdVbs)/Argl1);
 468       d->cbsb = -(d->cbgb + d->cbdb + cbbb);
 469       trace4("", d->qbulk, d->cbgb, d->cbdb, d->cbsb);
 470
 471       if (m->xpart >= 1) {
 472         /*0/100 partitioning for drain/source chArges at saturation region*/
 473         double Argl9 = 0.125 * Argl5 * Argl5; //t
 474         d->qdrn = -s->cgate * (0.5*Vgs_Vth - 0.75*Arg1 + 0.125*Arg1*Argl5);
 475         d->cdgb = -s->cgate * (0.5 - Argl9);
 476         d->cddb = s->cgate * (0.75*A - 0.25*A*Arg1/Ent + Argl9*dEntdVds);
 477         double cdbb = s->cgate * (0.5 * dVthdVbs + d->vds * dAdVbs *
 478                                   (0.75 - 0.25 * Argl5 ) + Argl9 * dEntdVbs);
 479         d->cdsb = -(d->cdgb + d->cddb + cdbb);
 480         trace2("", Argl9, cdbb);
 481         trace4("tri 0/100", d->qdrn, d->cdgb, d->cddb, d->cdsb);
 482       }else{
 483         /*40/60 partitioning for drain/source chArges at saturation region*/
 484         double Vgs_VthSquare = Vgs_Vth*Vgs_Vth;
 485         trace2("", Vgs_Vth, Vgs_VthSquare);
 486         double Arg5 = Arg1*Arg1;
 487         double Vcom = Vgs_Vth*Vgs_Vth/6.0-1.25e-1*Arg1*Vgs_Vth+2.5e-2*Arg5;
 488         double Argl4 = Vcom/Ent/EntSquare;
 489         double Argl6;
 490         if (Ent > 1.0e-8) {
 491           Argl6 = Vcom / EntSquare;
 492         }else{
 493           untested();
 494           Argl6 = 4.0 / 1.5e1;
 495         }
 496         d->qdrn = -s->cgate * (0.5 * (Vgs_Vth-Arg1) + Arg1 * Argl6);
 497         d->cdgb = -s->cgate
 498           * (0.5 + Arg1*(4.0*Vgs_Vth-1.5*Arg1)/Argl1 - 2.0*Arg1*Argl4);
 499         d->cddb = s->cgate*(0.5*A+2.0*Arg1*dEntdVds*Argl4-A*(2.0*Vgs_VthSquare
 500                                         -3.0*Arg1*Vgs_Vth+0.9*Arg5)/Argl1);
 501         double cdbb =s->cgate*(0.5*dVthdVbs+0.5*d->vds*dAdVbs+2.0*Arg1*dEntdVbs
 502              *Argl4-d->vds*(2.0*Vgs_VthSquare*dAdVbs-4.0*A*Vgs_Vth*dVthdVbs-3.0
 503              *Arg1*Vgs_Vth*dAdVbs+1.5*A*Arg1*dVthdVbs+0.9*Arg5*dAdVbs)
 504                                /Argl1);
 505         d->cdsb = -(d->cdgb + d->cddb + cdbb);
 506         trace4("", Vcom, Argl4, Argl6, cdbb);
 507         trace4("lin 40/60", d->qdrn, d->cdgb, d->cddb, d->cdsb);
 508       }
 509     }else{                              /* saturation region */
 510       assert(d->vds >= VdsPinchoff);
 511       double Args1 = 1.0 / (3.0 * A);
 512       trace2("sat", s->cgate, Args1);
 513
 514       d->qgate = s->cgate * (d->vgs - s->vfb - s->phi - Vgs_Vth * Args1);
 515       d->cggb = s->cgate * (1.0 - Args1);
 516       d->cgdb = 0.0;
 517       double cgbb = s->cgate * Args1 * (dVthdVbs + Vgs_Vth * dAdVbs / A);
 518       d->cgsb = -(d->cggb + d->cgdb + cgbb);
 519       trace4("", d->qgate, d->cggb, d->cgdb, d->cgsb);
 520
 521       d->qbulk = s->cgate * (s->vfb + s->phi - Vth0 + (1.0-A)*Vgs_Vth*Args1);
 522       d->cbgb = s->cgate * (Args1 - 1.0 / 3.0);
 523       d->cbdb = 0.0;
 524       double cbbb = -s->cgate * ((2.0 / 3.0 + Args1) * dVthdVbs
 525                                  + Vgs_Vth * Args1 * dAdVbs / A);
 526       d->cbsb = -(d->cbgb + d->cbdb + cbbb);
 527       trace4("", d->qbulk, d->cbgb, d->cbdb, d->cbsb);
 528
 529       if (m->xpart >= 1) {
 530         /*0/100 partitioning for drain/source chArges at saturation region*/
 531         d->qdrn = 0.0;
 532         d->cdgb = 0.0;
 533         d->cddb = 0.0;
 534         d->cdsb = 0.0;
 535         trace4("sat 0/100", d->qdrn, d->cdgb, d->cddb, d->cdsb);
 536       }else{
 537         /*40/60 partitioning for drain/source chArges at saturation region*/
 538         const double co4v15 = 4./15.;
 539         d->qdrn = -co4v15 * s->cgate * Vgs_Vth;
 540         d->cdgb = -co4v15 * s->cgate;
 541         d->cddb =  0.0;
 542         double cdbb = co4v15 * s->cgate * dVthdVbs;
 543         d->cdsb = -(d->cdgb + d->cddb + cdbb);
 544         trace4("sat 40/60", d->qdrn, d->cdgb, d->cddb, d->cdsb);
 545       }
 546     }
 547     if (d->reversed) {
 548       d->ids *= -1;
 549       d->gmr = d->gmf;
 550       d->gmbr = d->gmbf;
 551       d->gmf = d->gmbf = 0;
 552     }else{
 553       d->gmr = d->gmbr = 0.;
 554     }
 555   }
 556 }
 557 /*--------------------------------------------------------------------------*/
 558 /*--------------------------------------------------------------------------*/