013_bsim4_cap_table : BSIM4 device capacitance table
Requires: SmartSpice & Smartview
Minimum Versions: SMARTSPICE 4.6.5.R
The input deck is a single BSIM4 device .......
Input Files
bsim4.va
											/*****************************************************************/
/* Berkeley BSIM4.3.0 Verilog-A model                            */
/*****************************************************************/
//      UPDATED March 8,19 2004 
//      Contributed By:
//      Geoffrey Coram, Ph.D Senior CAD Engineer Analog Devices, Inc.
`define VOLTAGE_MAXDELTA 0.3
`include "discipline.h"
// The following line must be uncommented if RGATEMOD parameter is 3.
//`define RGATE3
// The following line must be uncommented if RBODYMOD parameter is set to 1.
//`define RBODY
//****** Physical constants ******//
`define EPS0          8.85418e-12
`define KboQ          8.617087e-5 
`define EPSSI         1.03594e-10
`define Charge_q      1.60219e-19
`define Charge        1.6021918e-19
//****** Mathematical constants and constants of limitation ******//
`define PI            3.141592654
`define EXP_THRESHOLD 34.0
`define MIN_EXP       1.713908431e-15
`define MAX_EXP       5.834617425e14
//****** Constants for the model ******//
`define MM 3
`define DELTA   1.0e-9
`define DELTA_1 0.02
`define DELTA_3 0.02
`define DELTA_4 0.02
   
//****** Beginning of the model ******//
module mosfet(drain, gate, source, bulk);
  inout drain, gate, source, bulk;
  electrical drain, gate, source, bulk;     // External nodes
  electrical drainp, sourcep;               // Internal nodes for rdsmod
  electrical gatep, gatem;                  // Internal nodes for rgatemod
  electrical drainb, sourceb, bulkp;        // Internal nodes for rbodymod
  
  //****** Definition of all instance and model parameters *****//
  //****** To customize your model to fit your device, just*****//
  //****** modify those parameters in your modelcard       *****//
  //****** -99.0 is used for parameters which are          *****//
  //****** calculated if you don't determine them          *****//
  
  //****** Minimum conductance ******//
  parameter GMIN  = 1e-12;
  
  //****** Geometrical Parameters ******//
  parameter PS = 0.0;
  parameter PD = 0.0;
  parameter AS = 0.0;
  parameter AD = 0.0;
  
  //****** Overlap Capacitance Parameters ******//
  parameter CGBO = -99.0;                     // Gate-bulk overlap capacitance per length
  parameter CGDO = -99.0;                     // Gate-drain overlap capacitance per width
  parameter CGSO = -99.0;                     // Gate-source overlap capacitance per width
  
  //****** Mosfet type ******//
  parameter TYPE = 1;                         // Define the type of the mosfet (NMOS or PMOS)
  
  //****** Parameters L and W ******//
  parameter L = -1.0;      // Length                            
  parameter W = -1.0;      // Width                            
  
  //****** Model Selectors/Controllers ******//
  parameter MOBMOD   = -99.0;     // Mobility model selector
  parameter RDSMOD   = -99.0;     // Bias-dependent source/drain resistance model selector
  parameter IGCMOD   = 0;         // Gate-to-channel tunneling current model selector
  parameter IGBMOD   = 0;         // Gate-to-substrate tunneling current model selector
  parameter CAPMOD   = 2;         // Capacitance model selector
  parameter RGATEMOD = 2;         // Gate resistance model selector
  parameter RBODYMOD = 0;         // Substrate resistance network model selector
  parameter DIOMOD   = 1;         // Source/drain junction diode IV model selctor
  parameter TEMPMOD  = -99.0;     // Temperature mode selector
  parameter GEOMOD   = 0;         // Geometry-dependent parasitics model selector
  parameter RGEOMOD  = 0;         // Source/drain diffusion resistance and contact model selector
  parameter PERMOD   = 1;         // 
  
  //****** Process Parameters  ******//
  parameter EPSROX   = 3.9;       // Gate dielctric constant relative to vacuum
  parameter TOXE     = -99.0;     // Electrical gate equivalent oxide thickness
  parameter TOXP     = TOXE;      // Physical gate equivalent oxide thickness
  parameter TOXM     = TOXE;      // Tox at which parameters are extracted
  parameter DTOX     = 0.0;       // Defined as TOXE-TOXP
  parameter XJ       = 1.5e-7;    // S/D junction depth
  parameter GAMMA1   = -99.0;     // Body-effect coefficient near the surface
  parameter GAMMA2   = -99.0;     // Body-effect coefficient in the bulk
  parameter NDEP     = -99.0;     // Channel doping concentration at depletion edge for zero body bias
  parameter NSUB     = 6.0e16;    // Substrate doping concentration
  parameter NGATE    = 0.0;       // Poly Si gate doping concentration
  parameter NSD      = 1.0e20;    // S/D doping concentration
  parameter VBX      = -99.0;     // Vbs at which the depletion region width equals XT
  parameter XT       = 1.55e-7;   // Doping depth
  parameter RSH      = 0.0;       // S/D sheet resistance
  parameter RSHG     = 0.0;       // Gate electrode sheet resistance
  
  //****** Basic Parameters  ******//
  parameter VTH0     = -99.0;     // Long-channel threshold voltage at Vbs=0
  parameter VFB      = -99.0;     // Flat-band voltage
  parameter PHIN     = 0.0;       // Non-uniform vertical doping effect on surface potential
  parameter K1       = -99.0;     // First-order body bias coefficient
  parameter K2       = -99.0;     // Second-order body bias coefficient
  parameter K3       = 80.0;      // Narrow width coefficient
  parameter K3B      = 0.0;       // Body effect coefficient of K3
  parameter W0       = 2.5e-6;    // Narrow width parameter
  parameter LPE0     = 1.74e-7;   // Lateral non-uniform doping parameter at Vbs=0
  parameter LPEB     = 0.0;       // Lateral non-uniform doping effect on K1
  parameter VBM      = -3.0;      // Maximum applied body bias in VTH0 calculation
  parameter DVT0     = 2.2;       // First coefficient of short-channel effect on Vth
  parameter DVT1     = 0.53;      // Second coefficient of short-channel effect on Vth
  parameter DVT2     = -0.032;    // Body-bias coefficient of short-channel effect on Vth
  parameter DVTP0    = 0.0;       // First coefficient of drain-induced Vth shift due to for long-channel pocket devices
  parameter DVTP1    = 0.0;       // Second coefficient of drain-induced Vth shift due to for long-channel pocket devices
  parameter DVT0W    = 0.0;       // First coefficient of narrow width effect on Vth for small channel length
  parameter DVT1W    = 5.3e6;     // Second coefficient of narrow width effect on Vth for small channel length
  parameter DVT2W    = -0.032;    // Body-bias coefficient of narrow width effect on Vth for small channel length
  parameter U0       = -99.0;     // Low-field mobility
  parameter UA       = -99.0;     // Coefficient of first-order mobility degradation due to vertical field
  parameter UB       = 1.0e-19;   // Coefficient of second-order mobility degradation due to vertical field
  parameter UC       = -99.0;     // Coefficient of mobility degradation due to body-bias effect
  parameter EU       = -99.0;     // Exponent for mobility degradation of MOBMOD=2
  parameter VSAT     = 8.0e4;     // Saturation velocity
  parameter A0       = 1.0;       // Coefficient of channel-length dependence of bulk charge effect
  parameter AGS      = 0.0;       // Coefficient of Vgs dependence of bulk charge effect
  parameter B0       = 0.0;       // Bulk charge effect coefficient for channel width
  parameter B1       = 0.0;       // Bulk charge effect width offset
  parameter KETA     = -0.047;    // Body-bias coefficient of bulk charge effect
  parameter A1       = 0.0;       // First non-saturation effect parameter
  parameter A2       = 1.0;       // Second non-saturation factor
  parameter WINT     = 0.0;       // Channel-width offset parameter
  parameter LINT     = 0.0;       // Channel-length offset parameter
  parameter DWG      = 0.0;       // Coefficient of gate bias dependence of Weff
  parameter DWB      = 0.0;       // Coefficient of body bias dependence of Weff
  parameter VOFF     = -0.08;     // Offset voltage in subthreshold region for large W and L
  parameter VOFFL    = 0.0;       // Channel-length dependence of VOFF
  parameter MINV     = 0.0;       // Vgsteff fitting parameter for moderate inversion condition
  parameter NFACTOR  = 1.0;       // Subthreshold swing factor
  parameter ETA0     = 0.08;      // DIBL coefficient in subthreshold region
  parameter ETAB     = -0.07;     // Body-bias coefficient for the subthreshold DIBL effect
  parameter DROUT    = 0.56;      // Channel-length dependence of DIBL effect on Rout
  parameter DSUB     = DROUT;     // DIBL coefficient exponent in subthreshold region
  parameter CIT      = 0.0;       // Interface trap capacitance
  parameter CDSC     = 2.4e-4;    // Coupling cpacitance between S/D and channel
  parameter CDSCB    = 0.0;       // Body-bias sensivity of CDSC
  parameter CDSCD    = 0.0;       // Drain-bias sensivity of CDSC
  parameter PCLM     = 1.3;       // Channel-length modulation parameter
  parameter PDIBL1   = 0.39;      // Parameter for DIBL effect on Rout
  parameter PDIBL2   = 0.0086;    // Parameter for DIBL effect on Rout 
  parameter PDIBLB   = 0.0;       // Body-bias coefficient of DIBL effect on Rout
  parameter PSCBE1   = 4.24e8;    // First substrate current induced body-effect parameter
  parameter PSCBE2   = 1.0e-5;    // Second substrate current induced body-effect parameter
  parameter PVAG     = 0.0;       // Gate-bias dependence of Early voltage
  parameter DELTA    = 0.01;      // Parameter for DC Vdseff
  parameter FPROUT   = 0.0;       // Effect of pocket implant on Rout degradation
  parameter PDITS    = 0.0;       // Impact of drain-induced Vth shift on Rout
  parameter PDITSD   = 0.0;       // Vds dependence of drain-induced Vth shift for Rout
  parameter PDITSL   = 0.0;       // Channel-length dependence of drain-induced Vth shift for Rout
  parameter LAMBDA   = -99.0;     // Velocity overshoot coefficient
  parameter VTL      = -99.0;     // Thermal velocity
  parameter LC       = 5.0e-9;    // Velocity back scattering coefficient
  parameter XN       = 3.0;       // Velocity back scattering coefficient
  
  //****** Assymetric and Bias-Dependent Rds Model Parameters  ******//
  parameter RDSW     = 200.0;     // Zero bias LDD resistance per unit width for RDSMOD=0
  parameter RDSWMIN  = 0.0;       // LDD resistance per unit width at high Vgs and zero Vbs for RDSMOD=0
  parameter RDW      = 100.0;     // Zero bias lightly-doped drain resistance Rd per unit width for RDSMOD=1
  parameter RDWMIN   = 0.0;       // Lightly-doped drain resistance Rd per unit width at high Vgs and zero Vbs for RDSMOD=1
  parameter RSW      = 100.0;     // Zero bias lightly-doped source resistance Rd per unit width for RDSMOD=1
  parameter RSWMIN   = 0.0;       // Lightly-doped source resistance Rd per unit width at high Vgs and zero Vbs for RDSMOD=1
  parameter PRWG     = 1.0;       // Gate-bias dependence of LDD resistance
  parameter PRWB     = 0.0;       // Body-bias dependence of LDD resistance
  parameter WR       = 1.0;       // Channel-width dependence parameter of LDD resistance
  parameter NRS      = -99.0;     // Number of source diffusion squares
  parameter NRD      = -99.0;     // Number of drain diffusion squares
  
  //****** Impact Ionization Current Model Parameters  ******//
  parameter ALPHA0   = 0.0;       // First parameter of impact ionization current
  parameter ALPHA1   = 0.0;       // Isub parameter length scaling
  parameter BETA0    = 30.0;      // Second parameter of impact ionization current
  
  //****** Gate Induced Drain Leakage Model Parameters  ******//
  parameter AGIDL    = 0.0;       // Pre-exponential coefficient for GIDL
  parameter BGIDL    = 2.3e9;     // Exponential coefficient for GIDL
  parameter CGIDL    = 0.5;       // Parameter for body-bias effect on GIDL
  parameter EGIDL    = 0.8;       // Fitting parameter for band bending for GIDL
  
  //****** Gate Dielectric Tunneling Current Model Parameters  ******//
  parameter AIGBACC  = 0.43;      // Parameter for Igb in accumulation
  parameter BIGBACC  = 0.054;     // Parameter for Igb in accumulation
  parameter CIGBACC  = 0.075;     // Parameter for Igb in accumulation
  parameter NIGBACC  = 1.0;       // Parameter for Igb in accumulation
  parameter AIGBINV  = 0.35;      // Parameter for Igb in inversion 
  parameter BIGBINV  = 0.03;      // Parameter for Igb in inversion
  parameter CIGBINV  = 0.006;     // Parameter for Igb in inversion
  parameter EIGBINV  = 1.1;       // Parameter for Igb in inversion
  parameter NIGBINV  = 3.0;       // Parameter for Igb in inversion
  parameter AIGC     = -99.0;     // Parameter for Igcs and Igcd
  parameter BIGC     = -99.0;     // Parameter for Igcs and Igcd
  parameter CIGC     = -99.0;     // Parameter for Igcs and Igcd
  parameter AIGSD    = -99.0;     // Parameter for Igs and Igd
  parameter BIGSD    = -99.0;     // Parameter for Igs and Igd
  parameter CIGSD    = -99.0;     // Parameter for Igs and Igd
  parameter DLCIG    = LINT;      // S/D overlap length for Igs and Igd
  parameter NIGC     = 1.0;       // Parameter for Igcs, Igcd, Igs and Igd
  parameter POXEDGE  = 1.0;       // Factor for the gate oxide thickness in S/D overlap regions
  parameter PIGCD    = 1.0;       // Vds dependence of Igcs and Igcd
  parameter NTOX     = 1.0;       // Exponent for the gate oxide ratio
  parameter TOXREF   = 3.0e-9;    // Nominal gate oxide thickness for gate dielectric tunneling current model only
  
  //****** Charge and Capacitance Model Parameters  ******//
  parameter XPART    = 0.0;       // Charge partition parameter
  parameter CGS0     = 0.0;       // Non LDD region source-gate overlap capacitance per unit channel width
  parameter CGD0     = 0.0;       // Non LDD region drain-gate overlap capacitance per unit channel width
  parameter CGB0     = 0.0;       // Gate-bulk overlap capcitance per unit channel length
  parameter CGSL     = 0.0;       // Overlap capacitance between gate and lightly-doped source region
  parameter CGDL     = 0.0;       // Overlap capacitance between gate and lightly-doped source region
  parameter CKAPPAS  = 0.6;       // Coefficient of bias-dependent overlap capacitance for the source side
  parameter CKAPPAD  = CKAPPAS;   // Coefficient of bias-dependent overlap capacitance for the drain side
  parameter CF       = -99.0;     // Fringing field capacitance
  parameter CLC      = 1.0e-7;    // Constant term for the short channel model
  parameter CLE      = 0.6;       // Exponential term for the short channel model
  parameter DLC      = LINT;      // Channel-length offset parameter for CV model
  parameter DWC      = WINT;      // Channel-width offset parameter for CV model
  parameter VFBCV    = -1.0;      // Flat-band voltage parameter
  parameter NOFF     = 1.0;       // CV parameter in Vgstett,CV for weak to strong inversion
  parameter VOFFCV   = 0.0;       // CV parameter in Vgstett,CV for weak to strong inversion
  parameter ACDE     = 1.0;       // Exponential coefficient for charge thickness in CAPMOD=2 for accumulation and depletion regions
  parameter MOIN     = 15.0;      // Coefficient for the gate-bias dependent surface potential
  
  //****** High-Speed/RF Model Parameters  ******//
  parameter XRCRG1   = 12.0;      // Parameter for distributed channel-resistance effect for intrinsic-input resistance
  parameter XRCRG2   = 1.0;       // Parameter to account for the excess channel diffusion resistance
  parameter RBPB     = 50.0;      // Resistance connected between bNodePrime and bNode
  parameter RBPD     = 50.0;      // Resistance connected between bNodePrime and dbNode
  parameter RBPS     = 50.0;      // Resistance connected between bNodePrime and sbNode
  parameter RBDB     = 50.0;      // Resistance connected between dbNodePrime and bNode
  parameter RBSB     = 50.0;      // Resistance connected between sbNodePrime and bNode
  parameter GBMIN    = 1.0e-12;   // Conductance in parallel with each of the five substrate resistances to avoid potential numerical instability
  
  //****** Layout-Dependent Parasistics Model Parameters  ******//
  parameter DMCG     = 0.0;       // Distance from S/D contact center to the gate edge
  parameter DMCI     = DMCG;      // Distance from S/D contact center to the isolation edge in the channel-length direction
  parameter DMDG     = 0.0;       // Same as DMCG but for merged device only
  parameter DMCGT    = 0.0;       // DMCG of test structures
  parameter NF       = 1.0;       // Number of device fingers
  parameter DWJ      = DWC;       // Offset of the S/D junction width
  parameter MIN      = 0.0;       // Whether to minimize the number of drain or source diffusions for even-number fingered device
  parameter XGW      = 0.0;       // Distance from the gate contact to the channel edge
  parameter XGL      = 0.0;       // Offset of the gate length due to variations in patterning
  parameter XL       = 0.0;       // Channel length offset due to mask/etch effect
  parameter XW       = 0.0;       // Channel width offset due to mask/etch effect
  parameter NGCON    = 1.0;       // Number of gate contacts
  
  //****** Assymetric Source/Drain Junction Diode Model Parameters  ******//
  parameter IJTHSREV = 0.1;       // Limiting current in reverse bias region
  parameter IJTHDREV = IJTHSREV;  // Idem
  parameter IJTHSFWD = 0.1;       // Limiting current in forward bias region
  parameter IJTHDFWD = IJTHSFWD;  // Idem
  parameter XJBVS    = 1.0;       // Fitting parameter for diode breakdown
  parameter XJBVD    = XJBVS;     // Idem
  parameter BVS      = 10.0;      // Breakdown voltage
  parameter BVD      = BVS;       // Idem
  parameter JSS      = 1.0e-4;    // Bottom junction reverse saturation current density
  parameter JSD      = JSS;       // Idem
  parameter JSWS     = 0.0;       // Isolation-edge sidewall reverse saturation current density
  parameter JSWD     = JSWS;      // Idem
  parameter JSWGS    = 0.0;       // Gate-edge sidewall reverse saturation current density
  parameter JSWGD    = JSWGS;     // Idem
  parameter CJS      = 5.0e-4;    // Bottom junction capacitance per unit area at zero bias
  parameter CJD      = CJS;       // Idem
  parameter MJS      = 0.5;       // Bottom junction capacitance grating coefficient
  parameter MJD      = MJS;       // Idem
  parameter MJSWS    = 0.33;      // Isolation-edge sidewall junction capacitance grading coefficient
  parameter MJSWD    = MJSWS;     // Idem
  parameter CJSWS    = 5.0e-10;   // Isolation-edge sidewall junction capacitance per unit area
  parameter CJSWD    = CJSWS;     // Idem
  parameter CJSWGS   = CJSWS;     // Gate-edge sidewall junction capacitance per unit length
  parameter CJSWGD   = CJSWS;     // Idem
  parameter MJSWGS   = MJSWS;     // Gate-edge sidewall junction capacitance grading coefficient
  parameter MJSWGD   = MJSWS;     // Idem
  parameter PBS      = 1.0;       // Bottom junction built-in potential
  parameter PBD      = PBS;       // Idem
  parameter PBSWS    = 1.0;       // Isoaltion-edge sidewall junction built-in potential
  parameter PBSWD    = PBSWS;     // Idem
  parameter PBSWGS   = PBSWS;     // Gare-edge sidewall junction built-in potential
  parameter PBSWGD   = PBSWS;     // Idem
  
  //****** Temperature Dependence Parameters  ******//
  parameter TNOM     = 27;        // Temperature at which parameters are extracted
  parameter UTE      = -1.5;      // Mobility temperature exponent
  parameter KT1      = -0.11;     // Tempertature coefficient for threshold voltage
  parameter KT1L     = 0.0;       // Channel length dependence of the temperature coefficient for threshold voltage
  parameter KT2      = 0.022;     // Body-bias coefficient of Vth temperature effect
  parameter UA1      = 1.0e-9;    // Temperature coefficient for UA
  parameter UB1      = -1.0e-18;  // Temperature coefficient for UB
  parameter UC1      = -99.0;     // Temperature coefficient for UC
  parameter AT       = 3.3e4;     // Temperature coefficient for saturation velocity
  parameter PRT      = 0.0;       // Temperature coefficient for Rdsw
  parameter NJS      = 1.0;       // Emission coefficients of junction for drain and source jonctions
  parameter NJD      = NJS;       // Idem
  parameter XTIS     = 3.0;       // Junction current temperature exponents for source and drain junctions
  parameter XTID     = XTIS;      // Idem
  parameter TPB      = 0.0;       // Temperature coefficient of PB
  parameter TPBSW    = 0.0;       // Temperature coefficient of PBSW
  parameter TPBSWG   = 0.0;       // Temperature coefficient of PBSWG
  parameter TCJ      = 0.0;       // Temperature coefficient of CJ
  parameter TCJSW    = 0.0;       // Temperature coefficient of CJSW
  parameter TCJSWG   = 0.0;       // Temperature coefficient of CJSWG
  
  //****** Stress Effect Model Parameters  ******//
  parameter SA       = 0.0;       // Distance between OD edge to poly from one side
  parameter SB       = 0.0;       // Distance between OD edge to poly from other side
  parameter SD       = 0.0;       // Distance between neighbouring fingers
  parameter SAREF    = 1e-6;      // Reference distance between OD and edge to poly of one side
  parameter SBREF    = 1e-6;      // Reference distance between OD and edge to poly of the other side
  parameter WLOD     = 0.0;       // Width parameter for stress effect
  parameter KU0      = 0.0;       // Mobility degradation/enhancement coefficient for stress effect
  parameter KVSAT    = 0.0;       // Saturation velocity degradation/enhancement parameter for stress effect
  parameter TKU0     = 0.0;       // Temperature coefficient of KU0
  parameter LKU0     = 0.0;       // Length dependence of KU0
  parameter WKU0     = 0.0;       // Width dependence of KU0
  parameter PKU0     = 0.0;       // 
  parameter LLODKU0  = 0.0;       // Length parameter for U0 stress effect
  parameter WLODKU0  = 0.0;       // Width parameter for U0 stress effect
  parameter KVTH0    = 0.0;       // Threshold shift parameter for stress effect
  parameter LKVTH0   = 0.0;       // Length dependence of KVTH0
  parameter WKVTH0   = 0.0;       // Width dependence of KVTH0
  parameter PKVTH0   = 0.0;       // Cross-term dependence of KVTH0
  parameter LLODVTH  = 0.0;       // Length parameter for Vth stress effect
  parameter WLODVTH  = 0.0;       // Width parameter for Vth stress effect
  parameter STK2     = 0.0;       // K2 shift factor related to VTH0 change
  parameter LODK2    = 1.0;       // K2 shift modification factor for stress effect
  parameter STETA0   = 0.0;       // ETA0 shift factor related to VTH0 change
  parameter LODETA0  = 1.0;       // ETA0 shift modification factor for stress effect
  
  //****** DW and DL Parameters  ******//
  parameter WL       = 0.0;       // Coefficient of length dependence for width offset
  parameter WLN      = 1.0;       // Power of length dependence of width offset
  parameter WW       = 0.0;       // Coefficient of width dependence for width offset
  parameter WWN      = 1.0;       // Power of width dependence of width offset
  parameter WWL      = 0.0;       // Coefficient of length and width cross-term dependence for width 
  parameter LL       = 0.0;       // Coefficient of length dependence for length offset 
  parameter LLN      = 1.0;       // Power of length dependence of length offset
  parameter LW       = 0.0;       // Coefficient of width dependence for length offset
  parameter LWN      = 1.0;       // Power of width dependence of length offset
  parameter LWL      = 0.0;       // Coefficient of length and width cross-term dependence for length 
  parameter LLC      = LL;        // Coefficient of length dependence for CV channel length offset
  parameter LWC      = LW;        // Coefficient of width dependence for CV channel length offset
  parameter LWLC     = LWL;       // Power of length and width cross-term dependence for CV channel length offset
  parameter WLC      = WL;        // Coefficient of length dependence for CV channel width offset
  parameter WWC      = WW;        // Coefficient of width dependence for CV channel width offset
  parameter WWLC     = WWL;       // Power of length and width cross-term dependence for CV channel width offset
  
  integer i;
  real type, mode, gmin;
  real t0, t1, t2, t3, t4, t5, t6, t7, t8, t9, t10, t11, t12, t13;
  real tmp, tmp1, tmp2, tmp3, tmp4;
  real mobmod, rdsmod, tempmod, rgatemod, permod, geomod, rgeomod;
  real igcmod, igbmod, diomod, capmod, rbodymod;
  real toxe, toxp, coxe, coxp, epsrox;
  real tnom, vtm0, vtm, eg0, eg, ni;
  real l, w, lnew, wnew, nf, xl, xw;
  real ll, lw, lwl, lln, lwn, lint, dl, leff;
  real wl, ww, wwl, wln, wwn, wint, dw, weff;
  real gamma1, gamma2, ndep, nsub;
  real phi, phin, sqrtphi, k3, vfbzb, w0, phis, sqrtphis;
  real xj, nsd, xdep0, litl, vbi, vfbsd, ngate, xdep;
  real cdep0, ntox, toxratio, toxratioedge, toxref, poxedge;
  real mstar, voff, voffl, voffcbn, minv, ldeb;
  real k1, k2, vbx, vbm, xt, vbsc, vfb, vth0, k1ox, k2ox, toxm;
  real vtfbphi1, vtfbphi2, thetarout, dsub, theta0vb0;
  real drout, pdibl1, pdibl2, factor1, dvt0w, dvt1w;
  real dvt1, dvt0, lpe0, lpeb, kt1, kt1l, tratio;
  real vgs, vds, vbs, Vds, Vbs, vses, vdes, vbd, vgd;
  real vsbs, vdbd, vgb, vded, vgeg, vgmg, vgegm, vgmb;
  real vbsb, vbdb, vbeb, vbes, vges, vgms, vdbs, vbesb, vbedb;
  real vbseff, v0, dvt2, dvt2w, lt1, ltw, theta0, thetavth;
  real delt_vth, kt2, vth_narroww, eta0, etab, ddibl_sft_dvd;
  real dibl_sft, lpe_vb, vth, k3b, dvtp0, dvtp1;
  real cit, nfactor, cdsc, cdscb, cdscd, n;
  real vgs_eff, vgd_eff, Vgs_eff, vgst, expvgst, vgsteff;
  real Weff, dwg, dwb, weffcj, powweffwr, wr, ua, ua1, ub, ub1, uc, uc1;
  real vsat, vsattemp, at, prt, rdw, rsw, rdwmin, rswmin;
  real rdsw, rdswmin, deltemp, u0, u0temp, ute, eu;
  real wlc, wwlc, wwc, dwj, rds0, rd0, rs0, rds;
  real prwb, prwg, a0, ags, b0, b1, ueff;
  real abulk, abulk0, dabulk_dvg, keta, denomi;
  real wvcox, wvcoxrds, esat, esatl, a1, a2, Lambda;
  real vgst2vtm, vdsat, delta, vdseff, diffvds;
  real lambda, vasat, tcen, coxeff, coxeffwovl, beta;
  real abovvgst2vtm, fgche1, fgche2, gche, idl;
  real fprout, fp, pvag, pvagterm, pclm, cclm, vaclm;
  real vadibl, pdiblb, pdits, pditsl, pditsd, vadits;
  real pscbe1, pscbe2, vascbe, idsa, ids;
  real alpha0, alpha1, beta0, isub, cdrain;
  real vtl, vs, lc, xn, tfactor, fsevl;
  real grgeltd, xgl, rshg, ngcon, xgw, xrcrg1, xrcrg2, gcrg;
  real llodku0, wlod, w_tmp, wlodku0, lku0, wku0, pku0, tku0;
  real llodvth, wlodvth, ku0, lkvth0, wkvth0, pkvth0, kvth0;
  real ku0temp, ldrn, inv_saref, inv_sbref, inv_od_ref, rho_ref;
  real sa, sb, saref, sbref, sd, lodk2, lodeta0, kvsat;
  real inv_sa, inv_sb, inv_odeff, rho, od_offset, stk2;
  real dvth0_lod, dk2_lod, steta0, deta0_lod;
  real pseff, pdeff, adeff, aseff, dmcg, dmcgt, dmcgeff;
  real dmcieff, dmci, dmdg, dmdgeff;
  real ps, pd, as, ad, imin;
  real nrs, nrd, rsh, gsdiff, gddiff, gstot, gdtot, rs, rd;
  real agidl, bgidl, cgidl, egidl, igidl, igisl;
  real v3, v4, vfbeff, voxacc, voxdepinv, vxnvt, expvxnvt, vaux;
  real dlc, dlcig, aechvb, bechvb, aechvbedge, bechvbedge;
  real aigc, bigc, cigc, aigsd, bigsd, cigsd, aigbacc, bigbacc;
  real cigbacc, nigbacc, aigbinv, bigbinv, cigbinv, nigbinv;
  real nigc, pigcd, eigbinv, igc, modif_pigcd, igcs, igcd;
  real igs, igd, igb, igbacc, igbinv;
  real llc, lwc, lwlc, pbs, pbsws, pbswgs, pbd, pbswd, pbswgd;
  real cgbo, param_cgdo, param_cgso, cgdo, cgso;
  real xtis, xtid, jss, jsd, jsws, jswd, jswgs, jswgd;
  real jss_temp, jsd_temp, jsws_temp, jswd_temp;
  real jswgs_temp, jswgd_temp, njs, njd, cgdl, cgsl;
  real dwc, tcj, tcjsw, tcjswg, cjs, cjd, cjsws, cjswd;
  real cjswgs, cjswgd, cjs_temp, cjd_temp, cjsws_temp, cjswd_temp;
  real cjswgs_temp, cjswgd_temp, tpb, tpbsw, tpbswg, phibs, phibd;
  real phibsws, phibswd, phibswgs, phibswgd;
  real ijthsfwd, ijthsrev, ijthdfwd, ijthdrev, xjbvd, bvd;
  real xjbvs, bvs, weffcv, leffcv, cf, acde;
  real gbmin, rbdb, grbdb, rbpb, grbpb, rbsb, grbsb, rbpd, grbpd;
  real rbps, grbps, nvtms, nvtmd, isbs, isbd, xexpbvs, xexpbvd;
  real vjsmfwd, vjdmfwd, ivjsmfwd, ivjdmfwd, sslpfwd, dslpfwd; 
  real vjsmrev, vjdmrev, ivjsmrev, ivjdmrev, sslprev, dslprev; 
  real cbs, cbd, evbs, evbd, vbs_jct, vbd_jct;
  real csub, clc, cle, abulkcvfactor, xpart, vfbcv;
  real coxwl, arg1, qdrn, qbulk, qgate, qsrc, ccn;
  real abulkcv, alphaz, vdsatcv, two_third_coxwl;
  real vbseffcv, noff, voffcv, vgstnvt, qac0, qsub0;
  real vdseffcv, tox, ccen, link, coxwlcen, moin, deltaphi;
  real arg, czbd, czbs, czbdsw, czbssw, czbdswg, czbsswg, sarg;
  real mjs, mjd, mjsws, mjswd, mjswgs, mjswgd;
  real ckappas, ckappad, qgso, qgdo, dt0_dvg;
  real qgmb, qgmid, qgb, qbd, qbs, qb, qd, qs;
  //****** Function to calculate geometrical parameters ******//
  //****** and intermediaries for rd an rs calculation  ******//
  analog function real get_nuintd;
    input nf, minsd;
    real  nf, minsd;
  
    begin
     if ((nf%2) != 0)
       get_nuintd = 2.0 * max((nf - 1.0) / 2.0, 0.0);
    else
      begin
        if (minsd == 1) /* minimize # of source */
          get_nuintd = 2.0 * max((nf / 2.0 - 1.0), 0.0);
        else
          get_nuintd = nf;
      end
    end
  endfunction
  
  analog function real get_nuendd;
    input nf, minsd;
    real  nf, minsd;
  
    begin
      if ((nf%2) != 0)
        get_nuendd = 1.0;
      else
        begin
          if (minsd == 1) /* minimize # of source */
            get_nuendd = 2.0;
          else
            get_nuendd = 0.0;
        end
    end
  endfunction
  
  analog function real get_nuints;
    input nf, minsd;
    real  nf, minsd;
  
    begin
      if ((nf%2) != 0)
        get_nuints= 2.0 * max((nf - 1.0) / 2.0, 0.0);
      else
        begin
          if (minsd == 1) /* minimize # of source */
            get_nuints = nf;
          else
            get_nuints = 2.0 * max((nf / 2.0 - 1.0), 0.0);
        end
    end
  endfunction
  
  analog function real get_nuends;
    input nf, minsd;
    real  nf, minsd;
  
    begin
      if ((nf%2) != 0)
        get_nuends = 1.0;
      else
        begin
          if (minsd == 1) /* minimize # of source */
            get_nuends = 0.0;
          else
            get_nuends = 2.0;
        end
    end
  endfunction
  
  analog function real get_ps;
    input nf, geo, minsd, weffcj, dmcg, dmci, dmdg;
    real  nf, geo, minsd, weffcj, dmcg, dmci, dmdg;
    real  psiso, pssha, psmer;
    real  t0, t1, t2;
    real  nuints, nuends;
  
    begin
  
      nuints = 0.0;
      nuends = 0.0;
  
      if (geo < 9) 
        begin
          nuints = get_nuints(nf, minsd);
          nuends = get_nuends(nf, minsd);
        end 
  
      t0 = dmcg + dmci;
      t1 = dmcg + dmcg;
      t2 = dmdg + dmdg;
    
      psiso = t0 + t0 + weffcj;
      pssha = t1;
      psmer = t2;
    
      case(geo)
        0: get_ps = nuends * psiso + nuints * pssha;
        1: get_ps = nuends * psiso + nuints * pssha;
        2: get_ps = (nuends + nuints) * pssha;
        3: get_ps = (nuends + nuints) * pssha;
        4: get_ps = nuends * psiso + nuints * pssha;
        5: get_ps = (nuends + nuints) * pssha;
        6: get_ps = nuends * psmer + nuints * pssha;
        7: get_ps = nuends * psmer + nuints * pssha;
        8: get_ps = nuends * psmer + nuints * pssha;
        9: get_ps = psiso + (nf - 1.0) * pssha;
        10: get_ps = nf * pssha;
        default: $strobe("Warning: Specified GEO = %d not matched", geo); 
      endcase
    end
  endfunction
  
  analog function real get_pd;
    input nf, geo, minsd, weffcj, dmcg, dmci, dmdg;
    real  nf, geo, minsd, weffcj, dmcg, dmci, dmdg;
    real  pdiso, pdsha, pdmer;
    real  t0, t1, t2;
    real  nuintd, nuendd;
  
    begin
  
      nuintd = 0.0;
      nuendd = 0.0;
  
      if (geo < 9) 
        begin
          nuintd = get_nuintd(nf, minsd);
          nuendd = get_nuendd(nf, minsd);
        end 
  
      t0 = dmcg + dmci;
      t1 = dmcg + dmcg;
      t2 = dmdg + dmdg;
      
      pdiso = t0 + t0 + weffcj;
      pdsha = t1;
      pdmer = t2;
      
      case(geo)
        0: get_pd = nuendd * pdiso + nuintd * pdsha;
        1: get_pd = (nuendd + nuintd) * pdsha;
        2: get_pd = nuendd * pdiso + nuintd * pdsha;
        3: get_pd = (nuendd + nuintd) * pdsha;
        4: get_pd = nuendd * pdmer + nuintd * pdsha;
        5: get_pd = nuendd * pdmer + nuintd * pdsha;
        6: get_pd = nuendd * pdiso + nuintd * pdsha;
        7: get_pd = (nuendd + nuintd) * pdsha;
        8: get_pd = nuendd * pdmer + nuintd * pdsha;
        9: get_pd = nf * pdsha;
        10: get_pd = pdiso + (nf - 1.0) * pdsha;
        default: $strobe("Warning: Specified GEO = %d not matched", geo); 
      endcase
    end
  endfunction
  
  analog function real get_as;
    input nf, geo, minsd, weffcj, dmcg, dmci, dmdg;
    real  nf, geo, minsd, weffcj, dmcg, dmci, dmdg;
    real  asiso, assha, asmer;
    real  t0;
    real  nuints, nuends;
  
    begin
  
      nuints = 0.0;
      nuends = 0.0;
  
      if (geo < 9) 
        begin
          nuints = get_nuints(nf, minsd);
          nuends = get_nuends(nf, minsd);
        end 
  
      t0 = dmcg + dmci;
    
      asiso = t0 * weffcj;
      assha = dmcg * weffcj;
      asmer = dmdg * weffcj;
    
      case(geo)
        0: get_as = nuends * asiso + nuints * assha;
        1: get_as = nuends * asiso + nuints * assha;
        2: get_as = (nuends + nuints) * assha;
        3: get_as = (nuends + nuints) * assha;
        4: get_as = nuends * asiso + nuints * assha;
        5: get_as = (nuends + nuints) * assha;
        6: get_as = nuends * asmer + nuints * assha;
        7: get_as = nuends * asmer + nuints * assha;
        8: get_as = nuends * asmer + nuints * assha;
        9: get_as = asiso + (nf - 1.0) * assha;
        10: get_as = nf * assha;
        default: $strobe("Warning: Specified GEO = %d not matched", geo); 
      endcase
    end
  endfunction
  
  analog function real get_ad;
    input nf, geo, minsd, weffcj, dmcg, dmci, dmdg;
    real  nf, geo, minsd, weffcj, dmcg, dmci, dmdg;
    real  adiso, adsha, admer;
    real  t0;
    real  nuintd, nuendd;
  
    begin
  
      nuintd = 0.0;
      nuendd = 0.0;
  
      if (geo < 9) 
        begin
          nuintd = get_nuintd(nf, minsd);
          nuendd = get_nuendd(nf, minsd);
        end 
  
      t0 = dmcg + dmci;
    
      adiso = t0 * weffcj;
      adsha = dmcg * weffcj;
      admer = dmdg * weffcj;
    
      case(geo)
        0: get_ad = nuendd * adiso + nuintd * adsha;
        1: get_ad = (nuendd + nuintd) * adsha;
        2: get_ad = nuendd * adiso + nuintd * adsha;
        3: get_ad = (nuendd + nuintd) * adsha;
        4: get_ad = nuendd * admer + nuintd * adsha;
        5: get_ad = nuendd * admer + nuintd * adsha;
        6: get_ad = nuendd * adiso + nuintd * adsha;
        7: get_ad = (nuendd + nuintd) * adsha;
        8: get_ad = nuendd * admer + nuintd * adsha;
        9: get_ad = nf * adsha;
        10: get_ad = adiso + (nf - 1.0) * adsha;
        default: $strobe("Warning: Specified GEO = %d not matched", geo); 
      endcase
    end
  endfunction
  
  analog function real get_rendi;
    input weffcj, rsh, dmcg, dmci, dmdg, nuend, rgeo, type;
    real  weffcj, rsh, dmcg, dmci, dmdg, nuend, rgeo, type;
  
    begin	
      if (type == 1)
        begin
          case(rgeo)
            1, 2, 5:
              begin
                if (nuend == 0.0)
                  get_rendi = 0.0;
                else
                  get_rendi = rsh * dmcg / (weffcj * nuend);
              end
            3, 4, 6:
              begin
                if ((dmcg + dmci) == 0.0)
                  $strobe("(dmcg + dmci) can not be equal to zero");
                if (nuend == 0.0)
                  get_rendi = 0.0;
                else
                  get_rendi = rsh * weffcj / (3.0 * nuend * (dmcg + dmci));
              end
            default:
              $strobe("warning: specified rgeo = %d not matched", rgeo);
          endcase
        end
      else
        begin
          case(rgeo)
            1, 3, 7:
              begin
                if (nuend == 0.0)
                  get_rendi = 0.0;
                else
                  get_rendi = rsh * dmcg / (weffcj * nuend);
              end
            2, 4, 8:
              begin
                if ((dmcg + dmci) == 0.0)
                  $strobe("(dmcg + dmci) can not be equal to zero");
                if (nuend == 0.0)
                  get_rendi = 0.0;
                else
                  get_rendi = rsh * weffcj / (3.0 * nuend * (dmcg + dmci));
              end
            default:
              $strobe("warning: specified rgeo = %d not matched", rgeo);
          endcase
        end
    end
  endfunction
  
  analog function real get_renda;
    input weffcj, rsh, dmcg, dmci, dmdg, nuend, rgeo, type;
    real  weffcj, rsh, dmcg, dmci, dmdg, nuend, rgeo, type;
  
    begin
      if (type == 1)
        begin
          case(rgeo)
            1, 2, 5:
              begin
                if (nuend == 0.0)
                  get_renda = 0.0;
                else
                  get_renda = rsh * dmcg / (weffcj * nuend);
              end
            3, 4, 6:
              begin
                if (dmcg == 0.0)
                  $strobe("dmcg can not be equal to zero");
                if (nuend == 0.0)
                  get_renda = 0.0;
                else
                  get_renda = rsh * weffcj / (6.0 * nuend * dmcg);
              end
            default:
              $strobe("warning: specified rgeo = %d not matched", rgeo);
          endcase
        end
      else
        begin
          case(rgeo)
            1, 3, 7:
              begin
                if (nuend == 0.0)
                  get_renda = 0.0;
                else
                  get_renda = rsh * dmcg / (weffcj * nuend);
              end
            2, 4, 8:
              begin
                if (dmcg == 0.0)
                  $strobe("dmcg can not be equal to zero");
                if (nuend == 0.0)
                  get_renda = 0.0;
                else
                  get_renda = rsh * weffcj / (6.0 * nuend * dmcg);
              end
            default:
              $strobe("warning: specified rgeo = %d not matched", rgeo);
          endcase   
        end
    end
  endfunction
  
  analog function real get_rtot;
    input nf, geo, rgeo, minsd, weffcj, rsh, dmcg, dmci, dmdg, type;
    real  nf, geo, rgeo, minsd, weffcj, rsh, dmcg, dmci, dmdg, type;
    real  rint, rend;
    real  nuintd, nuendd, nuints, nuends;
  
    begin
      
      rend   = 0.0;
      nuintd = 0.0;
      nuendd = 0.0;
      nuints = 0.0;
      nuends = 0.0;
  
      if (geo < 9) /* since geo = 9 and 10 only happen when nf = even */
        begin
    	  nuintd = get_nuintd(nf, minsd);
    	  nuints = get_nuints(nf, minsd);
    	  nuendd = get_nuendd(nf, minsd);
    	  nuends = get_nuends(nf, minsd);
  
          /* internal s/d resistance -- assume shared s or d and all wide contacts */
          if (type == 1)
            begin
              if (nuints == 0.0)
                rint = 0.0;
              else
                rint = rsh * dmcg / ( weffcj * nuints); 
            end
          else
            begin
              if (nuintd == 0.0)
                rint = 0.0;
              else        
                rint = rsh * dmcg / ( weffcj * nuintd);
            end
        end
  
      /* end s/d resistance  -- geo dependent */
      case(geo)
        0:
          begin
            if (type == 1)
              rend = get_rendi(weffcj, rsh, dmcg, dmci, dmdg, nuends, rgeo, 1);
            else
              rend = get_rendi(weffcj, rsh, dmcg, dmci, dmdg, nuendd, rgeo, 0);
          end
        1:
          begin
            if (type == 1)
              rend = get_rendi(weffcj, rsh, dmcg, dmci, dmdg, nuends, rgeo, 1);
            else
              rend = get_renda(weffcj, rsh, dmcg, dmci, dmdg, nuendd, rgeo, 0);
          end
        2:
          begin
            if (type == 1)
              rend = get_renda(weffcj, rsh, dmcg, dmci, dmdg, nuends, rgeo, 1);
            else
              rend = get_rendi(weffcj, rsh, dmcg, dmci, dmdg, nuendd, rgeo, 0);
          end
        3:
          begin
            if (type == 1)
              rend = get_renda(weffcj, rsh, dmcg, dmci, dmdg, nuends, rgeo, 1);
            else
              rend = get_renda(weffcj, rsh, dmcg, dmci, dmdg, nuendd, rgeo, 0);
          end
        4:
          begin
            if (type == 1)
              rend = get_rendi(weffcj, rsh, dmcg, dmci, dmdg, nuends, rgeo, 1);
            else 
              rend = rsh * dmdg / weffcj;
          end
        5:
          begin
            if (type == 1)
              rend = get_renda(weffcj, rsh, dmcg, dmci, dmdg, nuends, rgeo, 1);
            else 
              rend = rsh * dmdg / (weffcj * nuendd);
          end
        6:
          begin
            if (type == 1)
              rend = rsh * dmdg / weffcj;
            else 
              rend = get_rendi(weffcj, rsh, dmcg, dmci, dmdg, nuendd, rgeo, 0);
          end
        7:
          begin
            if (type == 1)
              rend = rsh * dmdg / (weffcj * nuends);
            else 
              rend = get_renda(weffcj, rsh, dmcg, dmci, dmdg, nuendd, rgeo, 0);
          end
        8:
          begin
            rend = rsh * dmdg / weffcj;	
          end
        9:   /* all wide contacts assumed for geo = 9 and 10 */
          begin
            if (type == 1)
              begin
                rend = 0.5 * rsh * dmcg / weffcj;
                if (nf == 2.0)
                  rint = 0.0;
                else
                  rint = rsh * dmcg / (weffcj * (nf - 2.0));
              end
            else
              begin
                rend = 0.0;
                rint = rsh * dmcg / (weffcj * nf);
              end
          end
        10:
          begin
            if (type == 1)
              begin
                rend = 0.0;
                rint = rsh * dmcg / (weffcj * nf);
              end
            else
              begin
                rend = 0.5 * rsh * dmcg / weffcj;;
                if (nf == 2.0)
                  rint = 0.0;
                else
                  rint = rsh * dmcg / (weffcj * (nf - 2.0));
              end
          end
        default:
          $strobe("Warning: specified geo = %d not matched", geo);
      endcase  
  
      if (rint <= 0.0)
        get_rtot = rend;
      else if (rend <= 0.0)
        get_rtot = rint;
      else
        get_rtot = rint * rend / (rint + rend);
      if(get_rtot==0.0)
        $strobe("Warning: zero resistance returned from get_rtot");
    end
  endfunction
  
  analog function real get_vjm;
    input nvtm, ijth, isb, xexpbv;
    real  nvtm, ijth, isb, xexpbv;
    real  tb, tc, evjmovnv;
  
    begin
      tc = xexpbv;
      tb = 1.0 + ijth / isb - tc;
      evjmovnv = 0.5 * (tb + sqrt(tb * tb + 4.0 * tc));
      get_vjm = nvtm * ln(evjmovnv);
    end
  endfunction
  analog
    begin
       
      @(initial_step)
        begin
	  
          a0         = A0;
          a1         = A1;
          a2         = A2;
          acde       = ACDE;
          ad         = AD;
          agidl      = AGIDL;
          ags        = AGS;
          aigbacc    = AIGBACC;
          aigbinv    = AIGBINV;
          aigc       = AIGC;
          aigsd      = AIGSD;
          alpha0     = ALPHA0;
          alpha1     = ALPHA1;
          as         = AS;
          at         = AT;
          b0         = B0;
          b1         = B1;
          beta0      = BETA0;
          bgidl      = BGIDL;
          bigbacc    = BIGBACC;
          bigbinv    = BIGBINV;
          bigc       = BIGC;
          bigsd      = BIGSD;
          bvd        = BVD;
          bvs        = BVS;
          capmod     = CAPMOD;
          cdsc       = CDSC;
          cdscb      = CDSCB;
          cdscd      = CDSCD;
          cf         = CF;
          cgbo       = CGBO;
          cgdl       = CGDL;
          cgidl      = CGIDL;
          cgsl       = CGSL;
          cigbacc    = CIGBACC;
          cigbinv    = CIGBINV;
          cigc       = CIGC;
          cigsd      = CIGSD;
          cit        = CIT;
          cjd        = CJD;
          cjs        = CJS;
          cjswd      = CJSWD;
          cjswgd     = CJSWGD;
          cjswgs     = CJSWGS;
          cjsws      = CJSWS;
          ckappad    = CKAPPAD;
          ckappas    = CKAPPAS;
          clc        = CLC;
          cle        = CLE;
          delta      = DELTA;
          diomod     = DIOMOD;
          dlc        = DLC;
          dlcig      = DLCIG;
          dmcg       = DMCG;
          dmcgt      = DMCGT;
          dmci       = DMCI;
          dmdg       = DMDG;
          drout      = DROUT;
          dsub       = DSUB;
          dvt0       = DVT0;
          dvt0w      = DVT0W;
          dvt1       = DVT1;
          dvt1w      = DVT1W;
          dvt2       = DVT2;
          dvt2w      = DVT2W;
          dvtp0      = DVTP0;
          dvtp1      = DVTP1;
          dwb        = DWB;
          dwc        = DWC;
          dwg        = DWG;
          dwj        = DWJ;
          egidl      = EGIDL;
          eigbinv    = EIGBINV;
          epsrox     = EPSROX;
          eta0       = ETA0;
          etab       = ETAB;
          eu         = EU;
          fprout     = FPROUT;
          gamma1     = GAMMA1;
          gamma2     = GAMMA2;
          gbmin      = GBMIN;
          geomod     = GEOMOD;
          gmin       = GMIN;
          igbmod     = IGBMOD;
          igcmod     = IGCMOD;
          ijthdfwd   = IJTHDFWD;
          ijthdrev   = IJTHDREV;
          ijthsfwd   = IJTHSFWD;
          ijthsrev   = IJTHSREV;
          imin       = MIN;
          jsd        = JSD;
          jss        = JSS;
          jswd       = JSWD;
          jswgd      = JSWGD;
          jswgs      = JSWGS;
          jsws       = JSWS;
          k1         = K1;
          k2         = K2;
          k3         = K3;
          k3b        = K3B;
          keta       = KETA;
          kt1        = KT1;
          kt1l       = KT1L;
          kt2        = KT2;
          ku0        = KU0;
          kvsat      = KVSAT;
          kvth0      = KVTH0;
          l          = L;
          lambda     = LAMBDA;
          lc         = LC;
          lint       = LINT;
          lku0       = LKU0;
          lkvth0     = LKVTH0;
          ll         = LL;
	  llc        = LLC;
	  lln        = LLN;
          llodku0    = LLODKU0;
          llodvth    = LLODVTH;
          lodeta0    = LODETA0;
          lodk2      = LODK2;
          lpe0       = LPE0;
          lpeb       = LPEB;
          lw         = LW;
	  lwc        = LWC;
	  lwl        = LWL;
          lwn        = LWN;
	  lwlc       = LWLC;
	  minv       = MINV;
          mjd        = MJD;
          mjs        = MJS;
          mjswd      = MJSWD;
          mjswgd     = MJSWGD;
          mjswgs     = MJSWGS;
          mjsws      = MJSWS;
          mobmod     = MOBMOD;
          moin       = MOIN;
          ndep       = NDEP;
          nf         = NF;
          nfactor    = NFACTOR;
          ngate      = NGATE;
          ngcon      = NGCON;
          nigbacc    = NIGBACC;
          nigbinv    = NIGBINV;
          nigc       = NIGC;
          njd        = NJD;
          njs        = NJS;
          noff       = NOFF;
          nrd        = NRD;
          nrs        = NRS;
          nsd        = NSD;
          nsub       = NSUB;
          ntox       = NTOX;
          param_cgdo = CGDO;
          param_cgso = CGSO;
          pbd        = PBD;
          pbs        = PBS;
          pbswd      = PBSWD;
          pbswgd     = PBSWGD;
          pbswgs     = PBSWGS;
          pbsws      = PBSWS;
          pclm       = PCLM;
          pd         = PD;
          pdibl1     = PDIBL1;
          pdibl2     = PDIBL2;
          pdiblb     = PDIBLB;
          pdits      = PDITS;
          pditsd     = PDITSD;
          pditsl     = PDITSL;
          permod     = PERMOD;
          phin       = PHIN;
          pigcd      = PIGCD;
          pku0       = PKU0;
          pkvth0     = PKVTH0;
          poxedge    = POXEDGE;
          prt        = PRT;
          prwb       = PRWB;
          prwg       = PRWG;
          ps         = PS;
          pscbe1     = PSCBE1;
          pscbe2     = PSCBE2;
          pvag       = PVAG;
          rbdb       = RBDB;
          rbodymod   = RBODYMOD;
          rbpb       = RBPB;
          rbpd       = RBPD;
          rbps       = RBPS;
          rbsb       = RBSB;
          rdsmod     = RDSMOD;	
          rdsw       = RDSW;
          rdswmin    = RDSWMIN;
          rdw        = RDW;
          rdwmin     = RDWMIN;
          rgatemod   = RGATEMOD;
          rgeomod    = RGEOMOD;
          rsh        = RSH;
          rshg       = RSHG;
          rsw        = RSW;
          rswmin     = RSWMIN;
          sa         = SA;
          saref      = SAREF;
          sb         = SB;
          sbref      = SBREF;
          sd         = SD;
          steta0     = STETA0;
          stk2       = STK2;
          tcj        = TCJ;
          tcjsw      = TCJSW;
          tcjswg     = TCJSWG;
          tempmod    = TEMPMOD;
          tku0       = TKU0;
          toxe       = TOXE;
          toxm       = TOXM;
          toxp       = TOXP;
          toxref     = TOXREF;
          tpb        = TPB;
          tpbsw      = TPBSW;
          tpbswg     = TPBSWG;
          type       = TYPE;
          u0         = U0;
          ua         = UA;
          ua1        = UA1;
          ub         = UB;
          ub1        = UB1;
          uc         = UC;
          uc1        = UC1;
          ute        = UTE;
          vbm        = VBM;
          vbx        = VBX;
          vfb        = VFB;
          vfbcv      = VFBCV;
          voff       = VOFF;
          voffcv     = VOFFCV;
          voffl      = VOFFL;
          vsat       = VSAT;
          vth0       = VTH0;
          vtl        = VTL;
          w          = W;
          w0         = W0;
          wint       = WINT;
          wku0       = WKU0;
          wkvth0     = WKVTH0;
          wl         = WL;
          wlc        = WLC;
          wln        = WLN;
          wlod       = WLOD;
          wlodku0    = WLODKU0;
          wlodvth    = WLODVTH;
          wr         = WR;
          ww         = WW;
          wwc        = WWC;
          wwl        = WWL;
          wwlc       = WWLC;
          wwn        = WWN;
          xgl        = XGL;
          xgw        = XGW;
          xj         = XJ;
          xjbvd      = XJBVD;
          xjbvs      = XJBVS;
          xl         = XL;
          xn         = XN;
	  xpart      = XPART;
	  xrcrg1     = XRCRG1;
          xrcrg2     = XRCRG2;
          xt         = XT;
          xtid       = XTID;
          xtis       = XTIS;
          xw         = XW;
          
          if (mobmod == -99.0) 
            mobmod = 0;
          else if ((mobmod != 0) && (mobmod != 1) && (mobmod != 2))
            begin
              mobmod = 0;
              $strobe("Warning: MOBMOD has been set to its default value: 0.");
            end
          
          if (rdsmod == -99.0)
            rdsmod = 0;
          else if ((rdsmod != 0) && (rdsmod != 1))
            begin
              rdsmod = 0;
          	  $strobe("Warning: RDSMOD has been set to its default value: 0.");
            end
          
          if (tempmod == -99.0)
            tempmod = 0;
          else if ((tempmod != 0) && (tempmod != 1))
            begin
              tempmod = 0;
              $strobe("Warning: TEMPMOD has been set to its default value: 0.");
            end
          
          if ((diomod != 0) && (diomod != 1) && (diomod != 2))
            begin
              diomod = 1;
              $strobe("Warning: DIOMOD has been set to its default value: 1.");
            end
          
          if ((capmod != 0) && (capmod != 1) && (capmod != 2))
            begin
              capmod = 2;
              $strobe("Warning: CAPMOD has been set to its default value: 2.");
            end
          
          if ((permod != 0) && (permod != 1))
            begin
              permod = 1;
              $strobe("Warning: PERMOD has been set to its default value: 1.");
            end
  
          if ((igcmod != 0) && (igcmod != 1))
            begin
              igcmod = 0;
              $strobe("Warning: IGCMOD has been set to its default value: 0.");
            end
          
          if ((igbmod != 0) && (igbmod != 1))
            begin
              igbmod = 0;
              $strobe("Warning: IGBMOD has been set to its default value: 0.");
            end
          
          if ((rbodymod != 0) && (rbodymod != 1))
            begin
              rbodymod = 0;
              $strobe("Warning: RBODYMOD has been set to its default value: 0.");
            end
          
          if (toxref <= 0.0)
            begin
              $strobe("Fatal: TOXREF = %e is not positive.", TOXREF);
              $finish(1);
            end
          
          if (toxe != -99.0 && toxp != -99.0 && DTOX != 0.0 && (toxe != (toxp + DTOX)))
            $strobe("Warning: TOXE, TOXP and DTOX all given and TOXE != TOXP + DTOX. DTOX ignored.");
          else if (toxe != -99.0 && toxp == -99.0)
            toxp = toxe - DTOX;
          else if (toxe == -99.0 && toxp != -99.0)
            toxe = toxp + DTOX;
          else if (toxp == -99.0 && toxe == -99.0)
            begin
              toxe = 3.0e-9;
              toxp = toxe;
            end
          if (toxe < 1.0e-10)
            $strobe("Warning: TOXE = %e is less than 1A. Recommended TOXE >= 5A", toxe);
          if (toxp < 1.0e-10)
            $strobe("Warning: TOXP = %e is less than 1A. Recommended TOXP >= 5A", toxp);
          
          if (toxm == -99.0)
            toxm = toxe;
          if (toxm <= 0.0)
            begin
              $strobe("Fatal: TOXM = %e is not positive.", toxm);
              $finish(1);
            end
            
          if (toxm < 1.0e-10)
            $strobe("Warning: TOXM = %e is less than 1A. Recommended TOXM >= 5A", toxm);
          
          if (epsrox <= 0.0)
            begin
              $strobe("Warning: EPSROX is not positive. Default value taken.");
              epsrox = 3.9;
            end
          
          if (TNOM < -273.15)
            begin
              $strobe("Warning: TNOM is not physically possible. Default value taken.");
              tnom = 300.15;
            end
          else
            tnom = TNOM + 273.15;
          
          if (l <= 0.0 )
            begin
              $strobe("FATAL : L is not positive.");
              $finish(1);
            end
          
          if (w <= 0.0 )
            begin
              $strobe("FATAL : W is not positive.");
              $finish(1);
            end
          
          if (nf < 1.0)
            begin
              $strobe("Warning : NF must be at least equal to 1.Default value taken");
              nf = 1.0;
            end
          
          if (phin < -0.4)
            begin
              $strobe("Fatal: phin = %e is less than -0.4.", PHIN);
              $finish(1);
            end
          else
          
          if (nsub <= 0.0)
            begin
              $strobe("Fatal: NSUB = %e is not positive.", NSUB);
              $finish(1);
            end
          else if (NSUB <= 1.0e14)
            $strobe("Warning: NSUB = %e may be too small.", NSUB);
          else if (NSUB >= 1.0e21)
            $strobe("Warning: NSUB = %e may be too large.", NSUB);
          
          if (xj <= 0.0)
            $strobe("Fatal: XJ = %e is not positive.", XJ);
          
          if (ngate < 0.0)
            begin
              $strobe("Fatal: NGATE = %e is not positive.", NGATE);
              $finish(1);
            end
          if (ngate > 1.0e25)
            begin
              $strobe("Fatal: NGATE = %e is too high.", NGATE);
              $finish(1);
            end
          if ((ngate > 0.0) && (ngate <= 1.0e18))
            $strobe("Warning: NGATE = %e is less than 1.E18cm^-3.", NGATE);
          
          if (poxedge <= 0.0)
            begin
              $strobe("Fatal: POXEDGE = %e is non-positive.", POXEDGE);
              $finish(1);
            end
          
          if (dsub < 0.0)
            begin
              $strobe("Fatal: DSUB = %e is negative.", DSUB);
              $finish(1);
            end
          
          if (drout < 0.0)
            begin
              $strobe("Fatal: DROUT = %e is negative.", DROUT);
              $finish(1);
            end
          
          if (pdibl1 < 0.0)
            $strobe("Warning: PDIBL1 = %e is negative.", PDIBL1);
          
          if (pdibl2 < 0.0)
            $strobe("Warning: PDIBL2 = %e is negative.", PDIBL2);
          
          if (dvt1w < 0.0)
            begin
              $strobe("Fatal: DVT1W = %e is negative.", DVT1W);
              $finish(1);
            end
          
          if (dvt1 < 0.0)
            begin
              $strobe("Fatal: DVT1 = %e is negative.", DVT1);   
              $finish(1);
            end
          
          if (dvt0 < 0.0)
            $strobe("Warning: DVT0 = %e is negative.", DVT0);   
          
          if (lpe0 < -leff)
            begin
              $strobe("Fatal: LPE0 = %e is less than -leff.", LPE0);
              $finish(1);
            end
          
          if (w0 == -weff)
            begin
              $strobe("Fatal: (W0 + Weff) = 0 causing divided-by-zero.");
              $finish(1);
            end
          if (abs(1.0e-6 / (w0 + weff)) > 10.0)
            $strobe("Warning: (W0 + Weff) may be too small.");
          
          if (eta0 < 0.0)
            $strobe("Warning: ETA0 = %e is negative.", ETA0); 
          
          if (lpeb < -leff)
            begin
              $strobe("Fatal: LPEB = %e is less than -leff.", LPEB);
              $finish(1);
            end
          
          if (nfactor < 0.0)
            $strobe("Warning: NFACTOR = %e is negative.", NFACTOR);
          
          if (cdsc < 0.0)
            $strobe("Warning: CDSC = %e is negative.", CDSC);
          
          if (cdscd < 0.0)
            $strobe("Warning: CDSCD = %e is negative.", CDSCD);
          
          if (u0 == -99.0)
            u0 = (type == 1) ? 0.067 : 0.025;
          
          if (ua == -99.0)
            ua = (mobmod == 2) ? 1.0e-15 : 1.0e-9;
          
          if (uc == -99.0)
            uc = (mobmod == 1) ? -0.0465 : -0.0465e-9;   
          
          if (uc1 == -99.0)
            uc1 = (mobmod == 1) ? -0.056 : -0.056e-9;   
          
          if (eu == -99.0)
            eu = (type == 1) ? 1.67 : 1.0;;
          
          if (prwg < 0.0)
            begin
              $strobe("Warning: PRWG = %e is negative. Set to zero.", PRWG);
              prwg = 0.0;
            end
          
          if (a2 < 0.01)
            begin
              $strobe("Warning: A2 = %e is too small. Set to 0.01.", a2);
              a2 = 0.01;
            end
          else if (a2 > 1.0)
            begin
              $strobe("Warning: A2 = %e is larger than 1. A2 is set to 1 and A1 is set to 0.", a2);
              a2 = 1.0;
              a1 = 0.0;
            end
          
          if (delta < 0.0)
            begin
              $strobe("Fatal: DELTA = %e is less than zero.", delta);
              $finish(1);
            end
          
          if ((lambda != -99.0) && (lambda > 0.0)) 
            begin  
              if (lambda > 1.0e-9)
                $strobe("Warning: LAMBDA = %e may be too large.", LAMBDA);
            end
          
          if (fprout < 0.0)
            begin
              $strobe("Fatal: FPROUT = %e is negative.", FPROUT);
              $finish(1);
            end
          
          if (pclm <= 0.0)
            begin
              $strobe("Fatal: PCLM = %e is not positive.", PCLM);
              $finish(1);
            end
          
          if (pdits < 0.0)
            begin
              $strobe("Fatal: PDITS = %e is negative.", pdits);
              $finish(1);
            end
          
          if (pditsl < 0.0)
            begin
              $strobe("Fatal: PDITSL = %e is negative.", pditsl);
              $finish(1);
            end
          
          if (pscbe2 <= 0.0)
              $strobe("Warning: PSCBE2 = %e is not positive.", PSCBE2);
          
          if ((vtl != -99.0) && (vtl > 0.0))
            begin  
              if (vtl < 6.0e4)
                $strobe("Warning: Thermal velocity VTL = %e may be too small.", vtl);
            end
          
          if (xn < 3.0)
            begin
              $strobe("Warning: back scattering coeff XN = %e is too small. Reset to 3.0", xn);
              xn = 3.0;
            end
          
          if (rgatemod == 1)
            begin
              if (rshg <= 0.0)
                $strobe("Warning: RSHG should be positive for RGATEMOD = 1.");
            end    
          else if (rgatemod == 2)
            begin
              if (rshg <= 0.0)
                $strobe("Warning: RSHG <= 0.0 for rgateMod = 2.");
              else if (xrcrg1 <= 0.0)
                $strobe("Warning: XRCRG1 <= 0.0 for rgateMod = 2.");
            end
          if (rgatemod == 3)
            begin
              if (rshg <= 0.0)
                $strobe("Warning: RSHG should be positive for RGATEMOD = 3.");
              else if (xrcrg1 <= 0.0)
                $strobe("Warning: XRCRG1 should be positive for RGATEMOD = 3.");
            end  
          
          if (ngcon < 1.0)
            begin
              $strobe("Fatal: The parameter NGCON cannot be smaller than one.");
              $finish(1);
            end
          
          if ((l + xl) <= xgl)
            begin
              $strobe("Fatal: The parameter XGL must be smaller than Ldrawn+XL.");
              $finish(1);
            end
          
          if((sa > 0.0) && (sb > 0.0) && ((nf == 1.0) || ((nf > 1.0) && (sd > 0.0))) )
            begin
              if (saref <= 0.0)
              begin
                $strobe("Fatal: SAREF = %e is not positive.",saref);
                $finish(1);
              end
              if (sbref <= 0.0)
                begin
                  $strobe("Fatal: SBREF = %e is not positive.",sbref);
                  $finish(1); 
                end
              if (wlod < 0.0)
                begin
                  $strobe("Warning: WLOD = %e is less than 0.",wlod);
                  wlod = 0.0;
                end
              if (kvsat < -1.0 )
                begin
                  $strobe("Warning: KVSAT = %e is is too small; Reset to -1.0.",kvsat);
                  kvsat = -1.0;
                end
              if (kvsat > 1.0)
                begin
                  $strobe("Warning: KVSAT = %e is too big; Reset to 1.0.",kvsat);
                  kvsat = 1.0;
                end
              if (lodk2 <= 0.0)
                  $strobe("Warning: LODK2 = %e is not positive.",lodk2);
              if (lodeta0 <= 0.0)
                  $strobe("Warning: LODETA0 = %e ih not positive.",lodeta0);
            end
          
          if (aigc == -99.0)
            aigc = (type == 1) ? 0.43 : 0.31;
          if (bigc  == -99.0)
            bigc = (type == 1) ? 0.054 : 0.024;
          if (cigc == -99.0)
            cigc = (type == 1) ? 0.075 : 0.03;
          if (aigsd == -99.0)
            aigsd = (type == 1) ? 0.43 : 0.31;
          if (bigsd == -99.0)
            bigsd = (type == 1) ? 0.054 : 0.024;
          if (cigsd == -99.0)
            cigsd = (type == 1) ? 0.075 : 0.03;
          
          if (nigbinv <= 0.0)
            begin
              $strobe("Fatal: NIGBINV = %e is non-positive.", nigbinv);
              $finish(1);
            end
          
          if (nigbacc <= 0.0)
            begin
              $strobe("Fatal: NIGBACC = %e is non-positive.", nigbacc);
              $finish(1);
            end
          
          if (nigc <= 0.0)
            begin
              $strobe("Fatal: NIGC = %e is non-positive.", nigc);
              $finish(1); 
            end
          
          if (pigcd <= 0.0)
            begin
              $strobe("Fatal: PIGCD = %e is non-positive.", pigcd);
              $finish(1);
            end
          
          if (pbs < 0.1)  
            begin
              pbs = 0.1;
              $strobe("Given PBS is less than 0.1. PBS is set to 0.1.");
            end
          
          if (pbsws < 0.1)
            begin
              pbsws = 0.1;
              $strobe("Given PBSWS is less than 0.1. PBSWS is set to 0.1.");
            end
          
          if (pbswgs < 0.1)
            begin
              pbswgs = 0.1;
              $strobe("Given PBSWGS is less than 0.1. PBSWGS is set to 0.1.");
            end
          
          if (pbd < 0.1) 
            begin
              pbd = 0.1;
              $strobe("Given PBD is less than 0.1. PBD is set to 0.1.");
            end
          
          if (pbswd < 0.1)
            begin
              pbswd = 0.1;
              $strobe("Given PBSWD is less than 0.1. PBSWD is set to 0.1.");
            end
          
          if (pbswgd < 0.1)
            begin
              pbswgd = 0.1;
              $strobe("Given PBSWGD is less than 0.1. PBSWGD is set to 0.1.");
            end
          
          if (ijthdfwd <= 0.0)
            begin
              ijthdfwd = 0.1;
              $strobe("IJTHDFWD reset to %e.", ijthdfwd);
            end 
          
          if (ijthsfwd <= 0.0)
            begin
              ijthsfwd = 0.1;
              $strobe("IJTHSFWD reset to %e.", ijthsfwd);
            end 
          
          if (ijthdrev <= 0.0)
            begin
              ijthdrev = 0.1;
              $strobe("IJTHDREV reset to %e.", ijthdrev);
            end 
          
          if (ijthsrev <= 0.0)
            begin
              ijthsrev = 0.1;
              $strobe("IJTHSREV reset to %e.", ijthsrev);
            end 
          
          if ((xjbvd <= 0.0) && (diomod == 2))
            begin
              xjbvd = 1.0;
              $strobe("XJBVD reset to %e.", xjbvd);
            end 
          else if ((xjbvd < 0.0) && (diomod == 0))
            begin
              xjbvd = 1.0;
              $strobe("XJBVD reset to %e.", xjbvd);
            end 
          
          if (bvd <= 0.0)
            begin
              bvd = 10.0;
              $strobe("BVD reset to %e.\n", bvd);
            end 
          
          if ((xjbvs <= 0.0) && (diomod == 2))
            begin
              xjbvs = 1.0;
              $strobe("XJBVS reset to %e.\n", xjbvs);
            end 
          else if ((xjbvs < 0.0) && (diomod == 0))
            begin
              xjbvs = 1.0;
              $strobe("XJBVS reset to %e.\n", xjbvs);
            end 
          
          if (bvs <= 0.0)
            begin
              bvs = 10.0;
              $strobe("BVS reset to %g.\n", bvs);
            end 
          
          if (gbmin < 1.0e-20)
            $strobe("Warning: GBMIN = %e is too small.", gbmin);
          
          if (clc < 0.0)
            begin
              $strobe("Fatal: CLC = %e is negative.", CLC);
              $finish(1); 
            end
          
          if (noff < 0.1)
            $strobe("Warning: NOFF = %e is too small.", noff);
           
          if (voffcv < -0.5)
            $strobe("Warning: VOFFCV = %e is too small.", voffcv);
          
          if (moin < 5.0)
            $strobe("Warning: MOIN = %e is too small.", moin);
          if (moin > 25.0)
            $strobe("Warning: MOIN = %e is too large.", moin);
          
          if (ckappas < 0.02)
            begin
              $strobe("Warning: CKAPPAS = %e is too small.", ckappas);
              ckappas = 0.02;
            end
          
          if (ckappad < 0.02)
            begin
              $strobe("Warning: CKAPPAD = %e is too small.", ckappad);
              ckappad = 0.02;
            end
          
          //***** Oxide capacitances (line 110-111, file b4temp.c) *****//
          coxe = epsrox * `EPS0 / toxe;
          coxp = epsrox * `EPS0 / toxp;
  
          //***** Overlap capacitances (line 113-129, file b4temp.c) *****//
          if (param_cgdo == -99.0)
            begin
              if (dlc > 0.0)
                param_cgdo = dlc * coxe - cgdl ;
              else
                param_cgdo = 0.6 * xj * coxe;
            end
          
          if (param_cgso == -99.0)
            begin
              if (dlc > 0.0)
                param_cgso = dlc * coxe - cgsl ;
              else
                param_cgso = 0.6 * xj * coxe;
            end
          
          if (cgbo == -99.0)
            cgbo = 2.0 * dwc * coxe;
  
          tratio = $temperature / tnom;
          factor1 = sqrt(`EPSSI / (epsrox * `EPS0) * toxe);
  
          //***** Intrinsic carrier concentration (line 139-141, file b4temp.c) *****//
          vtm0 = `KboQ * tnom;
          eg0  = 1.16 - 7.02e-4 * tnom * tnom / (tnom + 1108.0);
          ni   = 1.45e10 * (tnom / 300.15) * sqrt(tnom / 300.15) * exp(21.5565981 - eg0 / (2.0 * vtm0));
          vtm  = `KboQ * $temperature;
  
          //***** Energy gap (line 145, file b4temp.c) *****//
          eg = 1.16 - 7.02e-4 * $temperature * $temperature / ($temperature + 1108.0);
  
          //***** Temperture dependance of Junction diode IV (line 149-191, file b4temp.c) *****//
          if ($temperature != tnom)
            begin
              t0         = eg0 / vtm0 - eg / vtm;
              t1         = ln($temperature / tnom);
              t2         = t0 + xtis * t1;
              t3         = exp(t2 / njs);
              jss_temp   = jss * t3;
              jsws_temp  = jsws * t3;
              jswgs_temp = jswgs * t3;
              t2         = t0 + xtid * t1;
              t3         = exp(t2 / njd);
              jsd_temp   = jsd * t3;
              jswd_temp  = jswd * t3;
              jswgd_temp = jswgd * t3;
            end
          else
            begin
              jss_temp   = jss;
              jsws_temp  = jsws;
              jswgs_temp = jswgs;
              jsd_temp   = jsd;
              jswd_temp  = jswd;
              jswgd_temp = jswgd;
            end
          
          if (jss_temp < 0.0)
            jss_temp   = 0.0;
          if (jsws_temp < 0.0)
            jsws_temp  = 0.0;
          if (jswgs_temp < 0.0)
            jswgs_temp = 0.0;
          if (jsd_temp < 0.0)
            jsd_temp   = 0.0;
          if (jswd_temp < 0.0)
            jswd_temp  = 0.0;
          if (jswgd_temp < 0.0)
            jswgd_temp = 0.0;
  
          //***** Temperature dependence of D/B and S/B diode capacitance (line 193-278, file b4temp.c) *****//
          deltemp = $temperature - tnom;
          t0 = tcj * deltemp;
  
          if (t0 >= -1.0)
            begin
              cjs_temp = cjs *(1.0 + t0);
              cjd_temp = cjd *(1.0 + t0);
            end
          else
            begin
              if (cjs > 0.0)
                begin
                  cjs_temp = 0.0;
                  $strobe("Temperature effect has caused CJS to be negative. CJS is clamped to zero.");
                end    
              if (cjd > 0.0)
                begin
                  cjd_temp = 0.0;
                  $strobe("Temperature effect has caused CJD to be negative. CJD is clamped to zero.\n");
                end    
            end
  
          t0 = tcjsw * deltemp;
  
          if (t0 >= -1.0)
            begin
              cjsws_temp = cjsws *(1.0 + t0);
              cjswd_temp = cjswd *(1.0 + t0);
            end
          else
            begin
              if (cjsws > 0.0)
                begin
                  cjsws_temp = 0.0;
                  $strobe("Temperature effect has caused CJSWS to be negative. CJSWS is clamped to zero.");
                end 
              if (cjswd > 0.0)
                begin
                  cjswd_temp = 0.0;
                  $strobe("Temperature effect has caused CJSWD to be negative. CJSWD is clamped to zero.");
                end    	
            end
  
          t0 = tcjswg * deltemp;
  
          if (t0 >= -1.0)
            begin
              cjswgs_temp = cjswgs *(1.0 + t0);
              cjswgd_temp = cjswgd *(1.0 + t0);
            end
          else
            begin
              if (cjswgs > 0.0)
                begin
                  cjswgs_temp = 0.0;
                  $strobe("Temperature effect has caused CJSWGS to be negative. CJSWGS is clamped to zero.");
                end
              if (cjswgd > 0.0)
                begin
                  cjswgd_temp = 0.0;
                  $strobe("Temperature effect has caused CJSWGD to be negative. CJSWGD is clamped to zero.");
                end    
            end
  
          phibs = pbs - tpb * deltemp;
  
          if (phibs < 0.01)
            begin
              phibs = 0.01;
              $strobe("Temperature effect has caused PBS to be less than 0.01. PBS is clamped to 0.01.");
            end
          
          phibd = pbd - tpb * deltemp;
  
          if (phibd < 0.01)
            begin
              phibd = 0.01;
              $strobe("Temperature effect has caused PBD to be less than 0.01. PBD is clamped to 0.01.");
            end    
  
          phibsws = pbsws - tpbsw * deltemp;
  
          if (phibsws <= 0.01)
            begin
              phibsws = 0.01;
              $strobe("Temperature effect has caused PBSWS to be less than 0.01. PBSWS is clamped to 0.01.");
            end 
  
          phibswd = pbswd - tpbsw * deltemp;
  
          if (phibswd <= 0.01)
            begin
              phibswd = 0.01;
              $strobe("Temperature effect has caused PBSWD to be less than 0.01. PBSWD is clamped to 0.01.");
            end    
  
          phibswgs = pbswgs - tpbswg * deltemp;
  
          if (phibswgs <= 0.01)
            begin
              phibswgs = 0.01;
              $strobe("Temperature effect has caused PBSWGS to be less than 0.01. PBSWGS is clamped to 0.01.");
            end 
  
          phibswgd = pbswgd - tpbswg * deltemp;
  
          if (phibswgd <= 0.01)
            begin
              phibswgd = 0.01;
              $strobe("Temperature effect has caused PBSWGD to be less than 0.01. PBSWGD is clamped to 0.01.");
            end
  
          //***** Effective length and width (line 362-396, file b4temp.c) *****//
          lnew  = l + xl ;
          wnew  = w / nf + xw;
          t0    = pow(lnew, lln);
          t1    = pow(wnew, lwn);
          tmp1  = ll / t0 + lw / t1 + lwl / (t0 * t1);
          dl    = lint + tmp1;
          tmp2  = llc / t0 + lwc / t1 + lwlc / (t0 * t1);
          dlc   = dlc + tmp2;
          dlcig = dlcig + tmp2;
          t2    = pow(lnew, wln);
          t3    = pow(wnew, wwn);
          tmp1  = wl / t2 + ww / t3 + wwl / (t2 * t3);
          dw    = wint + tmp1;
          tmp2  = wlc / t2 + wwc / t3 + wwlc / (t2 * t3); 
          dwj   = dwj + tmp2;
          leff  = lnew - 2.0 * dl;
          weff  = wnew - 2.0 * dw;
          leffcv = lnew - 2.0 * dlc;
          
          if (leffcv <= 0.0)
            begin
              $strobe("Fatal: Effective channel length for C-V <= 0");
              $finish(1); 
            end
          
          weffcv = wnew - 2.0 * dwc;
  
          if (weffcv <= 0.0)
            begin
              $strobe("Fatal: Effective channel width for C-V <= 0");
              $finish(1); 
            end
  
          if (leff <= 1.0e-9)
            $strobe("Warning: leff = %e <= 1.0e-9. Recommended leff >= 1e-8.", leff);
  
          if (weff <= 1.0e-9)
            $strobe("Warning: weff = %e <= 1.0e-9. Recommended weff >= 1e-7.", weff);
  
          if (leffcv <= 1.0e-9)
            $strobe("warning: leff for CV = %e <= 1.0e-9. recommended leffcv >=1e-8 ", leffcv);
        
          if (weffcv <= 1.0e-9)
            $strobe("warning: weff for CV = %e <= 1.0e-9. recommended weffcv >= 1e-7 ", weffcv);
  
          //***** weffcj (line 429-437, file b4temp.c) *****//
          weffcj = wnew - 2.0 * dwj;
  
          if (weffcj <= 0.0)
            begin
              $strobe("Fatal: Effective channel width for S/D junctions <= 0.");
              $finish(1);
            end
  
          //***** Temperature model (line 955-1026, file b4temp.c) *****//
  
          abulkcvfactor = 1.0 + pow((clc / leffcv), cle);
          t0 = (tratio - 1.0);
          powweffwr = pow(weffcj * 1.0e6, wr) * nf;
          t1 = 0.0;
          t2 = 0.0;
          t3 = 0.0; 
          t4 = 0.0;
          
          if (tempmod == 0)
            begin
              ua = ua + ua1 * t0;
              ub = ub + ub1 * t0;
              uc = uc + uc1 * t0;
              vsattemp = vsat - at * t0;
              t10 = prt * t0;
          
              if(rdsmod == 1) 
                begin
                  /* External Rd(V) */
          	  t1 = rdw + t10;
                  t2 = rdwmin + t10;
          
          	  /* External Rs(V) */
          	  t3 = rsw + t10;
                  t4 = rswmin + t10;
                end  
          
              /* Internal Rds(V) in IV */
              rds0 = (rdsw + t10) * nf / powweffwr;
              rdswmin = (rdswmin + t10) * nf / powweffwr;
            end
          else /* TEMPMOD = 1 */
            begin
       	      ua = ua * (1.0 + ua1 * deltemp);
              ub = ub * (1.0 + ub1 * deltemp);
              uc = uc * (1.0 + uc1 * deltemp);
              vsattemp = vsat * (1.0 - at * deltemp);
              t10 = 1.0 + prt * deltemp;
          
              if(rdsmod == 1) 
                begin
                  /* External Rd(V) */
          	  t1 = rdw * t10;
                  t2 = rdwmin * t10;
          
          	  /* External Rs(V) */
          	  t3 = rsw * t10;
                  t4 = rswmin * t10;
                end  
          
              /* Internal Rds(V) in IV */
              rds0 = rdsw * t10 * nf / powweffwr;
              rdswmin = rdswmin * t10 * nf / powweffwr;
            end
                
          if (t1 < 0.0)
            begin
              t1 = 0.0;
              $strobe("Warning: rdw at current temperature is negative; set to 0.");
            end
          
          if (t2 < 0.0)
            begin
              t2 = 0.0;
              $strobe("Warning: rdwmin at current temperature is negative; set to 0.");
            end
           
          rd0 = t1 / powweffwr;
          rdwmin = t2 / powweffwr;
          
          if (t3 < 0.0)
            begin
              t3 = 0.0;
              $strobe("Warning: rsw at current temperature is negative; set to 0.");
            end
           
          if (t4 < 0.0)
            begin
              t4 = 0.0;
              $strobe("Warning: rswmin at current temperature is negative; set to 0.");
            end
          
          rs0 = t3 / powweffwr;
          rswmin = t4 / powweffwr;
          
          if (u0 > 1.0) 
            u0 = u0 / 1.0e4;
          
          u0temp = u0 * pow(tratio, ute); 
          
          if (u0temp <= 0.0)
            begin
              $strobe("Fatal: U0 at current temperature = %e is not positive.", u0temp);
              $finish(1);
            end
          
          if (eu < 0.0)
            begin
              eu = 0.0;
              $strobe("Warning: EU has been negative; reset to 0.0.");
            end
          
          if (vsattemp <= 0.0)
            begin
              $strobe("Fatal: VSAT at current temperature = %e is not positive.", vsattemp);
              $finish(1);
            end
          
          if (vsattemp < 1.0e3)
            $strobe("Warning: VSAT at current temperature = %e may be too small.", vsattemp);
          
          if (rds0 < 0.0)
            begin
              $strobe("Warning: Rds at current temperature = %e is negative. Set to zero.", rds0);
              rds0 = 0.0;
            end
          
          if (rdsw < 0.0)
            begin
              $strobe("Warning: rdsw = %e is negative. Set to zero.", rdsw);
              rdsw = 0.0;
              rds0 = 0.0;
            end
          
          if (rdswmin < 0.0)
            begin
              $strobe("Warning: rdswmin at current temperature = %e is negative. Set to zero.", rdswmin);
              rdswmin = 0.0;
            end
  
          if (b1 == -weff)
            begin
              $strobe("Fatal: (B1 + weff) = 0 causing divided-by-zero.");
              $finish(1);
            end
          
          if (abs(1.0e-6 / (b1 + weff)) > 10.0)
             $strobe("Warning: (B1 + weff) may be too small.");
  
          //***** Source End Velocity Limit (line 1024-1034, file b4temp.c) *****//
          if ((vtl != -99.0) && (vtl > 0.0))
            begin  
              if (lc < 0.0) 
                lc = 0.0;
                
              t0 = leff / (xn * leff + lc);
              tfactor = (1.0 - t0) / (1.0 + t0 );
            end
  
          //***** Overlap capacitances (line 1034-1040, file b4temp.c) *****//
          if (cf == -99.0)
            cf = 2.0 * epsrox * `EPS0 / `PI * ln(1.0 + 0.4e-6 / toxe);
  
          param_cgdo = (param_cgdo + cf) * weffcv;
          param_cgso = (param_cgso + cf) * weffcv;
          cgbo = cgbo * leffcv * nf;
  
          //***** Test on ndep and  gamma1 (line 1043-1045, file b4temp.c) *****//
          if (ndep == -99.0 && gamma1 != -99.0)
            begin
              t0     = gamma1 * coxe;
              ndep   = 3.01248e22 * t0 * t0;
            end
          else if (ndep != -99.0 && gamma1 == -99.0)
            gamma1   = 5.753e-12 * sqrt(ndep) / coxe;
          else if (ndep == -99.0 && gamma1 == -99.0)
            begin
              ndep   = 1.0e17;
              gamma1 = 5.753e-12 * sqrt(ndep) / coxe;
            end
           
       	  if (ndep <= 0.0)
            begin
              $strobe("Fatal: NDEP = %e is not positive.", ndep);
              $finish(1); 
            end
         
          if (ndep <= 1.0e12)
            $strobe("Warning: NDEP = %e may be too small.", ndep);
          else if (ndep >= 1.0e21)
            $strobe("Warning: NDEP = %e may be too large.", ndep);
  
          if (gamma2 == -99.0)
            gamma2   = 5.753e-12 * sqrt(nsub) / coxe;
          
          //***** Potential surface (line 1048-1049, file b4temp.c) *****//
          phi     = vtm0 * ln(ndep / ni) + phin + 0.4;
          sqrtphi = sqrt(phi);
  
          //***** Calculation of some intermediaries (line 1054-1091, file b4temp.c) *****//
          xdep0        = sqrt(2.0 * `EPSSI / (`Charge_q * ndep * 1.0e6)) * sqrtphi; 
          litl         = sqrt(3.0 * xj * toxe);
          vbi          = vtm0 * ln(nsd * ndep / (ni * ni));
  
          if (ngate > 0.0)
            vfbsd      = vtm0 * ln(ngate / nsd);
          else
            vfbsd      = 0.0;
  
          cdep0        = sqrt(`Charge_q * `EPSSI * ndep * 1.0e6 / 2.0 / phi);
          toxratio     = exp(ntox * ln(toxref / toxe)) / toxe / toxe;
          toxratioedge = exp(ntox * ln(toxref / (toxe * poxedge))) / toxe / toxe / poxedge / poxedge;
          mstar        = 0.5 + atan(minv) / `PI;
          voffcbn      =  voff + voffl / leff;
          ldeb         = sqrt(`EPSSI * vtm0 / (`Charge_q * ndep * 1.0e6)) / 3.0;
          acde         = acde * pow((ndep / 2.0e16), -0.25);
  
          if (capmod ==2)
            begin
              if (acde < 0.1)
               	$strobe("Warning: ACDE = %e is too small.", acde);
              if (acde > 1.6)
                $strobe("Warning: ACDE = %e is too large.", acde);
            end	
  
          //***** Calculation of K1 and K2 (line 1101-1149, file b4temp.c) *****//
          if (k1 != -99.0 || k2 != -99.0)
            begin
              if (k1 == -99.0)
                begin
                  $strobe("Warning: K1 should be specified with K2.");
                  k1 = 0.53;
                end  
              if (k2 == -99.0)
                begin
                  $strobe("Warning: K2 should be specified with K1.");
                  k2 = -0.0186;
                end
            end
          else
            begin
              if (vbx == -99.0)
                vbx = phi - 7.7348e-4 * ndep * xt * xt;
              if (vbx > 0.0)
                vbx = -vbx;
              if (vbm > 0.0)
                vbm = -vbm;
  
              t0 = gamma1 - gamma2;
              t1 = sqrt(phi - vbx) - sqrtphi;
              t2 = sqrt(phi * (phi - vbm)) - phi;
              k2 = t0 * t1 / (2.0 * t2 + vbm);
              k1 = gamma2 - 2.0 * k2 * sqrt(phi - vbm);
            end
  
          //***** Calculation of vbsc (line 1151-1161, file b4temp.c) *****//
          if (k2 < 0.0)
            begin
              t0 = 0.5 * k1 / k2;
              vbsc = 0.9 * (phi - t0 * t0);
              if (vbsc > -3.0)
                vbsc = -3.0;
              else if (vbsc < -30.0)
                vbsc = -30.0;
            end
          else
            vbsc = -30.0;
  
          if (vbsc > vbm)
            vbsc = vbm;
  
          //***** Flat-band voltage (line 1165-1179, file b4temp.c) *****//
          if (vfb == -99.0)
            begin
              if (vth0 != -99.0)
                vfb = type * vth0 - phi - k1 * sqrtphi;
              else
                vfb = -1.0;
            end
  
          //***** Flat-band voltage (line 1165-1179, file b4temp.c) *****//
           if (vth0 == -99.0)
             vth0 = type * (vfb + phi + k1 * sqrtphi);
  
          //***** Calculation of intermediaries (line 1181-1185, file b4temp.c) *****//
          k1ox = k1 * toxe / toxm;
          k2ox = k2 * toxe / toxm;
  
          //***** Calculation of vfbzb (line 1186-1265, file b4temp.c) *****//
          t3 = type * vth0 - vfb - phi;
          t4 = t3 + t3;
          t5 = 2.5 * t3;
          vtfbphi1 = (type == 1) ? t4 : t5; 
          
          if (vtfbphi1 < 0.0)
            vtfbphi1 = 0.0;
          
          vtfbphi2 = 4.0 * t3;
          
          if (vtfbphi2 < 0.0)
            vtfbphi2 = 0.0;
          
          tmp = sqrt(`EPSSI / (epsrox * `EPS0) * toxe * xdep0);
          t0 = dsub * leff / tmp;
          
          if (t0 < `EXP_THRESHOLD)
            begin
              t1 = exp(t0);
              t2 = t1 - 1.0;
              t3 = t2 * t2;
              t4 = t3 + 2.0 * t1 * `MIN_EXP;
              theta0vb0 = t1 / t4;
            end 
          else
            theta0vb0 = 1.0 / (`MAX_EXP - 2.0);
          
          t0 = drout * leff / tmp;
          
          if (t0 < `EXP_THRESHOLD)
            begin
              t1 = exp(t0);
              t2 = t1 - 1.0;
              t3 = t2 * t2;
              t4 = t3 + 2.0 * t1 * `MIN_EXP;
              t5 = t1 / t4;
            end
          else
            t5 = 1.0 / (`MAX_EXP - 2.0); /* 3.0 * `MIN_EXP omitted */
          
          thetarout = pdibl1 * t5 + pdibl2;
          tmp = sqrt(xdep0);
          tmp1 = vbi - phi;
          tmp2 = factor1 * tmp;
          t0 = dvt1w * weff * leff / tmp2;
          
          if (t0 < `EXP_THRESHOLD)
            begin
              t1 = exp(t0);
              t2 = t1 - 1.0;
              t3 = t2 * t2;
              t4 = t3 + 2.0 * t1 * `MIN_EXP;
              t8 = t1 / t4;
            end
          else
            t8 = 1.0 / (`MAX_EXP - 2.0);
          
          t0 = dvt0w * t8;
          t8 = t0 * tmp1;
          t0 = dvt1 * leff / tmp2;
          if (t0 < `EXP_THRESHOLD)
            begin
              t1 = exp(t0);
              t2 = t1 - 1.0;
              t3 = t2 * t2;
              t4 = t3 + 2.0 * t1 * `MIN_EXP;
              t9 = t1 / t4;
            end   
          else
            t9 = 1.0 / (`MAX_EXP - 2.0);
          
          t9 = dvt0 * t9 * tmp1;
          t4 = toxe * phi / (weff + w0);
          t0 = sqrt(1.0 + lpe0 / leff);
          t5 = k1ox * (t0 - 1.0) * sqrtphi + (kt1 + kt1l / leff) * (tratio - 1.0);
          tmp3 = type * vth0 - t8 - t9 + k3 * t4 + t5;
          vfbzb = tmp3 - phi - k1 * sqrtphi; 
  
          //***** Stress Effect (line 1267-1366, file b4temp.c) *****//
          ldrn = l;
          t0 = pow(lnew, llodku0);
          w_tmp = wnew + wlod;
          t1 = pow(w_tmp, wlodku0);
          tmp1 = lku0 / t0 + wku0 / t1 + pku0 / (t0 * t1);
          ku0 = 1.0 + tmp1;
          t0 = pow(lnew, llodvth);
          t1 = pow(w_tmp, wlodvth);
          tmp1 = lkvth0 / t0 + wkvth0 / t1 + pkvth0 / (t0 * t1);
          kvth0 = 1.0 + tmp1;
          kvth0 = sqrt(kvth0*kvth0 + `DELTA);
          t0 = (tratio - 1.0);
          ku0temp = ku0 * (1.0 + tku0 *t0) + `DELTA;
          inv_saref = 1.0/(saref + 0.5 * ldrn);
          inv_sbref = 1.0/(sbref + 0.5 * ldrn);
          inv_od_ref = inv_saref + inv_sbref;
          rho_ref =  KU0 / ku0temp * inv_od_ref;
  
          if ((sa > 0.0) && (sb > 0.0) && ((nf == 1.0) || ((nf > 1.0) && (sd > 0.0))))
            begin
              inv_sa = 0;
              inv_sb = 0;
          
              for(i = 0; i < nf; i = i+1)
                begin
                  t0 = 1.0 / nf / (sa + 0.5*ldrn + i * (sd +ldrn));
                  t1 = 1.0 / nf / (sb + 0.5*ldrn + i * (sd +ldrn));
                  inv_sa = inv_sa + t0;
                  inv_sb = inv_sb + t1;
                end
          
              inv_odeff = inv_sa + inv_sb; 
              rho = KU0 / ku0temp * inv_odeff;
              t0 = (1.0 + rho)/(1.0 + rho_ref);
              u0temp = u0temp * t0;
              t1 = (1.0 + kvsat * rho)/(1.0 + kvsat * rho_ref);
              vsattemp = vsattemp * t1;
              od_offset = inv_odeff - inv_od_ref;
              dvth0_lod = KVTH0 / kvth0 * od_offset;
              dk2_lod = stk2 / pow(kvth0, lodk2) * od_offset;
              deta0_lod = steta0 / pow(kvth0, lodeta0) * od_offset;
              vth0 = vth0 + dvth0_lod;
          
              if (VFB == -99.0 && VTH0 == -99.0)
                vfb = -1.0;
              else  
                vfb = vfb + type * dvth0_lod;
          
              vfbzb = vfbzb + type * dvth0_lod;
              t3 = type * vth0 - vfb - phi;
              t4 = t3 + t3;
              t5 = 2.5 * t3;
              vtfbphi1 = (type == 1) ? t4 : t5;
              
              if (vtfbphi1 < 0.0)
                vtfbphi1 = 0.0;
          
              vtfbphi2 = 4.0 * t3;
              
              if (vtfbphi2 < 0.0)
                vtfbphi2 = 0.0;
          
              k2 = k2 + dk2_lod;
              
              if (k2 < 0.0)
                begin
                  t0 = 0.5 * k1 / k2;
                  vbsc = 0.9 * (phi - t0 * t0);
                  
                  if (vbsc > -3.0)
                    vbsc = -3.0;
                  else if (vbsc < -30.0)
                    vbsc = -30.0;
                end
              else
                vbsc = -30.0;
              
              if (vbsc > vbm)
                vbsc = vbm;
          
              k2ox = k2 * toxe / toxm;
              eta0 = eta0 + deta0_lod;
            end
        
          //*********** HF model parameters (line 1371-1392, file b4temp.c) *****************//
          if (rbodymod == 1.0)
            begin
              if (rbdb < 1.0e-3)
                grbdb = 1.0e3;
              else
                grbdb = gbmin + 1.0 / rbdb;
                
              if (rbpb < 1.0e-3)
                grbpb = 1.0e3;
              else
                grbpb = gbmin + 1.0 / rbpb;
                
              if (rbps < 1.0e-3)
                grbps = 1.0e3;
              else
                grbps = gbmin + 1.0 / rbps;
                
              if (rbsb < 1.0e-3)
                grbsb = 1.0e3;
              else
                grbsb = gbmin + 1.0 / rbsb;
                
              if (rbpd < 1.0e-3)
                grbpd = 1.0e3;
              else
                grbpd = gbmin + 1.0 / rbpd;
            end
  
  
          //*********** Process geomertry dependent parasitics (line 1396-1452, file b4temp.c) *****************//
          grgeltd = rshg * (xgw + weffcj / 3.0 / ngcon) / (ngcon * nf * (lnew - xgl));
  
          if (grgeltd > 0.0)
            grgeltd = 1.0 / grgeltd;
          else
            begin
              grgeltd = 1.0e3;
              
              if (rgatemod != 0)
                $strobe("Warning: The gate conductance reset to 1.0e3 ohms.");
            end
  
          dmcgeff = dmcg - dmcgt;
          dmcieff = dmci;
          dmdgeff = dmdg - dmcgt;
            
          if (ps > 0.0)
            begin
              if (permod == 0)
                pseff = ps;
              else
                pseff = ps - weffcj * nf;
            end
          else
            pseff = get_ps(nf, geomod, imin, weffcj, dmcgeff, dmcieff, dmdgeff);
          
          if (pd > 0.0)
            begin
              if (permod == 0)
                pdeff = pd;
              else
                pdeff = pd - weffcj * nf;
            end
          else
            pdeff = get_pd(nf, geomod, imin, weffcj, dmcgeff, dmcieff, dmdgeff);
          
          if (as > 0.0)
            aseff = as;
          else
            aseff = get_as(nf, geomod, imin, weffcj, dmcgeff, dmcieff, dmdgeff);
          
          if (ad > 0.0)
            adeff = ad;
          else
            adeff = get_ad(nf, geomod, imin, weffcj, dmcgeff, dmcieff, dmdgeff);
  
          //*********** Processing S/D resistance and conductance below (line 1453-1516, file b4temp.c) *****************//
          if (nrs != -99.0)
            gsdiff = rsh * nrs;
          else if (rgeomod > 0)
            gsdiff = get_rtot(nf, geomod, rgeomod, imin, weffcj, rsh, dmcgeff, dmcieff, dmdgeff, 1);
          else
            gsdiff = 0.0;
          
          if (gsdiff > 0.0)
            gsdiff = 1.0 / gsdiff;
          else
            begin
              gsdiff = 1.0e3; /* mho */
              $strobe ("Warning: source conductance reset to 1.0e3 mho.");
            end
          
          if (nrd != -99.0)
            gddiff = rsh * nrd;
          else if (rgeomod > 0)
            gddiff = get_rtot(nf, geomod, rgeomod, imin, weffcj, rsh, dmcgeff, dmcieff, dmdgeff, 0);
          else
            gddiff = 0.0;
          
          if (gddiff > 0.0)
            gddiff = 1.0 / gddiff;
          else
            begin                 
              gddiff = 1.0e3; /* mho */
              $strobe ("Warning: drain conductance reset to 1.0e3 mho.");
            end
  
          aechvb = (type == 1) ? 4.97232e-7 : 3.42537e-7;
          bechvb = (type == 1) ? 7.45669e11 : 1.16645e12;
          aechvbedge = aechvb * weff * dlcig * toxratioedge;
          bechvbedge = -bechvb * toxe * poxedge;
          aechvb = aechvb * weff * leff * toxratio;
          bechvb = bechvb * -toxe;
  
          //*********** Diode model intermediaries calculation (line 1519-1635, file b4temp.c) *****************//
          nvtms = vtm * njs;
          
          if ((aseff <= 0.0) && (pseff <= 0.0))
            isbs = 1.0e-14;
          else
            isbs = aseff * jss_temp + pseff * jsws_temp + weffcj * nf * jswgs_temp;
          
          if (isbs > 0.0)
            begin
              case(diomod)
                0:
                  begin
                    if ((bvs / nvtms) > `EXP_THRESHOLD)
                      xexpbvs = xjbvs * `MIN_EXP;
                    else
                      xexpbvs = xjbvs * exp(-bvs / nvtms);	
                  end
                1:
                  begin
                    vjsmfwd = get_vjm(nvtms, ijthsfwd, isbs, 0.0);
                    ivjsmfwd = isbs * exp(vjsmfwd / nvtms);
                  end
                2:
                  begin
                    if ((bvs / nvtms) > `EXP_THRESHOLD)
                      begin
                        xexpbvs = xjbvs * `MIN_EXP;
                        tmp = `MIN_EXP;
                      end
                    else
                      begin
                        xexpbvs = exp(-bvs / nvtms);
                        tmp = xexpbvs;
                        xexpbvs = xexpbvs * xjbvs;	
                      end
          
                    vjsmfwd = get_vjm(nvtms, ijthsfwd, isbs, xexpbvs);
                    t0 = exp(vjsmfwd / nvtms);
                    ivjsmfwd = isbs * (t0 - xexpbvs / t0 + xexpbvs - 1.0);
                    sslpfwd = isbs * (t0 + xexpbvs / t0) / nvtms;
          
                    t2 = ijthsrev / isbs;
          
                    if (t2 < 1.0)
                      begin
                        t2 = 10.0;
                        $strobe("Warning: ijthsrev too small and set to 10 times isbsat.\n");
                      end     
          
                    vjsmrev = -bvs - nvtms * ln((t2 - 1.0) / xjbvs);
                    t1 = xjbvs * exp(-(bvs + vjsmrev) / nvtms);
                    ivjsmrev = isbs * (1.0 + t1);
                    sslprev = -isbs * t1 / nvtms;
                  end 
                default: $strobe("Specified diomod = %d not matched", diomod);
              endcase 
            end
          
          nvtmd = vtm * njd;
          
          if ((adeff <= 0.0) && (pdeff <= 0.0))
            isbd = 1.0e-14;
          else
            isbd = adeff * jsd_temp + pdeff * jswd_temp + weffcj * nf * jswgd_temp;
          
          if (isbd > 0.0)
            begin
              case(diomod)
                0:
                  begin
                    if ((bvd / nvtmd) > `EXP_THRESHOLD)
                      xexpbvd = xjbvd * `MIN_EXP;
                    else
                      xexpbvd = xjbvd * exp(-bvd / nvtmd);	
                  end
                1:
                  begin
                    vjdmfwd = get_vjm(nvtmd, ijthdfwd, isbd, 0.0);
                    ivjdmfwd = isbd * exp(vjdmfwd / nvtmd);
                  end
                2:
                  begin
                    if ((bvd / nvtmd) > `EXP_THRESHOLD)
                      begin
                        xexpbvd = xjbvd * `MIN_EXP;
                        tmp = `MIN_EXP;
                      end
                    else
                      begin
                        xexpbvd = exp(-bvd / nvtmd);
                        tmp = xexpbvd;
                        xexpbvd = xexpbvd * xjbvd;	
                      end
          
                    vjsmfwd = get_vjm(nvtmd, ijthdfwd, isbd, xexpbvd);
                    t0 = exp(vjdmfwd / nvtmd);
                    ivjdmfwd = isbd * (t0 - xexpbvd / t0 + xexpbvd - 1.0);
                    dslpfwd = isbd * (t0 + xexpbvd / t0) / nvtmd;
          
                    t2 = ijthdrev / isbd;
          
                    if (t2 < 1.0)
                      begin
                        t2 = 10.0;
                        $strobe("Warning: ijthdrev too small and set to 10 times idbsat.\n");
                      end     
          
                    vjdmrev = -bvd - nvtmd * ln((t2 - 1.0) / xjbvd);
                    t1 = xjbvd * exp(-(bvd + vjdmrev) / nvtmd);
                    ivjdmrev = isbd * (1.0 + t1);
                    dslprev = -isbd * t1 / nvtmd;
                  end 
                default: $strobe("Specified diomod = %d not matched", diomod);
              endcase 
            end
        end
      //*********************************//
      //****** End of initial_step ******//
      //*********************************//
   
      //****** Calculation of all equations to define all currents ******//
      //****** Definition of the tensions ******//
      vds   = type * V(drainp, sourcep);
      vgs   = type * V(gatep, sourcep);
      vbs   = type * V(bulkp, sourcep);
      vges  = type * V(gate, sourcep);
      vgms  = type * V(gatem, sourcep);
      vsbs  = type * V(sourceb, sourcep);
      vdbs  = type * V(drainb, sourcep);
      vses  = type * V(source, sourcep);
      vdes  = type * V(drain, sourcep);
      vbes  = type * V(bulk, sourcep);
      vgd   = vgs - vds;
      vbd   = vbs - vds;
      vgb   = vgs - vbs;
      vded  = vdes - vds;
      vgeg  = vges - vgs;
      vgmg  = vgms - vgs;
      vgmb  = vgms - vbs;
      vgegm = vgeg - vgmg;
      vbeb  = vbes - vbs;
      vdbd  = vdbs - vds;
      vbesb = vbes - vsbs;
      vbedb = vbes - vdbs;
      vbsb  = vbs - vsbs;
      vbdb  = vbs - vdbs;
      //***** Source/drain junction diode DC model (line 634-829, file b4ld.c) *****//
      vbs_jct = (rbodymod == 0) ? vbs : vsbs; 
      vbd_jct = (rbodymod == 0) ? vbd : vdbd;
      nvtms = vtm * njs;
      
      if ((aseff <= 0.0) && (pseff <= 0.0))
        isbs = 1.0e-14;
      else
        isbs = aseff * jss_temp + pseff * jsws_temp + weffcj * nf * jswgs_temp;
      
      if (isbs <= 0.0) 
        cbs = gmin * vbs_jct;
      else
        begin
          case(diomod)
            0:
              begin
                evbs = exp(vbs_jct / nvtms);
                t1 = xjbvs * exp(-(bvs + vbs_jct) / nvtms);
                cbs = isbs * (evbs + xexpbvs - t1 - 1.0) + gmin * vbs_jct;
              end
            1:
              begin
                t2 = vbs_jct / nvtms;
      
                if (t2 < -`EXP_THRESHOLD)
                  cbs = isbs * (`MIN_EXP - 1.0) + gmin * vbs_jct;
                else if (vbs_jct <= vjsmfwd)
                  begin
                    evbs = exp(t2);
                    cbs = isbs * (evbs - 1.0) + gmin * vbs_jct;
                  end
                else
                  begin
                    t0 = ivjsmfwd / nvtms;
                    cbs = ivjsmfwd - isbs + t0 * (vbs_jct - vjsmfwd) + gmin * vbs_jct;
                  end   	
              end
            2:
              begin
                if (vbs_jct < vjsmrev)
                  begin
                    t0 = vbs_jct / nvtms;
      
                    if (t0 < -`EXP_THRESHOLD)
                      evbs = `MIN_EXP;
                    else
                      evbs = exp(t0);
      
                    t1 = evbs - 1.0;
                    t2 = ivjsmrev + sslprev * (vbs_jct - vjsmrev);
                    cbs = t1 * t2 + gmin * vbs_jct;
                  end         
                else if (vbs_jct <= vjsmfwd)
                  begin
                    t0 = vbs_jct / nvtms;
      
                    if (t0 < -`EXP_THRESHOLD)
                      evbs = `MIN_EXP;
                    else
                      evbs = exp(t0);
      
                    t3 = (bvs + vbs_jct) / nvtms;
      
                    if (t1 > `EXP_THRESHOLD)
                      t2 = `MIN_EXP;
                    else
                      t2 = exp(-t1);
      
                    cbs = isbs * (evbs + xexpbvs - 1.0 - xjbvs * t2) + gmin * vbs_jct;
                  end 
                else
                  cbs = ivjsmfwd + sslpfwd * (vbs_jct - vjsmfwd) + gmin * vbs_jct;
              end
          endcase 
        end
      
      nvtmd = vtm * njd;
      
      if ((adeff <= 0.0) && (pdeff <= 0.0))
        isbd = 1.0e-14;
      else
        isbd = adeff * jsd_temp + pdeff * jswd_temp + weffcj * nf * jswgd_temp;
      
      if (isbd <= 0.0)
	cbd = gmin * vbd_jct;
      else
        begin
          case(diomod)
            0:
              begin
                evbd = exp(vbd_jct / nvtmd);
                t1 = xjbvd * exp(-(bvd + vbd_jct) / nvtmd);
                cbd = isbd * (evbd + xexpbvd - t1 - 1.0) + gmin * vbd_jct;
              end
            1:
              begin
                t2 = vbd_jct / nvtmd;
      
                if (t2 < -`EXP_THRESHOLD)
                begin
                  cbd = isbd * (`MIN_EXP - 1.0) + gmin * vbd_jct;
                    end
                else if (vbd_jct <= vjdmfwd)
                  begin
                    evbd = exp(t2);
                    cbd = isbd * (evbd - 1.0) + gmin * vbd_jct;
                  end
                else
                  begin
                    t0 = ivjdmfwd / nvtmd;
                    cbd = ivjdmfwd - isbd + t0 * (vbd_jct - vjdmfwd) + gmin * vbd_jct;
                  end   	
              end
            2:
              begin
                if (vbd_jct < vjdmrev)
                  begin
                    t0 = vbd_jct / nvtmd;
      
                    if (t0 < -`EXP_THRESHOLD)
                      evbd = `MIN_EXP;
                    else
                      evbd = exp(t0);
      
                    t1 = evbd - 1.0;
                    t2 = ivjdmrev + dslprev * (vbd_jct - vjdmrev);
                    cbd = t1 * t2 + gmin * vbd_jct;
                  end         
                else if (vbd_jct <= vjdmfwd)
                  begin
                    t0 = vbd_jct / nvtmd;
      
                    if (t0 < -`EXP_THRESHOLD)
                      evbd = `MIN_EXP;
                    else
                      evbd = exp(t0);
      
                    t1 = (bvd + vbd_jct) / nvtmd;
      
                    if (t1 > `EXP_THRESHOLD)
                      t2 = `MIN_EXP;
                    else
                      t2 = exp(-t1);
      
                    cbd = isbd * (evbd + xexpbvd - 1.0 - xjbvd * t2) + gmin * vbd_jct;
                  end 
                else
                  cbd = ivjdmfwd + dslpfwd * (vbd_jct - vjdmfwd) + gmin * vbd_jct;
              end
          endcase 
        end
      
      //***** Mode choice (line 831-844, file b4ld.c) *****//
      if (vds >= 0.0)
        begin
          mode = 1;
          Vds  = vds;
          Vbs  = vbs;
        end
      else
        begin
          mode = -1;
	  Vds  = -vds;
          Vbs  = vbd;
        end
      //***** Effective Vbs (line 846-863, file b4ld.c) *****//
      t0 = Vbs - vbsc - 0.001;
      t1 = sqrt(t0 * t0 - 0.004 * vbsc);
      
      if (t0 >= 0.0)
        vbseff = vbsc + 0.5 * (t0 + t1);
      else
        begin
          t2 = -0.002 / (t1 - t0);
          vbseff = vbsc * (1.0 + t2);
        end
      
      // Correction to forward body bias 
      t9 = 0.95 * phi;
      t0 = t9 - vbseff - 0.001;
      t1 = sqrt(t0 * t0 + 0.004 * t9);
      vbseff = t9 - 0.5 * (t0 + t1);
      
      //***** Calculation of phis (line 865-867, file b4ld.c) *****//
      phis = phi - vbseff;
      sqrtphis = sqrt(phis);
      
      //***** Threshold Voltage (line 878-969, file b4ld.c) *****// 
      xdep = xdep0 * sqrtphis / sqrt(phi);
      t3 = sqrt(xdep);
      v0 = vbi - phi;
      t0 = dvt2 * vbseff;
      
      if (t0 >= - 0.5)
        begin
          t1 = 1.0 + t0;
          t2 = dvt2;
        end
      else
        begin
          t4 = 1.0 / (3.0 + 8.0 * t0);
          t1 = (1.0 + 3.0 * t0) * t4; 
          t2 = dvt2 * t4 * t4;
        end  
      
      lt1 = factor1 * t3 * t1;
      t0 = dvt2w * vbseff;
      
      if (t0 >= - 0.5)
        begin
          t1 = 1.0 + t0;
          t2 = dvt2w;
        end  
      else
        begin
          t4 = 1.0 / (3.0 + 8.0 * t0);
          t1 = (1.0 + 3.0 * t0) * t4; 
          t2 = dvt2w * t4 * t4;
        end
        
      ltw = factor1 * t3 * t1;
      t0 = dvt1 * leff / lt1;
      
      if (t0 < `EXP_THRESHOLD)
        begin
          t1 = exp(t0);
          t2 = t1 - 1.0;
          t3 = t2 * t2;
          t4 = t3 + 2.0 * t1 * `MIN_EXP;
          theta0 = t1 / t4;
        end
      else
        theta0 = 1.0 / (`MAX_EXP - 2.0); /* 3.0 * `MIN_EXP omitted */
      
      thetavth = dvt0 * theta0;
      delt_vth = thetavth * v0;
      t0 = dvt1w * weff * leff / ltw;
      
      if (t0 < `EXP_THRESHOLD)
        begin
          t1 = exp(t0);
          t2 = t1 - 1.0;
          t3 = t2 * t2;
          t4 = t3 + 2.0 * t1 * `MIN_EXP;
          t5 = t1 / t4;
        end
      else
        t5 = 1.0 / (`MAX_EXP - 2.0); /* 3.0 * `MIN_EXP omitted */
      
      t0 = dvt0w * t5;
      t2 = t0 * v0;
      tratio =  $temperature / tnom - 1.0;
      t0 = sqrt(1.0 + lpe0 / leff);
      t1 = k1ox * (t0 - 1.0) * sqrtphi + (kt1 + kt1l / leff + kt2 * vbseff) * tratio;
      vth_narroww = toxe * phi / (weff + w0);
      t3 = eta0 + etab * vbseff;
      
      if (t3 < 1.0e-4)
        begin
          t9 = 1.0 / (3.0 - 2.0e4 * t3);
          t3 = (2.0e-4 - t3) * t9;
          t4 = t9 * t9;
        end
      else
        t4 = 1.0;
      
      ddibl_sft_dvd = t3 * theta0vb0;
      dibl_sft = ddibl_sft_dvd * Vds;
      lpe_vb = sqrt(1.0 + lpeb / leff);
      vth = type * vth0 + (k1ox * sqrtphis - k1 * sqrt(phi)) * lpe_vb - k2ox 
            * vbseff - delt_vth - t2 + (k3 + k3b * vbseff) * vth_narroww + t1 - dibl_sft;
      
      //***** Swing factor (line 978-998, file b4ld.c) *****//
      tmp1 = `EPSSI / xdep;
      tmp2 = nfactor * tmp1;
      tmp3 = cdsc + cdscb * vbseff + cdscd * Vds;
      tmp4 = (tmp2 + tmp3 * theta0 + cit) / coxe;
      
      if (tmp4 >= -0.5)
        n = 1.0 + tmp4;
      else
        begin
          t0 = 1.0 / (3.0 + 8.0 * tmp4);
          n = (1.0 + 3.0 * tmp4) * t0;
        end
      
      //***** Vth correction for Pocket Implant (line 1002-1024, file b4ld.c) *****//
      if (dvtp0 > 0.0)
        begin
          t0 = -dvtp1 * Vds;
      
          if (t0 < -`EXP_THRESHOLD)
            t2 = `MIN_EXP;
          else
            t2 = exp(t0);
      
          t3 = leff + dvtp0 * (1.0 + t2);
          t4 = vtm * ln(leff / t3);
          vth = vth - n * t4;
        end 
      
      //***** Poly Gate Si Depletion Effect (line 1028 & 4584-4612, file b4ld.c) *****//
      t0 = vfb + phi;
      
      if ((ngate > 1.0e18) && (ngate < 1.0e25) && (vgs > t0)) 
        begin
          t1 = 1.0e6 * `Charge * `EPSSI * ngate / (coxe * coxe);
          t8 = vgs - t0;
          t4 = sqrt(1.0 + 2.0 * t8 / t1);
          t2 = 2.0 * t8 / (t4 + 1.0);
          t3 = 0.5 * t2 * t2 / t1;
          t7 = 1.12 - t3 - 0.05;
          t6 = sqrt(t7 * t7 + 0.224);
          t5 = 1.12 - 0.5 * (t7 + t6);
          vgs_eff = vgs - t5;
        end    
      else 
        vgs_eff = vgs;
      
      if ((ngate > 1.0e18) && (ngate < 1.0e25) && (vgd > t0)) 
        begin
          t1 = 1.0e6 * `Charge * `EPSSI * ngate / (coxe * coxe);
          t8 = vgd - t0;
          t4 = sqrt(1.0 + 2.0 * t8 / t1);
          t2 = 2.0 * t8 / (t4 + 1.0);
          t3 = 0.5 * t2 * t2 / t1;
          t7 = 1.12 - t3 - 0.05;
          t6 = sqrt(t7 * t7 + 0.224);
          t5 = 1.12 - 0.5 * (t7 + t6);
          vgd_eff = vgd - t5;
        end    
      else 
        vgd_eff = vgd;
      	  
      if(mode > 0) 
        Vgs_eff = vgs_eff;
      else
        Vgs_eff = vgd_eff;
      
      vgst = Vgs_eff - vth;
      
      //***** Calculation of vgsteff (line 1051-1109, file b4ld.c) *****//
      t0 = n * vtm;
      t1 = mstar * vgst;
      t2 = t1 / t0;
      
      if (t2 > `EXP_THRESHOLD)
        t10 = t1;
      else if (t2 < -`EXP_THRESHOLD)
        begin
          t10 = vtm * ln(1.0 + `MIN_EXP);
          t10 = t10 * n;
        end
      else
        begin
          expvgst = exp(t2);
          t3 = vtm * ln(1.0 + expvgst);
          t10 = n * t3;
        end
      
      t1 = voffcbn - (1.0 - mstar) * vgst;
      t2 = t1 / t0;
      
      if (t2 < -`EXP_THRESHOLD)
        begin
          t3 = coxe * `MIN_EXP / cdep0;
          t9 = mstar + t3 * n;
        end
      else if (t2 > `EXP_THRESHOLD)
        begin
          t3 = coxe * `MAX_EXP / cdep0;
          t9 = mstar + t3 * n;
        end  
      else
        begin
          expvgst = exp(t2);
          t3 = coxe / cdep0;
          t4 = t3 * expvgst;
          t5 = t1 * t4 / t0;
          t9 = mstar + n * t4;
        end  
      
      vgsteff = t10 / t9;
      
      //***** Effective Channel Geometry (line 1111-1123, file b4ld.c) *****//
      t9 = sqrtphis - sqrt(phi);
      Weff = weff - 2.0 * (dwg * vgsteff + dwb * t9); 
      
      if (Weff < 2.0e-8) /* to avoid the discontinuity problem due to Weff */
        begin
          t0 = 1.0 / (6.0e-8 - 2.0 * Weff);
          Weff = 2.0e-8 * (4.0e-8 - Weff) * t0;
        end
      
      //***** Source/Drain Resistance (line 1126-1149, file b4ld.c) *****//
      if (rdsmod == 1)
        rds = 0.0;
      else
        begin
          t0 = 1.0 + prwg * vgsteff;
          t1 = prwb * t9;
          t2 = 1.0 / t0 + t1;
          t3 = t2 + sqrt(t2 * t2 + 0.01);
          t4 = rds0 * 0.5;
          rds = rdswmin + t3 * t4;
        end
      
      //***** Bulk Charge Effect (line 1151-1205, file b4ld.c) *****//
      t9 = 0.5 * k1ox * lpe_vb / sqrtphis;
      t1 = t9 + k2ox - k3b * vth_narroww;
      t9 = sqrt(xj * xdep);
      tmp1 = leff + 2.0 * t9;
      t5 = leff / tmp1; 
      tmp2 = a0 * t5;
      tmp3 = weff + b1; 
      tmp4 = b0 / tmp3;
      t2 = tmp2 + tmp4;
      t6 = t5 * t5;
      t7 = t5 * t6;
      abulk0 = 1.0 + t1 * t2; 
      t8 = ags * a0 * t7;
      dabulk_dvg = -t1 * t8;
      abulk = abulk0 + dabulk_dvg * vgsteff; 
      
      if (abulk0 < 0.1) /* added to avoid the problems caused by abulk0 */
        begin
          t9 = 1.0 / (3.0 - 20.0 * abulk0);
          abulk0 = (0.2 - abulk0) * t9;
        end
      
      if (abulk < 0.1)
        begin
          t9 = 1.0 / (3.0 - 20.0 * abulk);
          abulk = (0.2 - abulk) * t9;
        end
      
      t2 = keta * vbseff;
      
      if (t2 >= -0.9)
        t0 = 1.0 / (1.0 + t2);
      else
        begin
          t1 = 1.0 / (0.8 + t2);
          t0 = (17.0 + 20.0 * t2) * t1;
        end
      
      abulk = abulk * t0;
      abulk0 = abulk0 * t0;
      
      //***** Effective Mobility (line 1207-1255, file b4ld.c) *****//
      if (mobmod == 0)
        begin
          t0 = vgsteff + vth + vth;
          t2 = ua + uc * vbseff;
          t3 = t0 / toxe;
          t5 = t3 * (t2 + ub * t3);
        end
      else if (mobmod == 1)
        begin
          t0 = vgsteff + vth + vth;
          t2 = 1.0 + uc * vbseff;
          t3 = t0 / toxe;
          t4 = t3 * (ua + ub * t3);
          t5 = t4 * t2;
        end
      else
        begin
          t0 = (vgsteff + vtfbphi1) / toxe;
          t1 = exp(eu * ln(t0));
          t2 = ua + uc * vbseff;
          t5 = t1 * t2;
        end
      
      if (t5 >= -0.8)
        denomi = 1.0 + t5;
      else
        begin
          t9 = 1.0 / (7.0 + 10.0 * t5);
          denomi = (0.6 + t5) * t9;
        end
      
      ueff = u0temp / denomi;
      
      //***** Saturation Voltage (line 1257-1357, file b4ld.c) *****// 
      wvcox = Weff * vsattemp * coxe;
      wvcoxrds = wvcox * rds; 
      esat = 2.0 * vsattemp / ueff;
      esatl = esat * leff;
      t0 = -esatl /ueff;
      
      if (a1 == 0.0)
        Lambda = a2;
      else if (a1 > 0.0)
        begin
          t0 = 1.0 - a2;
          t1 = t0 - a1 * vgsteff - 0.0001;
          t2 = sqrt(t1 * t1 + 0.0004 * t0);
          Lambda = a2 + t0 - 0.5 * (t1 + t2);
        end
      else
        begin
          t1 = a2 + a1 * vgsteff - 0.0001;
          t2 = sqrt(t1 * t1 + 0.0004 * a2);
          Lambda = 0.5 * (t1 + t2);
        end
      
      vgst2vtm = vgsteff + 2.0 * vtm;
      
      if ((rds == 0.0) && (Lambda == 1.0))
        begin
          t0 = 1.0 / (abulk * esatl + vgst2vtm);
          t1 = t0 * t0;
          t2 = vgst2vtm * t0;
          t3 = esatl * vgst2vtm;
          vdsat = t3 * t0;
        end                     
      else
        begin
          t9 = abulk * wvcoxrds;
          t8 = abulk * t9;
          t7 = vgst2vtm * t9;
          t6 = vgst2vtm * wvcoxrds;
          t0 = 2.0 * abulk * (t9 - 1.0 + 1.0 / Lambda); 
          t1 = vgst2vtm * (2.0 / Lambda - 1.0) + abulk * esatl + 3.0 * t7;
          t2 = vgst2vtm * (esatl + 2.0 * t6);
          t3 = sqrt(t1 * t1 - 2.0 * t0 * t2);
          vdsat = (t1 - t3) / t0;
        end
      
      //***** Effective Vds (line 1359-1398, file b4ld.c) *****//
      t1 = vdsat - Vds - delta;
      t2 = sqrt(t1 * t1 + 4.0 * delta * vdsat);
      t0 = t1 / t2;
      t9 = 2.0 * delta;
      
      if (t1 >= 0.0)
        vdseff = vdsat - 0.5 * (t1 + t2);
      else
        begin
          t4 = t9 / (t2 - t1);
          t5 = 1.0 - t4;
          vdseff = vdsat * t5;
        end
       
      if (Vds == 0.0)
        vdseff = Vds;
      
      if (vdseff > Vds)
        vdseff = Vds;
      
      diffvds = Vds - vdseff;
      
      //***** Velocity Overshoot (line 1400-1439, file b4ld.c) *****//
      if((lambda != -99.0) && (lambda > 0.0) )
        begin 
          t1 =  leff * ueff;
          t2 = lambda / t1;
          t5 = 1.0 / (esat * litl);
          t6 = 1.0 + diffvds  * t5;
          t7 = 2.0 / (t6 * t6 + 1.0);
          t8 = 1.0 - t7;
          t10 = 1.0 + t2 * t8;
          esatl = esatl * t10;
        end
      
      //***** Early Voltage at vdsat (line 1441-1464 , file b4ld.c) *****//
      tmp4 = 1.0 - 0.5 * abulk * vdsat / vgst2vtm;
      t9 = wvcoxrds * vgsteff;
      t0 = esatl + vdsat + 2.0 * t9 * tmp4;
      t9 = wvcoxrds * abulk; 
      t1 = 2.0 / Lambda - 1.0 + t9; 
      vasat = t0 / t1;
      
      //***** Drain Current for Triode Region (line 1466-1523 , file b4ld.c) *****//
      tmp1 = vtfbphi2;
      tmp2 = 2.0e8 * toxp;
      dt0_dvg = 1.0 / tmp2;
      t0 = (vgsteff + tmp1) * dt0_dvg;
      tmp3 = exp(0.7 * ln(t0));
      t1 = 1.0 + tmp3;
      tcen = 1.9e-9 / t1;
      coxeff = `EPSSI * coxp / (`EPSSI + coxp * tcen);
      coxeffwovl = coxeff * Weff / leff;
      beta = ueff * coxeffwovl;
      abovvgst2vtm = abulk / vgst2vtm;
      t0 = 1.0 - 0.5 * vdseff * abovvgst2vtm;
      fgche1 = vgsteff * t0;
      t9 = vdseff / esatl;
      fgche2 = 1.0 + t9;
      gche = beta * fgche1 / fgche2;
      t0 = 1.0 + gche * rds;
      idl = gche / t0; 
      
      //***** Degradation Factor due to Pocket Implant (line 1525-1535 , file b4ld.c) *****//
      if (fprout <= 0.0)
       fp = 1.0;
      else
        begin
          t9 = fprout * sqrt(leff) / vgst2vtm;
          fp = 1.0 / (1.0 + t9);
        end
      
      //***** Early Voltage with Channel Length Modulatiom (line 1537-1585, file b4ld.c) *****//
      t8 = pvag / esatl;
      t9 = t8 * vgsteff;
      
      if (t9 > -0.9)
        pvagterm = 1.0 + t9;
      else
        begin
          t4 = 1.0 / (17.0 + 20.0 * t9);
          pvagterm = (0.8 + t9) * t4;
        end
      
      if ((pclm > 0.0) && (diffvds > 1.0e-10))
        begin
          t0 = 1.0 + rds * idl;
          t2 = vdsat / esat;
          t1 = leff + t2;
          cclm = fp * pvagterm * t0 * t1 / (pclm * litl);
          vaclm = cclm * diffvds;
        end         
      else
        begin
          vaclm = `MAX_EXP;
          cclm = `MAX_EXP;
        end
      
      //***** Early Voltage with Drain-INduced Barrier Lowering (line 1587-1635, file b4ld.c) *****//
      if (thetarout > 0.0)
        begin
          t8 = abulk * vdsat;
          t0 = vgst2vtm * t8;
          t1 = vgst2vtm + t8;
          t2 = thetarout;
          vadibl = (vgst2vtm - t0 / t1) / t2;
          t7 = pdiblb * vbseff;
      
          if (t7 >= -0.9)
            begin
              t3 = 1.0 / (1.0 + t7);
              vadibl = vadibl * t3;
            end
          else
            begin 
              t4 = 1.0 / (0.8 + t7);
              t3 = (17.0 + 20.0 * t7) * t4;
              vadibl = vadibl * t3;
            end
      
          vadibl = vadibl * pvagterm;
        end    
      else
        vadibl = `MAX_EXP;
      
      //***** Early Voltage with Drain-INduced Threshold Shift (line 1643-1664 , file b4ld.c) *****//
      if ((pditsd * Vds) > `EXP_THRESHOLD)
        t1 = `MAX_EXP;
      else
        t1 = exp(pditsd * Vds);
      
      if (pdits > 0.0)
        vadits = (1.0 + (1.0 + pditsl * leff) * t1) / pdits * fp;
      else
        vadits = `MAX_EXP;
      
      //***** Early Voltage with Substrate Current Induced Body Effect (line 1666-1685 , file b4.c) *****//
      if (pscbe2 > 0.0)
        begin
          if (diffvds > (pscbe1 * litl / `EXP_THRESHOLD))
            begin
              t0 =  pscbe1 * litl / diffvds;
              vascbe = leff * exp(t0) / pscbe2;
            end
          else
            vascbe = `MAX_EXP * leff/pscbe2;
        end  
      else
        vascbe = `MAX_EXP;
      
      //***** Drain Current (line 1687-1719 , file b4ld.c) *****//
      t9 = diffvds / vadibl;
      t0 = 1.0 + t9;
      idsa = idl * t0;
      t9 = diffvds / vadits;
      t0 = 1.0 + t9;
      idsa = idsa * t0;
      t0 = ln((vasat + vaclm) / vasat);
      t1 = t0 / cclm;
      t9 = 1.0 + t1;
      idsa = idsa * t9;
      
      //***** Substrate current (line 1722-1760, file b4ld.c) *****//
      tmp = alpha0 + alpha1 * leff;
      
      if ((tmp <= 0.0) || (beta0 <= 0.0))
        isub = 0.0;
      else
        begin
          t2 = tmp / leff;
      
          if (diffvds > beta0 / `EXP_THRESHOLD)
            begin
              t0 = -beta0 / diffvds;
              t1 = t2 * diffvds * exp(t0);
            end
          else
            begin
              t3 = t2 * `MIN_EXP;
              t1 = t3 * diffvds;
            end
      
          t4 = idsa * vdseff;
          isub = t1 * t4;
        end
      
      csub = isub;
      
      //***** Current Drain (line 1761-1785, file b4ld.c) *****//
      t9 = diffvds / vascbe;
      t0 = 1.0 + t9;
      ids = idsa * t0;
      cdrain = ids * vdseff;
      
      //***** Source End Velocity Limit (line 1787-1817, file b4ld.c) *****//
      if ((vtl != -99.0) && (vtl > 0.0))
        begin
          t12 = 1.0 / leff / coxeffwovl;
          t11 = t12 / vgsteff;
          t10 = -t11 / vgsteff;
          vs = cdrain * t11;
          t0 = 2 * `MM;
          t1 = vs / (vtl * tfactor);
          
          if (t1 <= 0)
            t2 = 1.0;
          else 
            t2 = 1.0 + exp(t0 * ln(t1));
          
          fsevl = 1.0 / exp(ln(t2)/ t0);
          cdrain = cdrain * fsevl; 
        end 
      
      //***** Rg calculation (line 1824-1858, file b4ld.c) *****// 
      if (rgatemod > 1)
        begin
          t9 = xrcrg2 * vtm;
          t0 = t9 * beta;
          gcrg = xrcrg1 * (t0 + ids);
      
          if (nf != 1.0)
            gcrg = gcrg * nf; 
      
          if (rgatemod == 2)
            begin
              t10 = grgeltd * grgeltd;
              t11 = grgeltd + gcrg;
              gcrg = grgeltd * gcrg / t11;
            end
        end
      
      //*************** Calculate bias-dependent external S/D resistance (line 1861-1939, file b4ld.c) ****************//
      if (rdsmod == 1.0)
        begin   /* rs(v) */
          t0      = vgs - vfbsd;
          t1      = sqrt(t0 * t0 + 1.0e-4);
          vgs_eff = 0.5 * (t0 + t1);
          t0      = 1.0 + prwg * vgs_eff;
          t1      = -prwb * vbs;
          t2      = 1.0 / t0 + t1;
          t3      = t2 + sqrt(t2 * t2 + 0.01);
          t4      = rs0 * 0.5;
          rs      = rswmin + t3 * t4;
          t0      = 1.0 + gsdiff * rs;
          gstot   = gsdiff / t0;
          /* rd(v) */
          t0      = vgd - vfbsd;
          t1      = sqrt(t0 * t0 + 1.0e-4);
          vgd_eff = 0.5 * (t0 + t1);
          t0      = 1.0 + prwg * vgd_eff;
          t1      = -prwb * vbd;
          t2      = 1.0 / t0 + t1;
          t3      = t2 + sqrt(t2 * t2 + 0.01);
          t4      = rd0 * 0.5;
          rd      = rdwmin + t3 * t4;
          t0      = 1.0 + gddiff * rd;
          gdtot   = gddiff / t0;
        end
      else
        begin
          gstot  = 0.0;
          gdtot  = 0.0;
        end
      
      //************** Calculate GIDL and GISL current (line 1941-2021, file b4ld.c) *******************//
      t0 = 3.0 * toxe;
      t1 = (vds - vgs_eff - egidl ) / t0;
      
      if ((agidl <= 0.0) || (bgidl <= 0.0) || (t1 <= 0.0) || (cgidl <= 0.0) || (vbd > 0.0))
        igidl  = 0.0;
      else 
        begin
          t2 = bgidl / t1;
      
          if (t2 < 100.0)
            igidl  = agidl * weffcj * t1 * exp(-t2);
          else
            begin
              igidl = agidl * weffcj * 3.720075976e-44;
              igidl = igidl * t1;
            end
      
          t4 = vbd * vbd;
          t5 = -vbd * t4;
          t6 = cgidl + t5;
          t7 = t5 / t6;
          t8 = 3.0 * cgidl * t4 / t6 / t6;
          igidl = igidl * t7;
        end
      
      t1 = (-vds - vgd_eff - egidl ) / t0;
      
      if ((agidl <= 0.0) || (bgidl <= 0.0) || (t1 <= 0.0) || (cgidl <= 0.0) || (vbs > 0.0))
        igisl = 0.0;
      else 
        begin
          t2 = bgidl / t1;
      
          if (t2 < 100.0) 
            igisl = agidl * weffcj * t1 * exp(-t2);
          else 
            begin
              igisl = agidl * weffcj * 3.720075976e-44;
              igisl = igisl * t1;
            end
      
          t4 = vbs * vbs;
          t5 = -vbs * t4;
          t6 = cgidl + t5;
          t7 = t5 / t6;
          t8 = 3.0 * cgidl * t4 / t6 / t6;
          igisl = igisl * t7;
        end
      
      //***************** Gate tunneling current (line 2024-2396, file b4ld.c) ****************//
      if ((igcmod != 0.0) || (igbmod != 0.0))
        begin
          v3 = vfbzb - Vgs_eff + vbseff - `DELTA_3;
      
          if (vfbzb<= 0.0)
            t0 = sqrt(v3 * v3 - 4.0 * `DELTA_3 * vfbzb);
          else
            t0 = sqrt(v3 * v3 + 4.0 * `DELTA_3 * vfbzb);
      
          vfbeff = vfbzb- 0.5 * (v3 + t0);
          voxacc = vfbzb- vfbeff;
      
          if (voxacc < 0.0)
            voxacc = 0.0;
      
          t0 = 0.5 * k1ox;
          t3 = Vgs_eff - vfbeff - vbseff - vgsteff;
      
          if (k1ox == 0.0)
            voxdepinv = 0.0;
          else if (t3 < 0.0)
            voxdepinv = -t3;
          else
            begin
              t1 = sqrt(t0 * t0 + t3);
              t2 = t0 / t1;
              voxdepinv = k1ox * (t1 - t0);
            end 
      
          voxdepinv = voxdepinv + vgsteff;
        end
      
      if (igcmod == 1.0)
        begin
          t0 = vtm * nigc;
          vxnvt = (Vgs_eff - type * vth0) / t0;
      
          if (vxnvt > `EXP_THRESHOLD)
            vaux = Vgs_eff - type * vth0;
          else if (vxnvt < -`EXP_THRESHOLD)
            vaux = t0 * ln(1.0 + `MIN_EXP);
          else
            begin
              expvxnvt = exp(vxnvt);
              vaux = t0 * ln(1.0 + expvxnvt);
            end   
        
          t2 = Vgs_eff * vaux;
          t11 = aechvb;
          t12 = bechvb;
          t3 = aigc * cigc - bigc;
          t4 = bigc * cigc;
          t5 = t12 * (aigc + t3 * voxdepinv - t4 * voxdepinv * voxdepinv);
        
          if (t5 > `EXP_THRESHOLD)
            t6 = `MAX_EXP;
          else if (t5 < -`EXP_THRESHOLD)
            t6 = `MIN_EXP;
          else
            t6 = exp(t5);
        
          igc = t11 * t2 * t6;
        
          if (pigcd != -99.0)
            modif_pigcd = pigcd;
          else
            begin
              t11 = bechvb * toxe;
              t12 = vgsteff + 1.0e-20;
              t13 = t11 / t12 / t12;
              modif_pigcd = t13 * (1.0 - 0.5 * vdseff / t12); 
            end
        
          t7 = -modif_pigcd * vdseff;
          t8 = t7 * t7 + 2.0e-4;
      
          if (t7 > `EXP_THRESHOLD)
            t9 = `MAX_EXP;
          else if (t7 < -`EXP_THRESHOLD)
            t9 = `MIN_EXP;
          else
            t9 = exp(t7);
        
          t0 = t8 * t8;
          t1 = t9 - 1.0 + 1.0e-4;
          t10 = (t1 - t7) / t8;
          igcs = igc * t10;
          t10 = (t7 * t9 - t1) / t8;
          igcd = igc * t10;
          t0 = vgs - vfbsd;
          vgs_eff = sqrt(t0 * t0 + 1.0e-4);
          t2 = vgs * vgs_eff;
          t11 = aechvbedge;
          t12 = bechvbedge;
          t3 = aigsd * cigsd - bigsd;
          t4 = bigsd * cigsd;
          t5 = t12 * (aigsd + t3 * vgs_eff - t4 * vgs_eff * vgs_eff);
      
          if (t5 > `EXP_THRESHOLD)
            t6 = `MAX_EXP;
          else if (t5 < -`EXP_THRESHOLD)
            t6 = `MIN_EXP;
          else
            t6 = exp(t5);
          
          igs = t11 * t2 * t6;
          t0 = vgd - vfbsd;
          vgd_eff = sqrt(t0 * t0 + 1.0e-4);
          t2 = vgd * vgd_eff;
          t5 = t12 * (aigsd + t3 * vgd_eff - t4 * vgd_eff * vgd_eff);
      
          if (t5 > `EXP_THRESHOLD)
            t6 = `MAX_EXP;
          else if (t5 < -`EXP_THRESHOLD)
            t6 = `MIN_EXP;
          else
            t6 = exp(t5);
      
          igd = t11 * t2 * t6;
        end  
      else
        begin
          igcs = 0.0;
          igcd = 0.0;
          igs = 0.0;
          igd = 0.0;
        end
      
      if (igbmod  == 1.0)
        begin
          t0 = vtm * nigbacc;
          t1 = -Vgs_eff + vbseff + vfbzb;
          vxnvt = t1 / t0;
       
          if (vxnvt > `EXP_THRESHOLD)
            vaux = t1;
          else if (vxnvt < -(`EXP_THRESHOLD))
          begin
            vaux = t0 * ln(1.0 + `MIN_EXP);
            end
          else
            begin
              expvxnvt = exp(vxnvt);
              vaux = t0 * ln(1.0 + expvxnvt);
            end 
         
          t2 = (Vgs_eff - vbseff) * vaux;
          t11 = 4.97232e-7 * weff * leff * toxratio;
          t12 = -7.45669e11 * toxe;
          t3 = aigbacc * cigbacc - bigbacc;
          t4 = bigbacc * cigbacc;
          t5 = t12 * (aigbacc + t3 * voxacc - t4 * voxacc * voxacc);
         
          if (t5 > `EXP_THRESHOLD)
            t6 = `MAX_EXP;
          else if (t5 < -`EXP_THRESHOLD)
            t6 = `MIN_EXP;
          else
            t6 = exp(t5);
         
          igbacc = t11 * t2 * t6;
          t0 = vtm * nigbinv;
          t1 = voxdepinv - eigbinv;
          vxnvt = t1 / t0;
      
          if (vxnvt > `EXP_THRESHOLD)
            vaux = t1;
          else if (vxnvt < -`EXP_THRESHOLD)
            vaux = t0 * ln(1.0 + `MIN_EXP);
          else
            begin
              expvxnvt = exp(vxnvt);
              vaux = t0 * ln(1.0 + expvxnvt);
            end
         
          t2 = (Vgs_eff - vbseff) * vaux;
          t11 = t11 * 0.75610;
          t12 = t12 * 1.31724;
          t3 = aigbinv * cigbinv - bigbinv;
          t4 = bigbinv * cigbinv;
          t5 = t12 * (aigbinv + t3 * voxdepinv - t4 * voxdepinv * voxdepinv);
         
          if (t5 > `EXP_THRESHOLD)
            t6 = `MAX_EXP;
          else if (t5 < -`EXP_THRESHOLD)
            t6 =`MIN_EXP;
          else
            t6 = exp(t5);
         
          igbinv = t11 * t2 * t6;
          igb = igbinv + igbacc;
        end
      else
        igb   = 0.0;   
        
      //***** Accounting of device fingers (line 2396-2452, file b4ld.c) *****//
      if (nf != 1.0)
        begin
          cdrain = cdrain * nf;
          csub   = csub * nf;
          igidl  = igidl * nf;
          igisl  = igisl * nf;
          igcs   = igcs * nf;
          igcd   = igcd * nf;
          igs    = igs * nf;
          igd    = igd * nf;
          igb    = igb * nf;
        end 
      //***** CV model (line 2484-3288, file b4ld.c) *****//
	 ccn = 1;
      
      if ((xpart < 0) || (ccn == 0))
        begin
          qgate  = 0.0;
          qdrn   = 0.0;
          qsrc   = 0.0;
          qbulk  = 0.0;
          qgmid  = 0.0;
        end
      else if (capmod == 0)
        begin
          if (vbseff < 0.0)
            vbseff = Vbs;
          else
            vbseff = phi - phis;
      
          vfb = vfbcv;
          vth = vfb + phi + k1ox * sqrtphis; 
          vgst = Vgs_eff - vth;
          coxwl = coxe * weffcv * leffcv * nf;
          arg1 = Vgs_eff - vbseff - vfb;
      
          if (arg1 <= 0.0)
            begin
              qgate = coxwl * arg1;
              qbulk = -qgate;
              qdrn = 0.0;
            end
          else if (vgst <= 0.0)
            begin
              t1 = 0.5 * k1ox;
              t2 = sqrt(t1 * t1 + arg1);
              qgate = coxwl * k1ox * (t2 - t1);
              qbulk = -qgate;
              qdrn = 0.0;
            end 
          else
            begin
              two_third_coxwl = 2.0 * (coxe * weffcv * leffcv * nf) / 3.0;
              abulkcv = abulk0 * abulkcvfactor;
              vdsat = vgst / abulkcv;
      
              if (xpart > 0.5)
                begin
                  /* 0/100 Charge partition model */
                  if (vdsat <= Vds)
                    begin
                      /* saturation region */
                      t1 = vdsat / 3.0;
                      qgate = coxwl * (Vgs_eff - vfb - phi - t1);
                      t2 = -two_third_coxwl * vgst;
                      qbulk = -(qgate + t2);
                      qdrn = 0.0;
                    end
                  else
                    begin
                      /* linear region */
                      alphaz = vgst / vdsat;
                      t1 = 2.0 * vdsat - Vds;
                      t2 = Vds / (3.0 * t1);
                      t3 = t2 * Vds;
                      t9 = 0.25 * coxwl;
                      t4 = t9 * alphaz;
                      t7 = 2.0 * Vds - t1 - 3.0 * t3;
                      t8 = t3 - t1 - 2.0 * Vds;
                      qgate = coxwl * (Vgs_eff - vfb - phi - 0.5 * (Vds - t3));
                      t10 = t4 * t8;
                      qdrn = t4 * t7;
                      qbulk = -(qgate + qdrn + t10);
                    end
                end
              else if (xpart < 0.5)
                begin
                  /* 40/60 Charge partition model */
                  if (Vds >= vdsat)
                    begin
                      /* saturation region */
                      t1 = vdsat / 3.0;
                      qgate = coxwl * (Vgs_eff - vfb - phi - t1);
                      t2 = -two_third_coxwl * vgst;
                      qbulk = -(qgate + t2);
                      qdrn = 0.4 * t2;
                    end 
                  else
                    begin
                      /* linear region  */
                      alphaz = vgst / vdsat;
                      t1 = 2.0 * vdsat - Vds;
                      t2 = Vds / (3.0 * t1);
                      t3 = t2 * Vds;
                      t9 = 0.25 * coxwl;
                      t4 = t9 * alphaz;
                      qgate = coxwl * (Vgs_eff - vfb - phi - 0.5 * (Vds - t3));
                      t6 = 8.0 * vdsat * vdsat - 6.0 * vdsat * Vds + 1.2 * Vds * Vds;
                      t8 = t2 / t1;
                      t7 = Vds - t1 - t8 * t6;
                      qdrn = t4 * t7;
                      t7 = 2.0 * (t1 + t3);
                      qbulk = -(qgate - t4 * t7);
                    end 
                end      
              else
                begin
                  /* 50/50 partitioning */
                  if (Vds >= vdsat)
                    begin
                      /* saturation region */
                      t1 = vdsat / 3.0;
                      qgate = coxwl * (Vgs_eff - vfb - phi - t1);
                      t2 = -two_third_coxwl * vgst;
                      qbulk = -(qgate + t2);
                      qdrn = 0.5 * t2;
                    end 
                  else
                    begin
                      /* linear region */
                      alphaz = vgst / vdsat;
                      t1 = 2.0 * vdsat - Vds;
                      t2 = Vds / (3.0 * t1);
                      t3 = t2 * Vds;
                      t9 = 0.25 * coxwl;
                      t4 = t9 * alphaz;
                      qgate = coxwl * (Vgs_eff - vfb - phi - 0.5 * (Vds - t3));
                      t7 = t1 + t3;
                      qdrn = -t4 * t7;
                      qbulk = - (qgate + qdrn + qdrn);
                    end
                end
            end
        end
      else
        begin
          if (vbseff < 0.0)
            vbseffcv = vbseff;
          else
            vbseffcv = phi - phis;
      
          coxwl = coxe * weffcv * leffcv * nf;
          t0 = vtm * n * noff;
          vgstnvt = (vgst - voffcv) / t0;
      
          if (vgstnvt > `EXP_THRESHOLD)
            vgsteff = vgst - voffcv;
          else if (vgstnvt < -`EXP_THRESHOLD)
            vgsteff = t0 * ln(1.0 + `MIN_EXP);
          else
            vgsteff = t0 * ln(1.0 + exp(vgstnvt));
      
          if (capmod == 1)
            begin
              v3 = vfbzb - Vgs_eff + vbseffcv - `DELTA_3;
      
              if (vfbzb <= 0.0)
                t0 = sqrt(v3 * v3 - 4.0 * `DELTA_3 * vfbzb);
              else
                t0 = sqrt(v3 * v3 + 4.0 * `DELTA_3 * vfbzb);
      
              t1 = 0.5 * (1.0 + v3 / t0);
              vfbeff = vfbzb - 0.5 * (v3 + t0);
              qac0 = coxwl * (vfbeff - vfbzb);
              t0 = 0.5 * k1ox;
              t3 = Vgs_eff - vfbeff - vbseffcv - vgsteff;
      
              if (k1ox == 0.0)
                begin
                  t1 = 0.0;
                  t2 = 0.0;
                end
              else if (t3 < 0.0)
                begin
                  t1 = t0 + t3 / k1ox;
                  t2 = coxwl;
                end
              else
                begin
                  t1 = sqrt(t0 * t0 + t3);
                  t2 = coxwl * t0 / t1;
                end
      
              qsub0 = coxwl * k1ox * (t1 - t0);
              abulkcv = abulk0 * abulkcvfactor;
              vdsatcv = vgsteff / abulkcv;
              t0 = vdsatcv - Vds - `DELTA_4;
              t1 = sqrt(t0 * t0 + 4.0 * `DELTA_4 * vdsatcv);
      
              if (t0 >= 0.0)
                vdseffcv = vdsatcv - 0.5 * (t0 + t1);
              else
                begin
                  t3 = (`DELTA_4 + `DELTA_4) / (t1 - t0);
                  t4 = 1.0 - t3;
                  t5 = vdsatcv * t3 / (t1 - t0);
                  vdseffcv = vdsatcv * t4;
                end 
      
              if (Vds == 0.0)
                vdseffcv = 0.0;
      
              t0 = abulkcv * vdseffcv;
              t1 = 12.0 * (vgsteff - 0.5 * t0 + 1.0e-20);
              t2 = t0 / t1;
              t3 = t0 * t2;
              qgate = coxwl * (vgsteff - 0.5 * t0 + t3);
              t7 = 1.0 - abulkcv;
              qbulk = coxwl * t7 * (0.5 * vdseffcv - t3);
              
              if (xpart > 0.5)
                begin
                  /* 0/100 Charge petition model */
                  t1 = t1 + t1;
                  qsrc = -coxwl * (0.5 * vgsteff + 0.25 * t0 - t0 * t0 / t1);
                end
              else if (xpart < 0.5)
                begin
                  /* 40/60 Charge petition model */
                  t1 = t1 / 12.0;
                  t2 = 0.5 * coxwl / (t1 * t1);
                  t3 = vgsteff * (2.0 * t0 * t0 / 3.0 + vgsteff * (vgsteff - 4.0 * t0 / 3.0)) - 2.0 * t0 * t0 * t0 / 15.0;
                  qsrc = -t2 * t3;
                end 
              else
                /* 50/50 Charge petition model */
                qsrc = -0.5 * (qgate + qbulk);
      
              qgate = qgate + qac0 + qsub0;
              qbulk = qbulk - (qac0 + qsub0);
              qdrn  = -(qgate + qbulk + qsrc);
            end 
          else if (capmod == 2)
            begin
              v3 = vfbzb - Vgs_eff + vbseffcv - `DELTA_3;
      
              if (vfbzb <= 0.0)
                t0 = sqrt(v3 * v3 - 4.0 * `DELTA_3 * vfbzb);
              else
                t0 = sqrt(v3 * v3 + 4.0 * `DELTA_3 * vfbzb);
      
              t1 = 0.5 * (1.0 + v3 / t0);
              vfbeff = vfbzb - 0.5 * (v3 + t0);
              tox = 1.0e8 * toxp;
              t0 = (Vgs_eff - vbseffcv - vfbzb) / tox;
              tmp = t0 * acde;
      
              if ((-`EXP_THRESHOLD < tmp) && (tmp < `EXP_THRESHOLD))
                tcen = ldeb * exp(tmp);
              else if (tmp <= -`EXP_THRESHOLD)
                tcen = ldeb * `MIN_EXP;
              else
                tcen = ldeb * `MAX_EXP;
      
              link = 1.0e-3 * toxp;
              v3 = ldeb - tcen - link;
              v4 = sqrt(v3 * v3 + 4.0 * link * ldeb);
              tcen = ldeb - 0.5 * (v3 + v4);
              ccen = `EPSSI / tcen;
              t2 = coxp / (coxp + ccen);
              coxeff = t2 * ccen;
              coxwlcen = coxwl * coxeff / coxe;
              qac0 = coxwlcen * (vfbeff - vfbzb);
              t0 = 0.5 * k1ox;
              t3 = Vgs_eff - vfbeff - vbseffcv - vgsteff;
              if (k1ox == 0.0)
                begin
                  t1 = 0.0;
                  t2 = 0.0;
                end
              else if (t3 < 0.0)
                begin
                  t1 = t0 + t3 / k1ox;
                  t2 = coxwlcen;
                end
              else
                begin
                  t1 = sqrt(t0 * t0 + t3);
                  t2 = coxwlcen * t0 / t1;
                end
      
              qsub0 = coxwlcen * k1ox * (t1 - t0);
      
              if (k1ox <= 0.0)
                begin
                  denomi = 0.25 * moin * vtm;
                  t0 = 0.5 * sqrtphi;
                end
              else
                begin
                  denomi = moin * vtm * k1ox * k1ox;
                  t0 = k1ox * sqrtphi;
                end
      
              t1 = 2.0 * t0 + vgsteff;
              deltaphi = vtm * ln(1.0 + t1 * vgsteff / denomi);
              tox = tox + tox;
              t0 = (vgsteff + vtfbphi2) / tox;
              tmp = exp(0.7 * ln(t0));
              t1 = 1.0 + tmp;
              tcen = 1.9e-9 / t1;
              ccen = `EPSSI / tcen;
              t0 = coxp / (coxp + ccen);
              coxeff = t0 * ccen;
              coxwlcen = coxwl * coxeff / coxe;
              abulkcv = abulk0 * abulkcvfactor;
              vdsatcv = (vgsteff - deltaphi) / abulkcv;
              t0 = vdsatcv - Vds - `DELTA_4;
              t1 = sqrt(t0 * t0 + 4.0 * `DELTA_4 * vdsatcv);
      
              if (t0 >= 0.0)
                vdseffcv = vdsatcv - 0.5 * (t0 + t1);
              else
                begin
                  t3 = (`DELTA_4 + `DELTA_4) / (t1 - t0);
                  t4 = 1.0 - t3;
                  vdseffcv = vdsatcv * t4;
                end 
      
              if (Vds == 0.0)
                vdseffcv = 0.0;
      
              t0 = abulkcv * vdseffcv;
              t1 = vgsteff - deltaphi;
              t2 = 12.0 * (t1 - 0.5 * t0 + 1.0e-20);
              t3 = t0 / t2;
              qgate = coxwlcen * (t1 - t0 * (0.5 - t3));
              t7 = 1.0 - abulkcv;
              qbulk = coxwlcen * t7 * (0.5 * vdseffcv - t0 * vdseffcv / t2);
              
              if (xpart > 0.5)
                /* 0/100 partition */
                qsrc = -coxwlcen * (t1 / 2.0 + t0 / 4.0 - 0.5 * t0 * t0 / t2);
              else if (xpart < 0.5)
                begin
                  /* 40/60 partition */
                  t2 = t2 / 12.0;
                  t3 = 0.5 * coxwlcen / (t2 * t2);
                  t4 = t1 * (2.0 * t0 * t0 / 3.0 + t1 * (t1 - 4.0 * t0 / 3.0)) - 2.0 * t0 * t0 * t0 / 15.0;
                  qsrc = -t3 * t4;
                end 
              else
                /* 50/50 partition */
                qsrc = -0.5 * qgate;
      
              qgate = qgate + qac0 + qsub0 - qbulk;
              qbulk = qbulk - (qac0 + qsub0);
              qdrn  = -(qgate + qbulk + qsrc);
            end 
        end
      
      if (ccn == 1)
        qsrc = -(qgate + qbulk + qdrn);
      
      //***** Junction Diode CV Model (line 3333-3450, file b4ld.c) *****//
      if (ccn == 1)
        begin
          czbd = cjd_temp * adeff;
          czbs = cjs_temp * aseff;
          czbdsw = cjswd_temp * pdeff;
          czbdswg = cjswgd_temp * weffcj * nf;
          czbssw = cjsws_temp * pseff;
          czbsswg = cjswgs_temp * weffcj * nf;
      
          /* Source Bulk Junction */
          if (vbs_jct == 0.0)
            qbs   = 0.0;
          else if (vbs_jct < 0.0)
            begin
              if (czbs > 0.0)
                begin
                  arg = 1.0 - vbs_jct / phibs;
      
                  if (mjs == 0.5)
                    sarg = 1.0 / sqrt(arg);
                  else
                    sarg = exp(-mjs * ln(arg));
      
                  qbs = phibs * czbs * (1.0 - arg * sarg) / (1.0 - mjs);
                end
              else
      	        qbs   = 0.0;
      
              if (czbssw > 0.0)
                begin
                  arg = 1.0 - vbs_jct / phibsws;
      
                  if (mjsws == 0.5)
                    sarg = 1.0 / sqrt(arg);
                  else
                    sarg = exp(-mjsws * ln(arg));
      
                  qbs = qbs + phibsws * czbssw * (1.0 - arg * sarg) / (1.0 - mjsws);
                end 
      
              if (czbsswg > 0.0)
                begin
                  arg = 1.0 - vbs_jct / phibswgs;
      
                  if (mjswgs == 0.5)
                    sarg = 1.0 / sqrt(arg);
                  else
                    sarg = exp(-mjswgs * ln(arg));
      
                  qbs = qbs + phibswgs * czbsswg * (1.0 - arg * sarg) / (1.0 - mjswgs);
                end
            end
          else
            begin
              t0  = czbs + czbssw + czbsswg;
              t1  = vbs_jct * (czbs * mjs / phibs + czbssw * mjsws / phibsws + czbsswg * mjswgs / phibswgs);
              qbs = vbs_jct * (t0 + 0.5 * t1);
            end
      
          /* Drain Bulk Junction */
          if (vbd_jct == 0.0)
            qbd   = 0.0;
          else if (vbd_jct < 0.0)
            begin
              if (czbd > 0.0)
                begin
                  arg = 1.0 - vbd_jct / phibd;
      
                  if (mjd == 0.5)
                    sarg = 1.0 / sqrt(arg);
                  else
                    sarg = exp(-mjd * ln(arg));
      
                  qbd = phibd * czbd * (1.0 - arg * sarg) / (1.0 - mjd);
                end
              else
                qbd   = 0.0;
      
              if (czbdsw > 0.0)
                begin
                  arg = 1.0 - vbd_jct / phibswd;
      
                  if (mjswd == 0.5)
                    sarg = 1.0 / sqrt(arg);
                  else
                    sarg = exp(-mjswd * ln(arg));
      
                  qbd = qbd + phibswd * czbdsw * (1.0 - arg * sarg) / (1.0 - mjswd);
                end
      
              if (czbdswg > 0.0)
                begin
                  arg = 1.0 - vbd_jct / phibswgd;
      
                  if (mjswgd == 0.5)
                    sarg = 1.0 / sqrt(arg);
                  else
                    sarg = exp(-mjswgd * ln(arg));
      
                  qbd = qbd + phibswgd * czbdswg * (1.0 - arg * sarg) / (1.0 - mjswgd);
                end
            end
          else
            begin
              t0  = czbd + czbdsw + czbdswg;
              t1  = vbd_jct * (czbd * mjd / phibd + czbdsw * mjswd / phibswd + czbdswg * mjswgd / phibswgd);    
              qbd = vbd_jct * (t0 + 0.5 * t1);
            end 
        end
      
      //***** Overlap capacitances (line 3519-3558, file b4ld.c) *****//
      if (ccn == 1)
        begin
          if (capmod == 0)
            begin
              cgdo = param_cgdo;
              cgso = param_cgso;
              qgdo = param_cgdo * vgd;
              qgso = param_cgso * vgs;
            end
          else /* For both capMod == 1 and 2 */
            begin
              t0 = vgd + `DELTA_1;
              t1 = sqrt(t0 * t0 + 4.0 * `DELTA_1);
              t2 = 0.5 * (t0 - t1);
              t3 = weffcv * cgdl;
              t4 = sqrt(1.0 - 4.0 * t2 / ckappad);
              cgdo = param_cgdo + t3 - t3 * (1.0 - 1.0 / t4) * (0.5 - 0.5 * t0 / t1);
              qgdo = (param_cgdo + t3) * vgd - t3 * (t2 + 0.5 * ckappad * (t4 - 1.0));
              t0 = vgs + `DELTA_1;
              t1 = sqrt(t0 * t0 + 4.0 * `DELTA_1);
              t2 = 0.5 * (t0 - t1);
              t3 = weffcv * cgsl;
              t4 = sqrt(1.0 - 4.0 * t2 / ckappas); 
              cgso = param_cgso + t3 - t3 * (1.0 - 1.0 / t4) * (0.5 - 0.5 * t0 / t1);
              qgso = (param_cgso + t3) * vgs - t3 * (t2 + 0.5 * ckappas * (t4 - 1.0));
            end
      
          if (nf != 1.0)
            begin
              cgdo = cgdo * nf;
              cgso = cgso * nf;
              qgdo = qgdo * nf;
              qgso = qgso * nf;
            end
      
          if (mode > 0)
            begin
              qdrn  = qdrn - qgdo;
      
      	      if (rgatemod == 3)
                begin
                  qgmb  = cgbo * vgmb;
                  qgmid = qgdo + qgso + qgmb;
                  qbulk = qbulk - qgmb;
                  qsrc  = -(qgate + qgmid + qbulk + qdrn);
                end
      	      else
                begin
                  qgb   = cgbo * vgb;
                  qgate = qgate + qgdo + qgso + qgb;
                  qbulk = qbulk - qgb;
                  qsrc  = -(qgate + qbulk + qdrn);
                end
            end
          else
            begin
              qsrc  = qdrn - qgso;
      	      if (rgatemod == 3)
                begin
                  qgmb  = cgbo * vgmb;
                  qgmid = qgdo + qgso + qgmb;
                  qbulk = qbulk - qgmb;
                  qdrn  = -(qgate + qgmid + qbulk + qsrc);
                end
      	      else
                begin
                  qgb   = cgbo * vgb;
                  qgate = qgate + qgdo + qgso + qgb;
                  qbulk = qbulk - qgb;
                  qdrn  = -(qgate + qbulk + qsrc);
                end
            end
      
          qd = qdrn - qbd;
          qs = qsrc - qbs;
          
          if (rbodymod == 0)
            qb = qbulk + qbd + qbs;
          else
            qb = qbulk;
        end
      
      //****** Here all currents going through the device are taken in account ******//
      I(sourceb, sourcep) <+ type * cbs;
      I(drainb, drainp)   <+ type * cbd;
      
      if (igbmod == 1)
        I(bulk, gatep) <+ -type * igb;
      if (igcmod == 1)
        begin
          I(source, gatep) <+ -type * (igs + igcs);
          I(drain, gatep)  <+ -type * (igd + igcd);
        end
       
      if (mode == 1)
        begin
          I(drainp , sourcep) <+ type * cdrain;
          I(drainp , bulkp)   <+ type * (csub + igidl);
          I(sourcep, bulkp)   <+ type * igisl;
        end
     else
        begin
          I(drainp , sourcep) <+ -type * cdrain;
          I(drainp , bulkp)   <+  type * (csub + igisl);
          I(sourcep, bulkp)   <+  type * igidl;
        end
      
      if (rdsmod == 0)
        begin
          V(source, sourcep) <+ 0;
          V(drainp, drain)   <+ 0;
        end
      else
        begin
          I(drain, drainp)   <+ type * gdtot * vded;
          I(source, sourcep) <+ type * gstot * vses;
        end
      
      if (rgatemod == 0)
        begin
          V(gate, gatem)  <+ 0;
          V(gatem, gatep) <+ 0;
        end
      else if (rgatemod == 1)
        begin
          V(gatep, gatem) <+ 0;
          I(gate, gatem)  <+ type * grgeltd * vgeg;
        end
      else if (rgatemod == 2)
        begin
          V(gatep, gatem) <+ 0;
          I(gate, gatem)  <+ type * gcrg * vgeg;
        end
      else
        begin
          I(gatem, gatep)   <+ type * gcrg * vgmg;
          I(gate, gatem)    <+ type * grgeltd * vgegm;
        `ifdef RGATE3  
          I(gatem, sourcep) <+ type * ddt(qgmid);
        `endif
        end
      
      if (rbodymod == 0)
        begin
          V(bulk, bulkp)       <+ 0;
          V(bulkp, sourceb)    <+ 0;
          V(bulkp, drainb)     <+ 0;
        end
      else
        begin
          I(sourceb, bulkp)   <+ -type * grbps * vbsb;
          I(drainb, bulkp)    <+ -type * grbpd * vbdb;
          I(bulk, bulkp)      <+  type * grbpb * vbeb;
          I(bulk, sourceb)    <+  type * grbsb * vbesb;
          I(bulk, drainb)     <+  type * grbdb * vbedb;
        end
      
      I(gatep, sourcep)   <+ type * ddt(qgate);
      I(drainp, sourcep)  <+ type * ddt(qd);
      I(bulkp, sourcep)   <+ type * ddt(qb);
      
      `ifdef RBODY
      if (rbodymod == 1)
        begin
          I(drainb,  sourcep) <+ type * ddt(qbd);
          I(sourceb, sourcep) <+ type * ddt(qbs);
        end
      `endif
      
    end
  endmodule
											
										rubberband_example_13.in
Testing CAPTAB option with BSIM4 and CAP models * based on EXAMPLE7. TWO-BIT MOSFET ADDER * .option BRIEF + itl1=400 + itl4=50 .verilog bsim4.va .GLOBAL VSS .SUBCKT NAND 1 2 3 YVLG1 3 1 VSS VSS MOSP W=10U L=1.3U YVLG2 3 2 VSS VSS MOSP W=10U L=1.3U YVLG3 3 1 10 0 MOSN W=5U L=1.3U YVLG4 10 2 0 0 MOSN W=5U L=1.3U C1 3 0 55FF .ENDS NAND ************* INA INB CIN OUT COUT .SUBCKT ONEBIT 1 2 3 4 5 X1 1 2 7 NAND X2 1 7 8 NAND X3 2 7 9 NAND X4 8 9 10 NAND X5 3 10 11 NAND X6 3 11 12 NAND X7 10 11 13 NAND X8 12 13 4 NAND X9 11 7 5 NAND .ENDS ONEBIT .SUBCKT TWOBIT INA1 INB1 INA2 INB2 BIT1 BIT2 CIN COUT X1 INA1 INB1 CIN BIT1 10 ONEBIT X2 INA2 INB2 10 BIT2 COUT ONEBIT .ENDS TWOBIT VINA1 1 0 PULSE( 0.8 0.0 0 10NS 10NS 140NS 300NS) VINB1 2 0 PULSE( 0.8 0.0 0 10NS 10NS 290NS 600NS) VINA2 3 0 PULSE( 0.8 0.0 0 10NS 10NS 140NS 300NS) VINB2 4 0 PULSE( 0.8 0.0 0 10NS 10NS 290NS 600NS) VINCIN CIN 0 PULSE( 0.0 0.8 0 10NS 10NS 65NS 150NS) VSS VSS 0 DC 0.8 X1 1 2 3 4 BIT1 BIT2 CIN COUT TWOBIT * .TRAN 10NS 240NS .MEASURE TRAN MAX_BIT1 MAX V(BIT1) .MEASURE TRAN MAX_COUT1 MAX V(X1.10) .MEASURE TRAN DEL_BIT1_CIN DELAY + V(CIN) FALL=1 VAL=0.4 + TARG=V(BIT1) FALL=1 VAL0=0 VAL1=MAX_BIT1 .MEASURE TRAN FALL_BIT2 WAVE + V(BIT2) FROM=40NS FALL=1 VAL0=0 VAL1=0.8 .MEASURE TRAN FALL_BIT1 WAVE + V(BIT1) FROM=40NS FALL=1 VAL0=0 VAL1=0.8 .MEASURE TRAN WIDTH_COUT1 DELAY + V(X1.10) RISE=1 VAL0=0 VAL1=MAX_COUT1 + TARG=V(X1.10) FALL=1 VAL0=0 VAL1=MAX_COUT1 .OPTIONS ACCT NOMOD .MODEL MOSN VLG MODULE=mosfet TYPE=1 TNOM=27 + EPSROX=3.9 TOXE=3e-9 TOXP=3e-9 TOXM=3e-9 TOXREF=3e-9 DTOX=0 + MOBMOD=0 RDSMOD=0 IGCMOD=0 IGBMOD=0 DIOMOD=1 CAPMOD=2 + RGATEMOD=0 RBODYMOD=0 *+ TRNQSMOD=0 ACNQSMOD=0 FNOIMOD=1 TNOIMOD=0 + PERMOD=1 GEOMOD=1 *+ PARAMCHK=0 BINUNIT=1 + NDEP=1.7e17 NGATE=9e19 VTH0=0.42 EACH=10% PHIN=0 XJ=9e-8 NSD=1e20 + VFB=-0.76 K1=0.33 K2=-1.87e-2 K3=80 K3b=0 W0=2.5e-6 KT1=-0.11 + KT1L=0 KT2=2.2e-2 LPE0=5.75e-8 LPEB=2.3e-10 DVT0=2.89 DVT1=0.53 + DVT2=-3.2e-2 DVT0W=0 DVT1W=0 DVT2W=0 DVTP0=7.32e-7 DVTP1=0.12 + DSUB=0.058 ETA0=0.001 ETAB=0 + CDSC=2.58e-4 CDSCD=6.1e-8 CDSCB=0 CIT=0 NFACTOR=1.1 VOFF=-0.02 VOFFL=0 + MINV=0.05 VSAT=9.2e4 AT=3.3e4 U0=4.19e-2 UA=8.7e-16 UB=3.06e-18 + UC=4.6e-13 UTE=-1.5 UA1=4.31e-9 UB1=7.61e-18 UC1=-5.6e-11 + RDSW=369.4 RDW=184.7 RSW=184.7 PRWG=3.22e-8 PRWB=6.8e-11 WR=1 + RDSWMIN=0 RDWMIN=0 RSWMIN=0 PRT=0 A0=1.1 AGS=1e-20 B0=-1e-20 + B1=0 KETA=-0.047 A1=0 A2=1 DELTA=0.014 PVAG=1e-20 PCLM=6.28e-4 + PDITS=0.2 PDITSL=2.3e6 PDITSD=0.23 FPROUT=0.2 *+ PDIBLCB=3.4e-8 PDIBLC1=0.81 PDIBLC2=9.84e-06 + DROUT=0.56 PSCBE1=8.14e8 PSCBE2=9.58e-07 + LINT=5e-9 WINT=5e-9 DWG=0 DWB=0 WL=0 WW=0 WWL=0 + WLN=1 WWN=1 LL=0 LW=0 LWL=0 LLN=1 LWN=1 + XGW=3e-7 XGL=4e-8 NGCON=1 RSHG=0.4 XRCRG1=12 XRCRG2=5 + DMCG=5e-6 DMCI=5e-6 DMDG=5e-6 DMCGT=6e-7 DWJ=4.5e-8 RSH=6 + RBDB=15 RBPD=15 RBSB=15 RBPS=15 RBPB=5 GBMIN=1e-10 + CGSO=7.43e-10 CGDO=7.43e-10 CGBO=2.56e-11 CGSL=1e-14 CGDL=1e-14 + CKAPPAS=0.5 CKAPPAD=0.5 NOFF=0.9 VOFFCV=0.02 ACDE=1 MOIN=15 XPART=0 + JSS=1e-4 JSWS=1e-11 JSWGS=1e-10 CJS=5e-4 CJSWS=5e-10 CJSWGS=3e-10 + PBS=1 PBSWS=1 PBSWGS=1 NJS=1 XTIS=3 MJS=0.5 MJSWS=0.33 MJSWGS=0.33 + XJBVS=1 BVS=10 IJTHSFWD=1e-2 IJTHSREV=1e-3 + JSD=1e-4 JSWD=1e-11 JSWGD=1e-10 CJD=5e-4 CJSWD=5e-10 CJSWGD=5e-10 + PBD=1 PBSWD=1 PBSWGD=1 NJD=1 XTID=3 MJD=0.5 MJSWD=0.33 MJSWGD=0.33 + XJBVD=1 BVD=10 IJTHDFWD=1e-2 IJTHDREV=1e-3 + TCJ=1e-3 TCJSW=1e-3 TCJSWG=1e-3 TPB=0.005 TPBSW=0.005 TPBSWG=0.005 + ALPHA0=7.4e-2 ALPHA1=0.005 BETA0=30 AGIDL=2e-4 BGIDL=2.1e9 CGIDL=2e-4 + EGIDL=0.8 AIGC=0.043 BIGC=0.0054 CIGC=0.0075 AIGSD=0.043 BIGSD=0.0054 + CIGSD=0.0075 DLCIG=2.2e-8 AIGBACC=0.043 BIGBACC=0.0054 CIGBACC=0.0075 + AIGBINV=0.35 BIGBINV=0.03 CIGBINV=0.006 NIGC=1 NIGBACC=1 NIGBINV=3 + NTOX=1 EIGBINV=1.1 PIGCD=1.0 POXEDGE=1 * .MODEL MOSP VLG MODULE=mosfet TYPE=-1 TNOM=27 + EPSROX=3.9 TOXE=3e-9 TOXP=3e-9 TOXM=3e-9 TOXREF=3e-9 DTOX=0 + MOBMOD=0 RDSMOD=0 IGCMOD=0 IGBMOD=0 DIOMOD=1 CAPMOD=2 + RGATEMOD=0 RBODYMOD=0 *+ TRNQSMOD=0 ACNQSMOD=0 FNOIMOD=1 TNOIMOD=0 + PERMOD=1 GEOMOD=1 *+ PARAMCHK=0 BINUNIT=1 + NDEP=1.7e17 NGATE=9e19 VTH0=-0.42 EACH=15% PHIN=0 XJ=9e-8 NSD=1e20 + VFB=0.76 K1=0.33 K2=-1.87e-2 K3=80 K3b=0 W0=2.5e-6 KT1=-0.11 + KT1L=0 KT2=2.2e-2 LPE0=5.75e-8 LPEB=2.3e-10 DVT0=2.89 DVT1=0.53 + DVT2=-3.2e-2 DVT0W=0 DVT1W=0 DVT2W=0 DVTP0=7.32e-7 DVTP1=0.12 + DSUB=0.058 ETA0=0.001 ETAB=0 + CDSC=2.58e-4 CDSCD=6.1e-8 CDSCB=0 CIT=0 NFACTOR=1.1 VOFF=-0.02 VOFFL=0 + MINV=0.05 VSAT=9.2e4 AT=3.3e4 U0=4.19e-2 UA=8.7e-16 UB=3.06e-18 + UC=4.6e-13 UTE=-1.5 UA1=4.31e-9 UB1=7.61e-18 UC1=-5.6e-11 + RDSW=369.4 RDW=184.7 RSW=184.7 PRWG=3.22e-8 PRWB=6.8e-11 WR=1 + RDSWMIN=0 RDWMIN=0 RSWMIN=0 PRT=0 A0=1.1 AGS=1e-20 B0=-1e-20 + B1=0 KETA=-0.047 A1=0 A2=1 DELTA=0.014 PVAG=1e-20 PCLM=6.28e-4 + PDITS=0.2 PDITSL=2.3e6 PDITSD=0.23 FPROUT=0.2 *+ PDIBLCB=3.4e-8 PDIBLC1=0.81 PDIBLC2=9.84e-06 + DROUT=0.56 PSCBE1=8.14e8 PSCBE2=9.58e-07 + LINT=5e-9 WINT=5e-9 DWG=0 DWB=0 WL=0 WW=0 WWL=0 + WLN=1 WWN=1 LL=0 LW=0 LWL=0 LLN=1 LWN=1 + XGW=3e-7 XGL=4e-8 NGCON=1 RSHG=0.4 XRCRG1=12 XRCRG2=5 + DMCG=5e-6 DMCI=5e-6 DMDG=5e-6 DMCGT=6e-7 DWJ=4.5e-8 RSH=6 + RBDB=15 RBPD=15 RBSB=15 RBPS=15 RBPB=5 GBMIN=1e-10 + CGSO=7.43e-10 CGDO=7.43e-10 CGBO=2.56e-11 CGSL=1e-14 CGDL=1e-14 + CKAPPAS=0.5 CKAPPAD=0.5 NOFF=0.9 VOFFCV=0.02 ACDE=1 MOIN=15 XPART=0 + JSS=1e-4 JSWS=1e-11 JSWGS=1e-10 CJS=5e-4 CJSWS=5e-10 CJSWGS=3e-10 + PBS=1 PBSWS=1 PBSWGS=1 NJS=1 XTIS=3 MJS=0.5 MJSWS=0.33 MJSWGS=0.33 + XJBVS=1 BVS=10 IJTHSFWD=1e-2 IJTHSREV=1e-3 + JSD=1e-4 JSWD=1e-11 JSWGD=1e-10 CJD=5e-4 CJSWD=5e-10 CJSWGD=5e-10 + PBD=1 PBSWD=1 PBSWGD=1 NJD=1 XTID=3 MJD=0.5 MJSWD=0.33 MJSWGD=0.33 + XJBVD=1 BVD=10 IJTHDFWD=1e-2 IJTHDREV=1e-3 + TCJ=1e-3 TCJSW=1e-3 TCJSWG=1e-3 TPB=0.005 TPBSW=0.005 TPBSWG=0.005 + ALPHA0=7.4e-2 ALPHA1=0.005 BETA0=30 AGIDL=2e-4 BGIDL=2.1e9 CGIDL=2e-4 + EGIDL=0.8 AIGC=0.043 BIGC=0.0054 CIGC=0.0075 AIGSD=0.043 BIGSD=0.0054 + CIGSD=0.0075 DLCIG=2.2e-8 AIGBACC=0.043 BIGBACC=0.0054 CIGBACC=0.0075 + AIGBINV=0.35 BIGBINV=0.03 CIGBINV=0.006 NIGC=1 NIGBACC=1 NIGBINV=3 + NTOX=1 EIGBINV=1.1 PIGCD=1.0 POXEDGE=1 .END
Graphics


