// ********************************************************
// ********************************************************
// *** BSIM-CMG 110.0.0 released by Sourabh Khandelwal on 01/01/2016 ****/
// *  BSIM Common Multi-Gate Model Equations (Verilog-A)
// ********************************************************
//
// ********************************************************
// * Copyright 2016 Regents of the University of California.
// * All rights reserved.
// *
// * Project Director: Prof. Chenming Hu.
// * Authors: Sriramkumar V., Navid Paydavosi, Juan Duarte, Sourabh Khandelwal, Darsen Lu,
// *          Chung-Hsun Lin, Mohan Dunga, Shijing Yao,
// *          Ali Niknejad, Chenming Hu
// ********************************************************
// ********************************************************
// *   NONDISCLOSURE STATEMENT
// Software is distributed as is, completely without warranty or service
// support. The University of California and its employees are not liable
// for the condition or performance of the software.
// The University of California owns the copyright and grants users a perpetual,
// irrevocable, worldwide, non-exclusive, royalty-free license with
// respect to the software as set forth below.
// The University of California hereby disclaims all implied warranties.
// The University of California grants the users the right to modify, copy,
// and redistribute the software and documentation, both within the user's
// organization and externally, subject to the following restrictions
// 1. The users agree not to charge for the University of California code
// itself but may charge for additions, extensions, or support.
// 2. In any product based on the software, the users agree to acknowledge
// the University of California that developed the software. This
// acknowledgment shall appear in the product documentation.
// 3. The users agree to obey all U.S. Government restrictions governing
// redistribution or export of the software.
// 4. The users agree to reproduce any copyright notice which appears on
// the software on any copy or modification of such made available
// to others
// Agreed to on __Jan 01, 2016__________________
// By: ___University of California, Berkeley____
//     ___Chenming Hu_____________________
//     ___Professor in Graduate School _______
// ********************************************************

// Clamped Exponential Function
analog function real lexp;
    input x;
    real x;

    begin
        if (x > `EXPL_THRESHOLD) begin
            lexp  =  `MAX_EXPL * (1.0 + x - `EXPL_THRESHOLD);
        end else if (x < -`EXPL_THRESHOLD) begin
            lexp  =  `MIN_EXPL;
        end else begin
            lexp  =  exp(x);
        end
    end
endfunction

// Clamped log Function
analog function real lln;
    input x;
    real x;

    begin
        lln  =  ln(max(x, `N_MINLOG));
    end
endfunction

// Hyperbolic Smoothing Function
analog function real hypsmooth;
    input x, c;
    real x, c;

    begin
        hypsmooth  =  0.5 * (x + sqrt(x * x + 4.0 * c * c));
    end
endfunction

// Hyperbolic Smoothing max Function
analog function real hypmax;
    input x, xmin, c;
    real x, xmin, c;

    begin
        hypmax  =  xmin + 0.5 * (x - xmin - c + sqrt((x - xmin - c) * (x - xmin - c) - 4.0 * xmin * c));
    end
endfunction

// Temperature Dependence Type
analog function real Tempdep;
    input PARAML, PARAMT, DELTEMP, TEMPMOD;
    real PARAML, PARAMT, DELTEMP, TEMPMOD;

    begin
        if (TEMPMOD != 0) begin
            Tempdep  =  PARAML + hypmax(PARAMT * DELTEMP, -PARAML, 1.0e-6);
        end else begin
            Tempdep  =  PARAML * hypsmooth(1.0 + PARAMT * DELTEMP - 1.0e-6, 1.0e-3);
        end
    end
endfunction

// Node Definitions
`ifdef __RGATEMOD__
    `define GateEdgeNode ge
`else
    `define GateEdgeNode g
`endif
`ifdef __NQSMOD1__
    `define IntrinsicGate gi
`else
    `define IntrinsicGate `GateEdgeNode
`endif

// ***************************
// *   Instance Parameters   *
// ***************************

// Note: Some instance parameters are also model parameters. Please refer to the technical note for details.
`IPRco( L              ,3.0e-8         ,"m"           ,1.0e-9      ,inf         ,"Designed gate length" )
`IPRco( D              ,4.0e-8         ,"m"           ,1.0e-9      ,inf         ,"Diameter of the cylinder (GEOMOD=3)" )
`IPRco( TFIN           ,1.5e-8         ,"m"           ,1.0e-9      ,inf         ,"Body (fin) thickness" )
`IPRco( FPITCH         ,8.0e-8         ,"m"           ,TFIN        ,inf         ,"Fin pitch" )
`IPIco( NF             ,1              ,""            ,1           ,inf         ,"Number of fingers" )
`IPRoz( NFIN           ,1.0            ,""            ,"Number of fins per finger (real number enables optimization)" )
`IPIcc( NGCON          ,1              ,""            ,1           ,2           ,"Number of gate contact (1 or 2 sided)" )
`IPRcz( ASEO           ,0.0            ,"m^2"         ,"Source-to-substrate overlap area through oxide" )
`IPRcz( ADEO           ,0.0            ,"m^2"         ,"Drain-to-substrate overlap area through oxide" )
`IPRcz( PSEO           ,0.0            ,"m"           ,"Perimeter of source-to-substrate overlap region through oxide" )
`IPRcz( PDEO           ,0.0            ,"m"           ,"Perimeter of drain-to-substrate overlap region through oxide" )
`IPRcz( ASEJ           ,0.0            ,"m^2"         ,"Source junction area (BULKMOD=1 or 2)" )
`IPRcz( ADEJ           ,0.0            ,"m^2"         ,"Drain junction area (BULKMOD=1 or 2)" )
`IPRcz( PSEJ           ,0.0            ,"m"           ,"Source-to-substrate PN junction perimeter (BULKMOD=1 or 2)" )
`IPRcz( PDEJ           ,0.0            ,"m"           ,"Drain-to-substrate PN junction perimeter (BULKMOD=1 or 2)" )
`IPRcz( COVS           ,0.0            ,"F/m"         ,"Constant gate-to-source overlap capacitance (CGEOMOD=1)" )
`IPRcz( COVD           ,COVS           ,"F/m"         ,"Constant gate-to-drain overlap capacitance (CGEOMOD=1)" )
`IPRcz( CGSP           ,0.0            ,"F/m"         ,"Constant gate-to-source fringe capacitance (CGEOMOD=1)" )
`IPRcz( CGDP           ,0.0            ,"F/m"         ,"Constant gate-to-drain fringe capacitance (CGEOMOD=1)" )
`IPRcz( CDSP           ,0.0            ,"F"           ,"Constant drain-to-source fringe capacitance (all CGEOMOD)" )
`IPRcz( NRS            ,0.0            ,""            ,"Number of source diffusion squares" )
`IPRcz( NRD            ,0.0            ,""            ,"Number of source diffusion squares" )
`IPRoz( LRSD           ,L              ,"m"           ,"Length of the source/drain" )
`IPRoz( NFINNOM        ,1.0            ,""            ,"Nominal number of fins per finger" )

// Variability Handles
`MPRnb( XL             ,0.0            ,"m"           ,"L offset for channel length due to mask/etch effect" )
`MPRnb( DTEMP          ,0.0            ,"Celsius"     ,"Variability in device temperature" )
`MPRnb( DELVTRAND      ,0.0            ,"V"           ,"Variability in Vth" )
`MPRcz( U0MULT         ,1.0            ,""            ,"Variability in carrier mobility" )
`MPRcz( IDS0MULT       ,1.0            ,""            ,"Variability in drain current for miscellaneous reasons" )

// ************************
// *  Model Parameters    *
// ************************
`MPIcc( DEVTYPE        ,`ntype         ,""            ,`ptype      ,`ntype      ,"0: PMOS; 1: NMOS" )
`MPIcc( TYPE           ,DEVTYPE        ,""            ,`ptype      ,`ntype      ,"0: PMOS; 1: NMOS" )
`MPIcc( BULKMOD        ,0              ,""            ,0           ,2           ,"0: SOI multi-gate; 1: Bulk multi-gate; 2: for decoupled bulk multi-gate" )
`MPIcc( GEOMOD         ,0              ,""            ,0           ,4           ,"0: Double gate; 1: Triple gate; 2: Quadruple gate; 3: Cylindrical gate; 4: Unified fin Shape" )
`MPIcc( CGEO1SW        ,0              ,""            ,0           ,1           ,"For CGEOMOD=1 only, this switch enables the parameters COVS, COVD, CGSP, and CGDP to be in F per fin, per gate-finger, per unit channel width" )
`MPIcc( RDSMOD         ,0              ,""            ,0           ,2           ,"0: Internal S/D resistance model; 1: External S/D resistance model; 2: Both bias dependent and independent part of S/D resistance internal" )
`MPIcc( ASYMMOD        ,0              ,""            ,0           ,1           ,"0: Turn off asymmetry model - forward mode parameters used; 1: Turn on asymmetry model" )
`MPIcc( IGCMOD         ,0              ,""            ,0           ,1           ,"0: Turn off Igc, Igs and Igd; 1: Turn on Igc, Igs and Igd" )
`MPIcc( IGBMOD         ,0              ,""            ,0           ,1           ,"0: Turn off Igb; 1: Turn on Igb" )
`MPIcc( GIDLMOD        ,0              ,""            ,0           ,1           ,"0: Turn off GIDL/GISL current; 1: Turn on GIDL/GISL current" )
`MPIcc( IIMOD          ,0              ,""            ,0           ,2           ,"0: Turn off impact ionization current; 1: BSIM4-based model; 2: BSIMSOI-based model" )
`MPIcc( TNOIMOD        ,0              ,""            ,0           ,1           ,"0: Charge-based, 1: Correlated thermal noise model" )
`MPIcc( NQSMOD         ,0              ,""            ,0           ,2           ,"0: Turn off NQS model; 1: NQS gate resistance (with gi node); 2: NQS charge deficit model from BSIM4 (with q node)" )
`MPIcc( SHMOD          ,0              ,""            ,0           ,1           ,"0: Turn off self-heating; 1: Turn on self-heating" )
`MPIcc( TEMPMOD        ,0              ,""            ,0           ,1           ,"1: Change temperature dependence of specific parameters" )
`MPIcc( RGATEMOD       ,0              ,""            ,0           ,1           ,"0: Turn off gate electrode resistance (without ge node); 1: Turn on gate electrode resistance (with ge node)" )
`MPIcc( RGEOMOD        ,0              ,""            ,0           ,1           ,"Geometry-dependent source/drain resistance; 0: RSH-based; 1: Holistic" )
`MPIcc( CGEOMOD        ,0              ,""            ,0           ,2           ,"Geometry-dependent parasitic capacitance model selector" )
`MPIcc( SH_WARN        ,0              ,""            ,0           ,1           ,"0: Disable self-heating warnings; 1: Enable self-heating warnings" )
`MPIcc( IGCLAMP        ,1              ,""            ,0           ,1           ,"0: Disable gate current clamps; 1: Enable gate current clamps" )
`MPRnb( LINT           ,0.0            ,"m"           ,"Length reduction parameter (dopant diffusion effect)" )
`MPRnb( LL             ,0.0            ,"m^(LLN+1)"   ,"Length reduction parameter (dopant diffusion effect)" )
`MPRnb( LLN            ,1.0            ,""            ,"Length reduction parameter (dopant diffusion effect)" )
`MPRnb( DLC            ,0.0            ,"m"           ,"Delta L for C-V model" )
`MPRnb( DLCACC         ,0.0            ,"m"           ,"Delta L for C-V model in accumulation region (BULKMOD=1 or 2)" )
`MPRnb( DLBIN          ,0.0            ,"m"           ,"Delta L for binning" )
`MPRnb( LLC            ,0.0            ,"m^(LLN+1)"   ,"Length reduction parameter (dopant diffusion effect)" )
`MPRco( EOT            ,1.0e-9         ,"m"           ,1.0e-10     ,inf         ,"Equivalent oxide thickness" )
`MPRco( TOXP           ,1.2e-9         ,"m"           ,1.0e-10     ,inf         ,"Physical oxide thickness" )
`MPRco( EOTBOX         ,1.4e-7         ,"m"           ,1.0e-9      ,inf         ,"Equivalent oxide thickness of the buried oxide (SOI FinFET)" )
`MPRco( HFIN           ,3.0e-8         ,"m"           ,1.0e-9      ,inf         ,"Fin height" )
`MPRcz( FECH           ,1.0            ,""            ,"End-channel factor for different orientation/shape" )
`MPRnb( DELTAW         ,0.0            ,"m"           ,"Change of effective width due to shape of fin/cylinder" )
`MPRcz( FECHCV         ,1.0            ,""            ,"CV end-channel factor for different orientation/shape" )
`MPRnb( DELTAWCV       ,0.0            ,"m"           ,"CV change of effective width due to shape of fin/cylinder" )
`MPRnb( NBODY          ,1.0e22         ,"/m^3"        ,"Channel (body) doping" )
`MPRnb( NBODYN1        ,0.0            ,""            ,"NFIN dependence of channel (body) doping" )
`MPRex( NBODYN2        ,1.0e5          ,""            ,0.0         ,"NFIN dependence of channel (body) doping" )
`MPRcc( NSD            ,2.0e26         ,"/m^3"        ,2.0e25      ,1.0e27      ,"Source/drain active doping concentration" )
`MPRcz( PHIG           ,4.61           ,"eV"          ,"Gate workfunction" )
`MPRnb( PHIGL          ,0.0            ,"eV/m"        ,"Length dependence of gate workfunction" )
`MPRnb( PHIGLT         ,0.0            ,"/m"          ,"Coupled NFIN and length dependence of gate workfunction" )
`MPRnb( PHIGN1         ,0.0            ,""            ,"NFIN dependence of gate workfunction" )
`MPRex( PHIGN2         ,1.0e5          ,""            ,0.0         ,"NFIN dependence of gate workfunction" )
`MPRco( EPSROX         ,3.9            ,""            ,1.0         ,inf         ,"Relative dielectric constant of the gate dielectric" )
`MPRco( EPSRSUB        ,11.9           ,""            ,1.0         ,inf         ,"Relative dielectric constant of the channel material" )
`MPRcz( EASUB          ,4.05           ,"eV"          ,"Electron affinity of substrate" )
`MPRnb( NI0SUB         ,1.1e16         ,"/m^3"        ,"Intrinsic carrier constant at 300.15K" )
`MPRnb( BG0SUB         ,1.12           ,"eV"          ,"Bandgap of substrate at 300.15K" )
`MPRnb( NC0SUB         ,2.86e25        ,"/m^3"        ,"Conduction band density of states" )
`MPRnb( NGATE          ,0.0            ,"/m^3"        ,"Parameter for poly gate doping. For metal gate please set NGATE = 0" )
`MPRnb( Imin           ,1.0e-15        ,"A/m^2"       ,"Parameter for Vgs clamping for inversion region calculation in accumulation" )

// Short Channel Effects
`MPRnb( CIT            ,0.0            ,"F/m^2"       ,"Parameter for interface trap" )
`MPRnb( CITR           ,CIT            ,""            ,"Parameter for interface trap in reverse mode for asymmetric model" )
`MPRnb( CDSC           ,7.0e-3         ,"F/m^2"       ,"Coupling capacitance between S/D and channel" )
`MPRnb( CDSCN1         ,0.0            ,""            ,"NFIN dependence of CDSC" )
`MPRnb( CDSCN2         ,1.0e5          ,""            ,"NFIN dependence of CDSC" )
`MPRnb( CDSCD          ,7.0e-3         ,"F/m^2"       ,"Drain-bias sensitivity of CDSC" )
`MPRnb( CDSCDN1        ,0.0            ,""            ,"NFIN dependence of CDSCD" )
`MPRex( CDSCDN2        ,1.0e5          ,""            ,0.0         ,"NFIN dependence of CDSCD" )
`MPRnb( CDSCDR         ,CDSCD          ,"F/m^2"       ,"Reverse-mode drain-bias sensitivity of CDSC" )
`MPRnb( CDSCDRN1       ,CDSCDN1        ,""            ,"NFIN dependence of CDSCD" )
`MPRex( CDSCDRN2       ,CDSCDN2        ,""            ,0.0         ,"NFIN dependence of CDSCD" )
`MPRnb( DVT0           ,0.0            ,""            ,"SCE coefficient" )
`MPRnb( DVT1           ,0.6            ,""            ,"SCE exponent coefficient. After binning it should be within (0:inf)" )
`MPRnb( DVT1SS         ,DVT1           ,""            ,"Subthreshold swing exponent coefficient. After binning it should be within (0:inf)" )
`MPRnb( PHIN           ,0.05           ,"V"           ,"Nonuniform vertical doping effect on surface potential" )
`MPRnb( ETA0           ,0.6            ,""            ,"DIBL coefficient" )
`MPRnb( ETA0N1         ,0.0            ,""            ,"NFIN dependence of ETA0" )
`MPRco( ETA0N2         ,1.0e5          ,""            ,1.0e-5      ,inf         ,"NFIN dependence of ETA0" )
`MPRnb( ETA0LT         ,0.0            ,"/m"          ,"Coupled NFIN and length dependence of ETA0" )
`MPRnb( TETA0          ,0.0            ,"/K"          ,"Temperature dependence of DIBL coefficient" )
`MPRnb( ETA0R          ,ETA0           ,""            ,"Reverse-mode DIBL coefficient" )
`MPRnb( TETA0R         ,TETA0          ,"/K"          ,"Temperature dependence of reverse-mode DIBL coefficient" )
`MPRnb( DSUB           ,1.06           ,""            ,"DIBL exponent coefficient" )
`MPRnb( DVTP0          ,0.0            ,""            ,"Coefficient for drain-induced Vth shift (DITS)" )
`MPRnb( DVTP1          ,0.0            ,""            ,"DITS exponent coefficient" )
`MPRnb( ADVTP0         ,0.0            ,""            ,"Pre-exponential coefficient for DITS" )
`MPRex( BDVTP0         ,1.0e-7         ,""            ,0.0         ,"Exponential coefficient for DITS" )
`MPRnb( ADVTP1         ,0.0            ,""            ,"Pre-exponential coefficient for DVTP1" )
`MPRex( BDVTP1         ,1.0e-7         ,""            ,0.0         ,"Exponential coefficient for DVTP1" )
`MPRnb( DVTP2          ,0.0            ,""            ,"DITS model parameter" )
`MPRnb( K1RSCE         ,0.0            ,"/V^(0.5)"    ,"K1 for reverse short channel effect calculation" )
`MPRnb( LPE0           ,5.0e-9         ,"m"           ,"Equivalent length of pocket region at zero bias" )
`MPRnb( DVTSHIFT       ,0.0            ,"V"           ,"Vth shift handle" )
`MPRnb( DVTSHIFTR      ,DVTSHIFT       ,""            ,"Vth shift handle for asymmetric mode" )
`MPRnb( THETASCE       ,0.0            ,""            ,"Vth roll-off length dependence. If defined by user, it will overwrite Theta_SCE in the code")
`MPRnb( THETADIBL      ,0.0            ,""            ,"DIBL length dependence. If defined by user, will overwrite Theta_DIBL in the code" )
`MPRnb( THETASW        ,0.0            ,""            ,"Subthreshold swing length dependence. If defined by user, it will overwrite Theta_SW in the code" )
`MPRnb( NVTM           ,0.0            ,"V"           ,"Subthreshold swing factor multiplied by Vtm. If defined by user, it will overwrite nVtm in the code" )

// Lateral non-uniform doping effect (IV-CV Vth shift)
`MPRnb( K0             ,0.0            ,"V"           ,"Lateral NUD voltage parameter" )
`MPRnb( K01            ,0.0            ,"V/K"         ,"Temperature dependence of lateral NUD voltage parameter" )
`MPRnb( K0SI           ,1.0            ,""            ,"Correction factor for strong inversion used in Mnud. After binning it should be within (0:inf)" )
`MPRnb( K0SI1          ,0.0            ,"/K"          ,"Temperature dependence of K0SI" )
`MPRnb( K2SI           ,K0SI           ,""            ,"Correction factor for strong inversion used in Mob" )
`MPRnb( K2SI1          ,K0SI1          ,""            ,"Temperature dependence of K2SI" )
`MPRnb( K0SISAT        ,0.0            ,""            ,"Correction factor for strong inversion used in Mnud" )
`MPRnb( K0SISAT1       ,0.0            ,""            ,"Temperature dependence of K0SISAT" )
`MPRnb( K2SISAT        ,K0SISAT        ,""            ,"Correction factor for strong inversion used in Mob" )
`MPRnb( K2SISAT1       ,K0SISAT1       ,""            ,"Temperature dependence of K2SISAT" )

// Body Effect for MG Devices on Bulk Substrate (ex: FinFETs on BULK)
`MPRnb( PHIBE          ,0.7            ,"V"           ,"Body effect voltage parameter. After binning it should be within [0.2:1.2]" )
`MPRco( K1             ,1.0e-6         ,"V^(0.5)"     ,1.0e-6      ,inf         ,"Body effect coefficient for subthreshold region" )
`MPRnb( K11            ,0.0            ,"V^(0.5)/K"   ,"Temperature dependence of K1" )
`MPRnb( K2SAT          ,0.0            ,""            ,"Correction factor for K2 in saturation (high Vds)" )
`MPRnb( K2SAT1         ,0.0            ,""            ,"Temperature dependence of K2SAT" )
`MPRnb( K2             ,0.0            ,""            ,"Body effect coefficient for BULKMOD==2" )
`MPRnb( K21            ,0.0            ,""            ,"Temperature dependence of K2" )

// Quantum Mechanical Effect
`MPRnb( QMFACTOR       ,0.0            ,""            ,"Prefactor + switch for QM Vth correction" )
`MPRnb( QMTCENCV       ,0.0            ,""            ,"Prefactor + switch for QM Width and Toxeff correction for CV" )
`MPRnb( QMTCENCVA      ,0.0            ,""            ,"Prefactor + switch for QM Width and Toxeff correction for CV (accumulation region)" )
`MPRnb( AQMTCEN        ,0.0            ,""            ,"Parameter for geometric dependence of Tcen on R/TFIN/HFIN" )
`MPRex( BQMTCEN        ,1.2e-8         ,""            ,0.0         ,"Parameter for geometric dependence of Tcen on R/TFIN/HFIN" )
`MPRnb( ETAQM          ,0.54           ,""            ,"Bulk charge coefficient for Tcen" )
`MPRnb( QM0            ,1.0e-3         ,"V"           ,"Knee-point for Tcen in inversion (Charge normalized to Cox)" )
`MPRnb( PQM            ,0.66           ,""            ,"Slope of normalized Tcen in inversion" )
`MPRnb( QM0ACC         ,1.0e-3         ,"V"           ,"Knee-point for Tcen in accumulation (Charge normalized to Cox)" )
`MPRnb( PQMACC         ,0.66           ,""            ,"Slope of normalized Tcen in accumulation" )

// Velocity Saturation Model
`MPRnb( VSAT           ,8.5e4          ,"m/s"         ,"Saturation velocity for the saturation region" )
`MPRnb( VSATR          ,VSAT           ,"m/s"         ,"Saturation velocity for the saturation region in the reverse mode" )
`MPRnb( VSATN1         ,0.0            ,""            ,"NFIN dependence of VSAT" )
`MPRex( VSATN2         ,1.0e5          ,""            ,0.0         ,"NFIN dependence of VSAT" )
`MPRnb( VSATRN1        ,VSATN1         ,""            ,"NFIN dependence of VSATR" )
`MPRex( VSATRN2        ,VSATN2         ,""            ,0.0         ,"NFIN dependence of VSATR" )
`MPRnb( AVSAT          ,0.0            ,""            ,"Pre-exponential coefficient for VSAT" )
`MPRex( BVSAT          ,1.0e-7         ,""            ,0.0         ,"Exponential coefficient for VSAT" )
`MPRnb( VSAT1          ,VSAT           ,"m/s"         ,"Velocity saturation parameter for Ion degradation - forward mode" )
`MPRnb( VSAT1N1        ,VSATN1         ,""            ,"NFIN dependence of VSAT1" )
`MPRex( VSAT1N2        ,VSATN2         ,""            ,0.0         ,"NFIN dependence of VSAT1" )
`MPRnb( VSAT1R         ,VSAT1          ,"m/s"         ,"Velocity saturation parameter for Ion degradation - reverse mode" )
`MPRnb( VSAT1RN1       ,VSAT1N1        ,""            ,"NFIN dependence of VSAT1R" )
`MPRex( VSAT1RN2       ,VSAT1N2        ,""            ,0.0         ,"NFIN dependence of VSAT1R" )
`MPRnb( AVSAT1         ,AVSAT          ,""            ,"Pre-exponential coefficient for VSAT1" )
`MPRex( BVSAT1         ,BVSAT          ,""            ,0.0         ,"Exponential coefficient for VSAT1" )
`MPRnb( DELTAVSAT      ,1.0            ,""            ,"velocity saturation parameter in the linear region" )
`MPRnb( PSAT           ,2.0            ,""            ,"Velocity saturation exponent, after binnig should be from [2.0:inf)" )
`MPRnb( APSAT          ,0.0            ,""            ,"Pre-exponential coefficient for PSAT" )
`MPRex( BPSAT          ,1.0            ,""            ,0.0         ,"Exponential coefficient for PSAT" )
`MPRnb( KSATIV         ,1.0            ,""            ,"Parameter for long channel Vdsat" )
`MPRnb( KSATIVR        ,KSATIV         ,""            ,"KSATIV in asymmetric mode" )
`MPRnb( VSATCV         ,VSAT           ,"m/s"         ,"Velocity saturation parameter for CV" )
`MPRnb( AVSATCV        ,AVSAT          ,""            ,"Pre-exponential coefficient for VSATCV" )
`MPRex( BVSATCV        ,BVSAT          ,""            ,0.0         ,"Exponential coefficient for VSATCV" )
`MPRnb( DELTAVSATCV    ,DELTAVSAT      ,""            ,"Velocity saturation parameter in the linear region for the capacitance model" )
`MPRnb( PSATCV         ,PSAT           ,""            ,"Velocity saturation exponent for C-V" )
`MPRnb( APSATCV        ,APSAT          ,""            ,"Pre-exponential coefficient for PSATCV" )
`MPRex( BPSATCV        ,BPSAT          ,""            ,0.0         ,"Exponential coefficient for PSATCV" )
`MPRnb( MEXP           ,4.0            ,""            ,"Smoothing function factor for Vdsat" )
`MPRnb( AMEXP          ,0.0            ,""            ,"Pre-exponential coefficient for MEXP" )
`MPRnb( BMEXP          ,1.0            ,""            ,"Exponential coefficient for MEXP" )
`MPRnb( MEXPR          ,MEXP           ,""            ,"Reverse-mode smoothing function factor for Vdsat" )
`MPRnb( AMEXPR         ,AMEXP          ,""            ,"Pre-exponential coefficient for MEXPR" )
`MPRnb( BMEXPR         ,BMEXP          ,""            ,"Exponential coefficient for MEXPR" )
`MPRnb( PTWG           ,0.0            ,"/V^2"        ,"Gmsat degradation parameter - forward mode" )
`MPRnb( PTWGR          ,PTWG           ,"/V^2"        ,"Gmsat degradation parameter - reverse mode" )
`MPRnb( APTWG          ,0.0            ,""            ,"Pre-exponential coefficient for PTWG" )
`MPRex( BPTWG          ,1.0e-7         ,""            ,0.0         ,"Exponential coefficient for PTWG" )
`MPRnb( AT             ,-1.56e-3       ,"/K"          ,"Saturation velocity temperature coefficient" )
`MPRnb( ATR            ,AT             ,""            ,"Reverse-mode saturation velocity temperature coefficient" )
`MPRnb( ATCV           ,AT             ,"/K"          ,"Saturation velocity temperature coefficient for CV" )
`MPRnb( TMEXP          ,0.0            ,"/K"          ,"Temperature coefficient for Vdseff smoothing" )
`MPRnb( TMEXPR         ,TMEXP          ,"/K"          ,"Reverse-mode temperature coefficient for Vdseff smoothing" )
`MPRnb( PTWGT          ,4.0e-3         ,"/K"          ,"PTWG temperature coefficient" )

// Mobility Model
`MPRnb( U0             ,3.0e-2         ,"m^2/(V*s)"   ,"Low-field mobility" )
`MPRnb( U0R            ,U0             ,"m^2/(V*s)"   ,"Reverse-mode low-field mobility" )
`MPRnb( U0N1           ,0.0            ,""            ,"NFIN dependence of U0" )
`MPRnb( U0N1R          ,U0N1           ,""            ,"Reverse-mode NFIN dependence of U0" )
`MPRex( U0N2           ,1.0e5          ,""            ,0.0         ,"NFIN dependence of U0" )
`MPRex( U0N2R          ,U0N2           ,""            ,0.0         ,"Reverse-mode NFIN dependence of U0" )
`MPRnb( U0LT           ,0.0            ,"/m"          ,"Coupled NFIN and length dependence of U0" )
`MPRnb( ETAMOB         ,2.0            ,""            ,"Effective field parameter" )
`MPRnb( UP             ,0.0            ,"um^LPA"      ,"Mobility L coefficient" )
`MPRnb( LPA            ,1.0            ,""            ,"Mobility L power coefficient" )
`MPRnb( UPR            ,UP             ,"um^LPA"      ,"Reverse-mode mobility L coefficient" )
`MPRnb( LPAR           ,LPA            ,""            ,"Reverse-mode mobility L power coefficient" )
`MPRnb( UA             ,0.3            ,"(cm/MV)^EU"  ,"Phonon/surface roughness scattering parameter" )
`MPRnb( UAR            ,UA             ,"(cm/MV)^EU"  ,"Reverse-mode phonon/surface roughness scattering parameter" )
`MPRnb( AUA            ,0.0            ,""            ,"Pre-exponential coefficient for UA" )
`MPRnb( AUAR           ,AUA            ,""            ,"Reverse-mode pre-exponential coefficient for UA" )
`MPRex( BUA            ,1.0e-7         ,""            ,0.0         ,"Exponential coefficient for UA" )
`MPRex( BUAR           ,BUA            ,""            ,0.0         ,"Reverse-mode exponential coefficient for UAR" )
`MPRnb( UC             ,0.0            ,"(1e-6*cm/MV^2)^EU"        ,"Body effect for mobility degradation parameter - BULKMOD=1 or 2" )
`MPRnb( UCR            ,UC             ,""            ,"Reverse-mode body effect for mobility degradation parameter - BULKMOD=1 or 2" )
`MPRnb( EU             ,2.5            ,"cm/MV"       ,"Phonon/surface roughness scattering parameter" )
`MPRnb( EUR            ,EU             ,"cm/MV"       ,"Reverse-mode phonon/surface roughness scattering parameter" )
`MPRnb( AEU            ,0.0            ,""            ,"Pre-exponential coefficient for EU" )
`MPRnb( AEUR           ,AEU            ,""            ,"Reverse-mode pre-exponential coefficient for EU" )
`MPRex( BEU            ,1.0e-7         ,""            ,0.0         ,"Exponential coefficient for EU" )
`MPRex( BEUR           ,BEU            ,""            ,0.0         ,"Reverse-mode exponential coefficient for EU" )
`MPRnb( UD             ,0.0            ,"cm/MV"       ,"Columbic scattering parameter" )
`MPRnb( UDR            ,UD             ,"cm/MV"       ,"Reverse-mode columbic scattering parameter" )
`MPRnb( AUD            ,0.0            ,""            ,"Pre-exponential coefficient for UD" )
`MPRnb( AUDR           ,AUD            ,""            ,"Reverse-mode pre-exponential coefficient for UD" )
`MPRex( BUD            ,5.0e-8         ,""            ,0.0         ,"Exponential coefficient for UD" )
`MPRex( BUDR           ,BUD            ,""            ,0.0         ,"Reverse-mode exponential coefficient for UD" )
`MPRnb( UCS            ,1.0            ,""            ,"Columbic scattering parameter" )
`MPRnb( UTE            ,0.0            ,""            ,"Mobility temperature coefficient" )
`MPRnb( UTER           ,UTE            ,""            ,"Reverse-mode for mobility temperature coefficient" )
`MPRnb( UTL            ,-1.5e-3        ,""            ,"Mobility temperature coefficient" )
`MPRnb( UTLR           ,UTL            ,""            ,"Reverse-mode for mobility temperature coefficient" )
`MPRnb( EMOBT          ,0.0            ,""            ,"Temperature coefficient of ETAMOB" )
`MPRnb( UA1            ,1.032e-3       ,""            ,"Mobility temperature coefficient for UA" )
`MPRnb( UA1R           ,UA1            ,""            ,"Reverse-mode mobility temperature coefficient for UA" )
`MPRnb( UC1            ,5.6e-11        ,""            ,"Mobility temperature coefficient for UC" )
`MPRnb( UC1R           ,UC1            ,""            ,"Reverse-mode mobility temperature coefficient for UC" )
`MPRnb( UD1            ,0.0            ,""            ,"Mobility temperature coefficient for UC" )
`MPRnb( UD1R           ,UD1            ,""            ,"Reverse-mode mobility temperature coefficient for UD" )
`MPRnb( UCSTE          ,-4.775e-3      ,""            ,"Mobility temperature coefficient" )
`MPRcc( CHARGEWF       ,0.0            ,""            ,-1.0        ,1.0         ,"Average channel charge weighting factor, +1: source-side, 0: middle, -1: drain-side" )

// Access Resistance Model
`MPRnb( RDSWMIN        ,0.0            ,"ohm*(um^(WR))"            ,"RDSMOD = 0 S/D extension resistance per unit width at high Vgs" )
`MPRnb( RDSW           ,1.0e2          ,"ohm*(um^(WR))"            ,"RDSMOD = 0 zero bias S/D extension resistance per unit width" )
`MPRnb( ARDSW          ,0.0            ,""            ,"Pre-exponential coefficient for RDSW" )
`MPRex( BRDSW          ,1.0e-7         ,""            ,0.0         ,"exponential coefficient for RDSW" )
`MPRnb( RSWMIN         ,0.0            ,"ohm*(um^(WR))"            ,"RDSMOD = 1 source extension resistance per unit width at high Vgs" )
`MPRnb( RSW            ,5.0e1          ,"ohm*(um^(WR))"            ,"RDSMOD = 1 zero bias source extension resistance per unit width" )
`MPRnb( ARSW           ,0.0            ,""            ,"Pre-exponential coefficient for RSW" )
`MPRex( BRSW           ,1.0e-7         ,""            ,0.0         ,"Exponential coefficient for RSW" )
`MPRnb( RDWMIN         ,0.0            ,"ohm*(um^(WR))"            ,"RDSMOD = 1 drain extension resistance per unit width at high Vgs" )
`MPRnb( RDW            ,5.0e1          ,""            ,"RDSMOD = 1 zero bias drain extension resistance per unit width" )
`MPRnb( ARDW           ,0.0            ,""            ,"Pre-exponential coefficient for RDW" )
`MPRex( BRDW           ,1.0e-7         ,""            ,0.0         ,"Exponential coefficient for RDW" )
`MPRcz( RSDR           ,0.0            ,"V^(-PRSDR)"  ,"Source-side drift resistance parameter - forward mode" )
`MPRcz( RSDRR          ,RSDR           ,"V^(-PRSDR)"  ,"Source-side drift resistance parameter - reverse mode" )
`MPRcz( RDDR           ,RSDR           ,"V^(-PRDDR)"  ,"Drain-side drift resistance parameter - forward mode" )
`MPRcz( RDDRR          ,RDDR           ,"V^(-PRDDR)"  ,"Drain-side drift resistance parameter - reverse mode" )
`MPRnb( PRSDR          ,1.0            ,""            ,"Source-side quasi-saturation parameter" )
`MPRnb( PRDDR          ,PRSDR          ,""            ,"Drain-side quasi-saturation parameter" )
`MPRnb( PRWGS          ,0.0            ,"/V"          ,"Gate bias dependence of source extension resistance" )
`MPRnb( PRWGD          ,PRWGS          ,"/V"          ,"Gate bias dependence of drain extension resistance" )
`MPRnb( WR             ,1.0            ,""            ,"W dependence parameter of S/D extension resistance" )
`MPRnb( PRT            ,1.0e-3         ,"/K"          ,"Series resistance temperature coefficient" )
`MPRnb( TRSDR          ,0.0            ,"/K"          ,"Source-side drift resistance temperature coefficient" )
`MPRnb( TRDDR          ,TRSDR          ,"/K"          ,"Drain-side drift resistance temperature coefficient" )

// DIBL Model
`MPRnb( PDIBL1         ,1.3            ,""            ,"DIBL output conductance parameter - forward mode" )
`MPRnb( PDIBL1R        ,PDIBL1         ,""            ,"DIBL output conductance parameter - reverse mode" )
`MPRnb( PDIBL2         ,2.0e-4         ,""            ,"DIBL output conductance parameter" )
`MPRnb( PDIBL2R        ,PDIBL2         ,""            ,"DIBL output conductance parameter - reverse mode" )
`MPRnb( DROUT          ,1.06           ,""            ,"L dependence of DIBL effect on Rout" )
`MPRnb( PVAG           ,1.0            ,""            ,"Vgs dependence on early voltage" )

// Channel Length Modulation Effect
`MPRnb( PCLM           ,1.3e-2         ,""            ,"Channel length modulation (CLM) parameter" )
`MPRnb( PCLMR          ,PCLM           ,""            ,"Reverse model PCLM parameter" )
`MPRnb( APCLM          ,0.0            ,""            ,"Pre-exponential coefficient for PCLM" )
`MPRnb( APCLMR         ,APCLM          ,""            ,"Reverse-mode pre-exponential coefficient for PCLM" )
`MPRex( BPCLM          ,1.0e-7         ,""            ,0.0         ,"Exponential coefficient for PCLM" )
`MPRex( BPCLMR         ,BPCLM          ,""            ,0.0         ,"Reverse-mode exponential coefficient for PCLM" )
`MPRnb( PCLMG          ,0.0            ,""            ,"Gate bias dependence parameter for CLM" )
`MPRnb( PCLMCV         ,PCLM           ,""            ,"CLM parameter for short-channel CV" )

// Non-Saturation Effect
`MPRnb( A1             ,0.0            ,""            ,"Non-saturation effect parameter for strong inversion Region" )
`MPRnb( A11            ,0.0            ,"(V^-2)/K"    ,"Temperature dependence of A1" )
`MPRnb( A2             ,0.0            ,""            ,"Non-saturation effect parameter for moderate Inversion Region" )
`MPRnb( A21            ,0.0            ,"(V^-1)/K"    ,"Temperature dependence of A2" )

// Gate Electrode Resistance
`MPRcz( RGEXT          ,0.0            ,"ohm"         ,"Effective gate electrode external resistance" )
`MPRco( RGFIN          ,1.0e-3         ,"ohm"         ,1.0e-3      ,inf         ,"Effective gate electrode per finger per fin resistance" )

// Geometry Dependent Source/Drain Resistance of RGEOMOD = 0
`MPRnb( RSHS           ,0.0            ,"ohm"         ,"Source-side sheet resistance" )
`MPRnb( RSHD           ,RSHS           ,"ohm"         ,"Drain-side sheet resistance" )

// Geometry Dependent Source/Drain Resistance of RGEOMOD = 1 for variability modeling
// These parameters are shared with CGEOMOD = 2
`MPRnb( HEPI           ,1.0e-8         ,"m"           ,"Height of the raised source/drain on top of the fin" )
`MPRnb( TSILI          ,1.0e-8         ,"m"           ,"Thickness of the silicide on top of the raised source/drain" )
`MPRcc( RHOC           ,1.0e-12        ,"ohm*(m^2)"   ,1.0e-18     ,1.0e-9      ,"Contact resistivity at the silicon/silicide interface" )
`MPRoz( RHORSD         ,1.0            ,"ohm*(m)"     ,"Average resistivity of silicon in the raised source/drain region" )
`MPRcc( CRATIO         ,0.5            ,""            ,0.0         ,1.0         ,"Ratio of the corner area filled with silicon to the total corner area" )
`MPRoo( DELTAPRSD      ,0.0            ,"m"           ,-FPITCH     ,inf         ,"Change in silicon/silicide interface length due to non-rectangular epi" )
`MPIcc( SDTERM         ,0              ,""            ,0           ,1           ,"Indicator of whether the source/drain are terminated with silicide" )
`MPRnb( LSP            ,-1             ,"m"           ,"Thickness of the gate sidewall spacer" )
`MPRco( EPSRSP         ,3.9            ,""            ,1.0         ,inf         ,"Relative dielectric constant of the spacer" )
`MPRoz( TGATE          ,3.0e-8         ,"m"           ,"Gate height on top of the hard mask" )
`MPRcz( TMASK          ,3.0e-8         ,"m"           ,"Height of hard mask on top of the fin" )
`MPRcz( ASILIEND       ,0.0            ,"m^2"         ,"Extra silicide cross sectional area at the two ends of the FinFET" )
`MPRcz( ARSDEND        ,0.0            ,"m^2"         ,"Extra raised source/drain cross sectional areaat the two ends of the FinFET" )
`MPRcz( PRSDEND        ,0.0            ,"m"           ,"Extra silicon/silicide interface perimeter at the two ends of the FinFET" )
`MPRcc( NSDE           ,2.0e25         ,"/(m^3)"      ,1.0e25      ,1.0e26      ,"Source/drain active doping concentration at Leff edge" )
`MPRnb( RGEOA          ,1.0            ,""            ,"Fitting parameter for RGEOMOD=1" )
`MPRnb( RGEOB          ,0.0            ,"/m"          ,"Fitting parameter for RGEOMOD=1" )
`MPRnb( RGEOC          ,0.0            ,"/m"          ,"Fitting parameter for RGEOMOD=1" )
`MPRnb( RGEOD          ,0.0            ,"/m"          ,"Fitting parameter for RGEOMOD=1" )
`MPRnb( RGEOE          ,0.0            ,"/m"          ,"Fitting parameter for RGEOMOD=1" )
`MPRnb( CGEOA          ,1.0            ,""            ,"Fitting parameter for CGEOMOD=2" )
`MPRnb( CGEOB          ,0.0            ,"/m"          ,"Fitting parameter for CGEOMOD=2" )
`MPRnb( CGEOC          ,0.0            ,"/m"          ,"Fitting parameter for CGEOMOD=2" )
`MPRnb( CGEOD          ,0.0            ,"/m"          ,"Fitting parameter for CGEOMOD=2" )
`MPRcz( CGEOE          ,1.0            ,""            ,"Fitting parameter for CGEOMOD=2" )

// Gate Current
`MPRnb( AIGBINV        ,1.11e-2        ,"((F*s^2/g)^0.5)*m^-1"     ,"Parameter for Igb in inversion" )
`MPRnb( AIGBINV1       ,0.0            ,"((F*s^2/g)^0.5)*m^-1/K"   ,"Parameter for Igb in inversion" )
`MPRnb( BIGBINV        ,9.49e-4        ,"((F*s^2/g)^0.5)*(m*V)^-1" ,"Parameter for Igb in inversion" )
`MPRnb( CIGBINV        ,6.0e-3         ,"/V"          ,"Parameter for Igb in inversion" )
`MPRnb( EIGBINV        ,1.1            ,"V"           ,"Parameter for Igb in inversion" )
`MPRnb( NIGBINV        ,3.0            ,""            ,"Parameter for Igb in inversion" )
`MPRnb( AIGBACC        ,1.36e-2        ,"((F*s^2/g)^0.5)*m^-1"     ,"Parameter for Igb in accumulation" )
`MPRnb( AIGBACC1       ,0.0            ,"((F*s^2/g)^0.5)*m^-1/K"   ,"Parameter for Igb in accumulation" )
`MPRnb( BIGBACC        ,1.71e-3        ,"((F*s^2/g)^0.5)*(m*V)^-1" ,"Parameter for Igb in accumulation" )
`MPRnb( CIGBACC        ,7.5e-2         ,"/V"          ,"Parameter for Igb in accumulation" )
`MPRnb( NIGBACC        ,1.0            ,""            ,"Parameter for Igb in accumulation" )
`MPRnb( AIGC           ,1.36e-2        ,"((F*s^2/g)^0.5)*m^-1"     ,"Parameter for Igc in inversion" )
`MPRnb( AIGC1          ,0.0            ,"((F*s^2/g)^0.5)*m^-1/K"   ,"Parameter for Igc in inversion" )
`MPRnb( BIGC           ,1.71e-3        ,"((F*s^2/g)^0.5)*(m*V)^-1" ,"Parameter for Igc in inversion" )
`MPRnb( CIGC           ,7.5e-2         ,"/V"          ,"Parameter for Igc in inversion" )
`MPRnb( PIGCD          ,1.0            ,""            ,"Parameter for Igc partition" )
`MPRnb( DLCIGS         ,0.0            ,"m"           ,"Delta L for Igs model" )
`MPRnb( AIGS           ,1.36e-2        ,"((F*s^2/g)^0.5)*m^-1"     ,"Parameter for Igs in inversion" )
`MPRnb( AIGS1          ,0.0            ,"((F*s^2/g)^0.5)*m^-1/K"   ,"Parameter for Igs in inversion" )
`MPRnb( BIGS           ,1.71e-3        ,"((F*s^2/g)^0.5)*(m*V)^-1" ,"Parameter for Igs in inversion" )
`MPRnb( CIGS           ,7.5e-2         ,"/V"          ,"Parameter for Igs in inversion" )
`MPRnb( DLCIGD         ,DLCIGS         ,"m"           ,"Delta L for Igd model" )
`MPRnb( AIGD           ,AIGS           ,"((F*s^2/g)^0.5)*m^-1"     ,"Parameter for Igd in inversion" )
`MPRnb( AIGD1          ,AIGS1          ,"((F*s^2/g)^0.5)*m^-1/K"   ,"Parameter for Igd in inversion" )
`MPRnb( BIGD           ,BIGS           ,"((F*s^2/g)^0.5)*(m*V)^-1" ,"Parameter for Igd in inversion" )
`MPRnb( CIGD           ,CIGS           ,"/V"          ,"Parameter for Igd in inversion" )
`MPRnb( VFBSD          ,0.0            ,"V"           ,"Flatband voltage for S/D region" )
`MPRnb( VFBSDCV        ,VFBSD          ,"V"           ,"Flatband voltage for S/D region for C-V calculations" )
`MPRoz( TOXREF         ,1.2e-9         ,"m"           ,"Target tox value" )
`MPRnb( TOXG           ,TOXP           ,"m"           ,"Oxide thickness for gate current model" )
`MPRnb( NTOX           ,1.0            ,""            ,"Exponent for Tox ratio" )
`MPRnb( POXEDGE        ,1.0            ,""            ,"Factor for the gate edge Tox" )

// GIDL/GISL Current
`MPRnb( AGISL          ,6.055e-12      ,"mho"         ,"Pre-exponential coefficient for GISL" )
`MPRnb( BGISL          ,3.0e8          ,"V/m"         ,"Exponential coefficient for GISL" )
`MPRnb( CGISL          ,0.5            ,"V^3"         ,"Parameter for body-effect of GISL" )
`MPRnb( EGISL          ,0.2            ,"V"           ,"Band bending parameter for GISL" )
`MPRnb( PGISL          ,1.0            ,""            ,"Parameter for body-bias effect on GISL" )
`MPRnb( AGIDL          ,AGISL          ,"mho"         ,"Pre-exponential coefficient for GIDL" )
`MPRnb( BGIDL          ,BGISL          ,"V/m"         ,"Exponential coefficient for GIDL" )
`MPRnb( CGIDL          ,CGISL          ,"V^3"         ,"Parameter for body-effect of GIDL" )
`MPRnb( EGIDL          ,EGISL          ,"V"           ,"Band bending parameter for GIDL" )
`MPRnb( PGIDL          ,PGISL          ,""            ,"Parameter for body-bias effect on GIDL" )

// Impact Ionization Current
// IIMOD = 1
`MPRnb( ALPHA0         ,0.0            ,"m/V"         ,"First parameter of Iii" )
`MPRnb( ALPHA01        ,0.0            ,"m/V/K"       ,"Temperature dependence of ALPHA0" )
`MPRnb( ALPHA1         ,0.0            ,"/V"          ,"L scaling parameter of Iii" )
`MPRnb( ALPHA11        ,0.0            ,"/V/K"        ,"Temperature dependence ALPHA1" )
`MPRnb( BETA0          ,0.0            ,"/V"          ,"Vds dependence parameter of Iii" )

// IIMOD = 2
`MPRnb( ALPHAII0       ,0.0            ,"m/V"         ,"First parameter of Iii for IIMOD=2" )
`MPRnb( ALPHAII01      ,0.0            ,"m/V/K"       ,"Temperature dependence of ALPHAII0" )
`MPRnb( ALPHAII1       ,0.0            ,"/V"          ,"L scaling parameter of Iii for IIMOD=2" )
`MPRnb( ALPHAII11      ,0.0            ,"m/V/K"       ,"Temperature dependence of ALPHAII1" )
`MPRnb( BETAII0        ,0.0            ,"/V"          ,"Vds dependence parameter of Iii" )
`MPRnb( BETAII1        ,0.0            ,""            ,"Vds dependence parameter of Iii" )
`MPRnb( BETAII2        ,0.1            ,"V"           ,"Vds dependence parameter of Iii" )
`MPRnb( ESATII         ,1.0e7          ,"V/m"         ,"Saturation channel E-field for Iii" )
`MPRnb( LII            ,0.5e-9         ,"V*m"         ,"Channel length dependence parameter of Iii" )
`MPRnb( SII0           ,0.5            ,"/V"          ,"Vgs dependence parameter of Iii" )
`MPRnb( SII1           ,0.1            ,""            ,"1st Vgs dependence parameter of Iii" )
`MPRnb( SII2           ,0.0            ,"V"           ,"2nd Vgs dependence parameter of Iii" )
`MPRnb( SIID           ,0.0            ,"V"           ,"3rd Vds dependence parameter of Iii" )
`MPRoo( IIMOD2CLAMP1   ,0.1            ,"V"           ,0.0         ,inf         ,"Clamp1 of SII1*Vg term in IIMOD=2 model" )
`MPRoo( IIMOD2CLAMP2   ,0.1            ,"V"           ,0.0         ,inf         ,"Clamp2 of SII0*Vg term in IIMOD=2 model" )
`MPRoo( IIMOD2CLAMP3   ,0.1            ,"V"           ,0.0         ,inf         ,"Clamp3 of Ratio term in IIMOD=2 model" )

// Accumulation Capacitance
`MPRco( EOTACC         ,EOT            ,"m"           ,1.0e-10     ,inf         ,"Equivalent oxide thickness for accumulation region" )
`MPRnb( DELVFBACC      ,0.0            ,"V"           ,"Change in flatband voltage: Vfb_accumulation - Vfb_inversion" )

// Fringe Capacitance
// CGEOMOD=0
`MPRcz( CFS            ,2.5e-11        ,"F/m"         ,"Outer fringe capacitance at source side" )
`MPRcz( CFD            ,CFS            ,"F/m"         ,"Outer fringe capacitance at drain side" )

// Overlap Capacitance for CGEOMOD = 0 and 2
`MPRcz( CGSO           ,0.0            ,"F/m"         ,"Non LDD region source-gate overlap capacitance per unit channel width" )
`MPRcz( CGDO           ,CGSO           ,"F/m"         ,"Non LDD region drain-gate overlap capacitance per unit channel width" )
`MPRcz( CGSL           ,0.0            ,"F/m"         ,"Overlap capacitance between gate and lightly-doped source region (for CGEOMOD = 0, 2)" )
`MPRcz( CGDL           ,CGSL           ,"F/m"         ,"Overlap capacitance between gate and lightly-doped drain region (for CGEOMOD = 0, 2)" )
`MPRco( CKAPPAS        ,0.6            ,"V"           ,2.0e-2      ,inf         ,"Coefficient of bias-dependent overlap capacitance for the source side (for CGEOMOD = 0, 2)" )
`MPRco( CKAPPAD        ,CKAPPAS        ,"V"           ,2.0e-2      ,inf         ,"Coefficient of bias-dependent overlap capacitance for the drain side (for CGEOMOD = 0, 2)" )
`MPRcz( CGBO           ,0.0            ,"F/m"         ,"Gate-to-substrate overlap capacitance per unit channel length per finger per NGCON" )
`MPRcz( CGBN           ,0.0            ,"F/m"         ,"Gate-to-substrate overlap capacitance per unit channel length per fin per finger" )
`MPRcz( CGBL           ,0.0            ,"F/m"         ,"Bias dependent component of gate-to-substrate overlap capacitance per unit channel length per fin per finger" )
`MPRco( CKAPPAB        ,0.6            ,""            ,2.0e-2      ,inf         ,"Bias dependent gate-to-substrate parasitic capacitance" )

// Source/Drain-to-Substrate Sidewall Capacitance
`MPRcz( CSDESW         ,0.0            ,"F/m"         ,"Coefficient for source/drain-to-substrate sidewall capacitance" )

// Junction Current and Capacitance
// Junction Capacitance
`MPRnb( CJS            ,5.0e-4         ,"F/m^2"       ,"Unit area source-side junction capacitance at zero bias" )
`MPRnb( CJD            ,CJS            ,"F/m^2"       ,"Unit area drain-side junction capacitance at zero bias" )
`MPRnb( CJSWS          ,5.0e-10        ,"F/m"         ,"Unit length source-side sidewall junction capacitance at zero bias" )
`MPRnb( CJSWD          ,CJSWS          ,"F/m"         ,"Unit length drain-side sidewall junction capacitance at zero bias" )
`MPRnb( CJSWGS         ,0.0            ,"F/m"         ,"Unit length source-side gate sidewall junction capacitance at zero bias" )
`MPRnb( CJSWGD         ,CJSWGS         ,"F/m"         ,"Unit length drain-side gate sidewall junction capacitance at zero bias" )
`MPRnb( PBS            ,1.0            ,"V"           ,"Source-side bulk junction built-in potential" )
`MPRnb( PBD            ,PBS            ,"V"           ,"Drain-side bulk junction built-in potential" )
`MPRnb( PBSWS          ,1.0            ,"V"           ,"Built-in potential for Source-side sidewall junction capacitance" )
`MPRnb( PBSWD          ,PBSWS          ,"V"           ,"Built-in potential for Drain-side sidewall junction capacitance" )
`MPRnb( PBSWGS         ,PBSWS          ,"V"           ,"Built-in potential for Source-side gate sidewall junction capacitance" )
`MPRnb( PBSWGD         ,PBSWGS         ,"V"           ,"Built-in potential for Drain-side gate sidewall junction capacitance" )
`MPRnb( MJS            ,0.5            ,""            ,"Source bottom junction capacitance grading coefficient" )
`MPRnb( MJD            ,MJS            ,""            ,"Drain bottom junction capacitance grading coefficient" )
`MPRnb( MJSWS          ,0.33           ,""            ,"Source sidewall junction capacitance grading coefficient" )
`MPRnb( MJSWD          ,MJSWS          ,""            ,"Drain sidewall junction capacitance grading coefficient" )
`MPRnb( MJSWGS         ,MJSWS          ,""            ,"Source-side gate sidewall junction capacitance grading coefficient" )
`MPRnb( MJSWGD         ,MJSWGS         ,""            ,"Drain-side gate sidewall junction capacitance grading coefficient" )

// Second Junction for Two-Step Junction Capacitance
`MPRcz( SJS            ,0.0            ,""            ,"Constant for source-side two-step second junction" )
`MPRcz( SJD            ,SJS            ,""            ,"Constant for drain-side two-step second junction" )
`MPRcz( SJSWS          ,0.0            ,""            ,"Constant for source-side sidewall two-step second junction" )
`MPRcz( SJSWD          ,SJSWS          ,""            ,"Constant for drain-side sidewall two-step second junction" )
`MPRcz( SJSWGS         ,0.0            ,""            ,"Constant for source-side gate sidewall two-step second junction" )
`MPRcz( SJSWGD         ,SJSWGS         ,""            ,"Constant for source-side gate sidewall two-step second junction" )
`MPRnb( MJS2           ,1.25e-1        ,""            ,"Source bottom two-step second junction capacitance grading coefficient" )
`MPRnb( MJD2           ,MJS2           ,""            ,"Drain bottom two-step second junction capacitance grading coefficient" )
`MPRnb( MJSWS2         ,8.3e-2         ,""            ,"Source sidewall two-step second junction capacitance grading coefficient" )
`MPRnb( MJSWD2         ,MJSWS2         ,""            ,"Drain sidewall two-step second junction capacitance grading coefficient" )
`MPRnb( MJSWGS2        ,MJSWS2         ,""            ,"Source-side gate sidewall two-step second junction capacitance grading coefficient" )
`MPRnb( MJSWGD2        ,MJSWGS2        ,""            ,"Drain-side gate sidewall two-step second junction capacitance grading coefficient" )

// Junction Current
`MPRnb( JSS            ,1.0e-4         ,"A/m^2"       ,"Bottom source junction reverse saturation current density" )
`MPRnb( JSD            ,JSS            ,"A/m^2"       ,"Bottom drain junction reverse saturation current density" )
`MPRnb( JSWS           ,0.0            ,"A/m"         ,"Unit length reverse saturation current for sidewall source junction" )
`MPRnb( JSWD           ,JSWS           ,"A/m"         ,"Unit length reverse saturation current for sidewall drain junction" )
`MPRnb( JSWGS          ,0.0            ,"A/m"         ,"Unit length reverse saturation current for gate-edge sidewall source junction" )
`MPRnb( JSWGD          ,JSWGS          ,"A/m"         ,"Unit length reverse saturation current for gate-edge sidewall drain junction" )
`MPRex( NJS            ,1.0            ,""            ,0.0         ,"Source junction emission coefficient" )
`MPRex( NJD            ,NJS            ,""            ,0.0         ,"Drain junction emission coefficient" )
`MPRnb( IJTHSFWD       ,0.1            ,"A"           ,"Forward source diode breakdown limiting current" )
`MPRnb( IJTHDFWD       ,IJTHSFWD       ,"A"           ,"Forward drain diode breakdown limiting current" )
`MPRnb( IJTHSREV       ,0.1            ,"A"           ,"Reverse source diode breakdown limiting current" )
`MPRnb( IJTHDREV       ,IJTHSREV       ,"A"           ,"Reverse drain diode breakdown limiting current" )
`MPRnb( BVS            ,1.0e1          ,"V"           ,"Source diode breakdown voltage" )
`MPRnb( BVD            ,BVS            ,"V"           ,"Drain diode breakdown voltage" )
`MPRnb( XJBVS          ,1.0            ,""            ,"Fitting parameter for source diode breakdown current" )
`MPRnb( XJBVD          ,XJBVS          ,""            ,"Fitting parameter for drain diode breakdown current" )

// Tunneling Component of Junction Current
`MPRnb( JTSS           ,0.0            ,"A/m^2"       ,"Bottom source junction trap-assisted saturation current density" )
`MPRnb( JTSD           ,JTSS           ,"A/m^2"       ,"Bottom drain junction trap-assisted saturation current density" )
`MPRnb( JTSSWS         ,0.0            ,"A/m"         ,"Unit length trap-assisted saturation current for sidewall source junction" )
`MPRnb( JTSSWD         ,JTSSWS         ,"A/m"         ,"Unit length trap-assisted saturation current for sidewall drain junction" )
`MPRnb( JTSSWGS        ,0.0            ,"A/m"         ,"Unit length trap-assisted saturation current for gate-edge sidewall source junction" )
`MPRnb( JTSSWGD        ,JTSSWGS        ,"A/m"         ,"Unit length trap-assisted saturation current for gate-edge sidewall drain junction" )
`MPRnb( JTWEFF         ,0.0            ,"m"           ,"Trap-assisted tunneling current width dependence" )
`MPRnb( NJTS           ,2.0e1          ,""            ,"Non-ideality factor for JTSS" )
`MPRnb( NJTSD          ,NJTS           ,""            ,"Non-ideality factor for JTSD" )
`MPRnb( NJTSSW         ,2.0e1          ,""            ,"Non-ideality factor for JTSSWS" )
`MPRnb( NJTSSWD        ,NJTSSW         ,""            ,"Non-ideality factor for JTSSWD" )
`MPRnb( NJTSSWG        ,2.0e1          ,""            ,"Non-ideality factor for JTSSWGS" )
`MPRnb( NJTSSWGD       ,NJTSSWG        ,""            ,"Non-ideality factor for JTSSWGD" )
`MPRnb( VTSS           ,1.0e1          ,"V"           ,"Bottom source junction trap-assisted current voltage dependent parameter" )
`MPRnb( VTSD           ,VTSS           ,"V"           ,"Bottom drain junction trap-assisted current voltage dependent parameter" )
`MPRnb( VTSSWS         ,1.0e1          ,"V"           ,"Unit length trap-assisted current voltage dependent parameter for sidewall source junction" )
`MPRnb( VTSSWD         ,VTSSWS         ,"V"           ,"Unit length trap-assisted current voltage dependent parameter for sidewall drain junction" )
`MPRnb( VTSSWGS        ,1.0e1          ,"V"           ,"Unit length trap-assisted current voltage dependent parameter for gate-edge sidewall source junction" )
`MPRnb( VTSSWGD        ,VTSSWGS        ,"V"           ,"Unit length trap-assisted current voltage dependent parameter for gate-edge sidewall drain junction" )

// Recombination-Generation Current
`MPRnb( LINTIGEN       ,0.0            ,"m"           ,"Lint for thermal generation current" )
`MPRnb( NTGEN          ,1.0            ,""            ,"Thermal generation current parameter" )
`MPRnb( AIGEN          ,0.0            ,"(m^-3)*(V^-1)"            ,"Thermal generation current parameter" )
`MPRnb( BIGEN          ,0.0            ,"(m^-3)*(V^-3)"            ,"Thermal generation current parameter" )

// NQS Gate Resistance Model & NQS Charge Deficit Model
// For NQSMOD=1, Set XRCRG1=0 to turn off NQS gate resistance
`MPRnb( XRCRG1         ,1.2e1          ,""            ,"Parameter for non-quasistatic gate resistance (NQSMOD = 1) and NQSMOD = 2" )
`MPRnb( XRCRG2         ,1.0            ,""            ,"Parameter for non-quasistatic gate resistance (NQSMOD = 1) and NQSMOD = 2" )

// NQS Charge Segmentation Model
`MPIcc( NSEG           ,4              ,""            ,4           ,10          ,"Number of segments for NQSMOD=3 (3, 5 & 10 supported)" )

// Flicker Noise
`MPRnb( EF             ,1.0            ,""            ,"Flicker noise frequency exponent" )
`MPRnb( EM             ,4.1e7          ,"V/m"         ,"Flicker noise parameter" )
`MPRnb( NOIA           ,6.25e39        ,"(eV^-1)*(s^(1-EF))*(m^-3)"             ,"Flicker noise parameter" )
`MPRnb( NOIB           ,3.125e24       ,"(eV^-1)*(s^(1-EF))*(m^-1)"             ,"Flicker noise parameter" )
`MPRnb( NOIC           ,8.75e7         ,"(eV^-1)*(s^(1-EF))*(m)"                ,"Flicker noise parameter" )
`MPRnb( LINTNOI        ,0.0            ,"m^2"         ,"L offset for flicker noise calculation" )

// Thermal Noise
`MPRcz( NTNOI          ,1.0            ,""            ,"Thermal noise parameter" )
`MPRnb( TNOIA          ,1.5            ,"/m"          ,"Thermal noise parameter" )
`MPRnb( TNOIB          ,3.5            ,"/m"          ,"Thermal noise parameter" )
`MPRnb( RNOIA          ,5.77e-1        ,""            ,"Thermal noise coefficient" )
`MPRnb( RNOIB          ,3.7e-1         ,""            ,"Thermal noise coefficient" )

// Parameters Controlled by Correlated Thermal Noise Switch
`ifdef  __TNOIMOD1__
    `MPRnb( TNOIC          ,3.5            ,""            ,"Thermal noise parameter for TNOIMOD=1" )
    `MPRnb( RNOIC          ,3.95e-1        ,""            ,"Thermal noise coefficient for TNOIMOD=1" )
    `MPRex( SCALEN         ,1.0e5          ,""            ,0.0         ,"Noise scaling parameter for TNOIMOD=1" )
`endif

// Temperature Effects
`MPRco( TNOM           ,27.0           ,"Celsius"     ,-`P_CELSIUS0,inf         ,"Temperature at which the model is extracted" )
`MPRnb( TBGASUB        ,7.02e-4        ,"eV/K"        ,"Bandgap temperature coefficient" )
`MPRnb( TBGBSUB        ,1.108e3        ,"K"           ,"Bandgap temperature coefficient" )
`MPRnb( KT1            ,0.0            ,"V"           ,"Vth temperature coefficient" )
`MPRnb( KT1L           ,0.0            ,"V*m"         ,"Vth temperature L coefficient" )
`MPRnb( TSS            ,0.0            ,"/K"          ,"Swing temperature coefficient" )
`MPRnb( IIT            ,-0.5           ,""            ,"Impact ionization temperature dependence for IIMOD = 1" )
`MPRnb( TII            ,0.0            ,""            ,"Impact ionization temperature dependence for IIMOD = 2" )
`MPRnb( TGIDL          ,-3.0e-3        ,"/K"          ,"GIDL/GISL temperature dependence" )
`MPRnb( IGT            ,2.5            ,""            ,"Gate current temperature dependence" )
`MPRnb( TCJ            ,0.0            ,"/K"          ,"Temperature coefficient for CJS/CJD" )
`MPRnb( TCJSW          ,0.0            ,"/K"          ,"Temperature coefficient for CJSWS/CJSWD" )
`MPRnb( TCJSWG         ,0.0            ,"/K"          ,"Temperature coefficient for CJSWGS/CJSWGD" )
`MPRnb( TPB            ,0.0            ,"/K"          ,"Temperature coefficient for PBS/PBD" )
`MPRnb( TPBSW          ,0.0            ,"/K"          ,"Temperature coefficient for PBSWS/PBSWD" )
`MPRnb( TPBSWG         ,0.0            ,"/K"          ,"Temperature coefficient for PBSWGS/PBSWGD" )
`MPRnb( XTIS           ,3.0            ,""            ,"Source junction current temperature exponent" )
`MPRnb( XTID           ,XTIS           ,""            ,"Drain junction current temperature exponent" )
`MPRnb( XTSS           ,2.0e-2         ,""            ,"Power dependence of JTSS on temperature" )
`MPRnb( XTSD           ,XTSS           ,""            ,"Power dependence of JTSD on temperature" )
`MPRnb( XTSSWS         ,2.0e-2         ,""            ,"Power dependence of JTSSWS on temperature" )
`MPRnb( XTSSWD         ,XTSSWS         ,""            ,"Power dependence of JTSSWD on temperature" )
`MPRnb( XTSSWGS        ,2.0e-2         ,""            ,"Power dependence of JTSSWGS on temperature" )
`MPRnb( XTSSWGD        ,XTSSWGS        ,""            ,"Power dependence of JTSSWGD on temperature" )
`MPRnb( TNJTS          ,0.0            ,""            ,"Temperature coefficient for NJTS" )
`MPRnb( TNJTSD         ,TNJTS          ,""            ,"Temperature coefficient for NJTSD" )
`MPRnb( TNJTSSW        ,0.0            ,""            ,"Temperature coefficient for NJTSSW" )
`MPRnb( TNJTSSWD       ,TNJTSSW        ,""            ,"Temperature coefficient for NJTSSWD" )
`MPRnb( TNJTSSWG       ,0.0            ,""            ,"Temperature coefficient for NJTSSWG" )
`MPRnb( TNJTSSWGD      ,TNJTSSWG       ,""            ,"Temperature coefficient for NJTSSWGD" )

// Self Heating
`MPRcz( RTH0           ,1.0e-2         ,"ohm*m*K/W"   ,"Thermal resistance" )
`MPRcz( CTH0           ,1.0e-5         ,"W*s/m/K"     ,"Thermal capacitance" )
`MPRcz( WTH0           ,0.0            ,"m"           ,"Width dependence coefficient for Rth and Cth" )
`MPRcz( ASHEXP         ,1.0            ,""            ,"Exponent to tune RTH dependence of NFINTOTAL" )
`MPRcz( BSHEXP         ,1.0            ,""            ,"Exponent to tune RTH dependence of NF" )

// Unified Model
`MPRoz( ACH_UFCM       ,1.0            ,"m^2"         ,"Area of the channel for the unified Model" )
`MPRoz( CINS_UFCM      ,1.0            ,"F/m"         ,"Insulator capacitance for the unified Model" )
`MPRoz( W_UFCM         ,1.0            ,"m"           ,"Effective channel width for the unified Model" )
`MPRcz( TFIN_TOP       ,1.5e-8         ,"m"           ,"Top body (fin) thickness for trapezoidal triple gate" )
`MPRco( TFIN_BASE      ,1.5e-8         ,"m"           ,1.0e-9      ,inf         ,"Base body (fin) thickness for trapezoidal triple gate" )
`MPRcz( QMFACTORCV     ,0.0            ,""            ,"Charge dependence taking QM effects into account" )
`MPRcz( ALPHA_UFCM     ,0.5556         ,""            ,"Mobile charge scaling term taking QM effects into account" )

// Binning Parameters
`include "bsimcmg_binning_parameters.include"

// Output Variables
`ifdef __OPINFO__
    (* desc=  "WEFF" *)        real WEFF;
    (* desc=  "LEFF" *)        real LEFF;
    (* desc=  "WEFFCV" *)      real WEFFCV;
    (* desc=  "LEFFCV" *)      real LEFFCV;
    (* desc=  "IDS" *)         real IDS;
    (* desc=  "IDEFF" *)       real IDEFF;
    (* desc=  "ISEFF" *)       real ISEFF;
    (* desc=  "IGTOT" *)       real IGTOT;
    (* desc=  "IDSGEN" *)      real IDSGEN;
    (* desc=  "III" *)         real III;
    (* desc=  "IGS" *)         real IGS;
    (* desc=  "IGD" *)         real IGD;
    (* desc=  "IGCS" *)        real IGCS;
    (* desc=  "IGCD" *)        real IGCD;
    (* desc=  "IGBS" *)        real IGBS;
    (* desc=  "IGBD" *)        real IGBD;
    (* desc=  "IGIDL" *)       real IGIDL;
    (* desc=  "IGISL" *)       real IGISL;
    (* desc=  "IJSB" *)        real IJSB;
    (* desc=  "IJDB" *)        real IJDB;
    (* desc=  "ISUB" *)        real ISUB;
    (* desc=  "BETA" *)        real BETA;
    (* desc=  "VTH" *)         real VTH;
    (* desc=  "VDSSAT" *)      real VDSSAT;
    (* desc=  "VFB" *)         real VFB;
    (* desc=  "GM" *)          real GM;
    (* desc=  "GDS" *)         real GDS;
    (* desc=  "GMBS" *)        real GMBS;
    (* desc=  "QGI" *)         real QGI;
    (* desc=  "QDI" *)         real QDI;
    (* desc=  "QSI" *)         real QSI;
    (* desc=  "QBI" *)         real QBI;
    (* desc=  "QG" *)          real QG;
    (* desc=  "QD" *)          real QD;
    (* desc=  "QS" *)          real QS;
    (* desc=  "QB" *)          real QB;
    (* desc=  "CGGI" *)        real CGGI;
    (* desc=  "CGSI" *)        real CGSI;
    (* desc=  "CGDI" *)        real CGDI;
    (* desc=  "CGEI" *)        real CGEI;
    (* desc=  "CDGI" *)        real CDGI;
    (* desc=  "CDDI" *)        real CDDI;
    (* desc=  "CDSI" *)        real CDSI;
    (* desc=  "CDEI" *)        real CDEI;
    (* desc=  "CSGI" *)        real CSGI;
    (* desc=  "CSDI" *)        real CSDI;
    (* desc=  "CSSI" *)        real CSSI;
    (* desc=  "CSEI" *)        real CSEI;
    (* desc=  "CEGI" *)        real CEGI;
    (* desc=  "CEDI" *)        real CEDI;
    (* desc=  "CESI" *)        real CESI;
    (* desc=  "CEEI" *)        real CEEI;
    (* desc=  "CGG" *)         real CGG;
    (* desc=  "CGS" *)         real CGS;
    (* desc=  "CGD" *)         real CGD;
    (* desc=  "CGE" *)         real CGE;
    (* desc=  "CDG" *)         real CDG;
    (* desc=  "CDD" *)         real CDD;
    (* desc=  "CDS" *)         real CDS;
    (* desc=  "CDE" *)         real CDE;
    (* desc=  "CSG" *)         real CSG;
    (* desc=  "CSD" *)         real CSD;
    (* desc=  "CSS" *)         real CSS;
    (* desc=  "CSE" *)         real CSE;
    (* desc=  "CEG" *)         real CEG;
    (* desc=  "CED" *)         real CED;
    (* desc=  "CES" *)         real CES;
    (* desc=  "CEE" *)         real CEE;
    (* desc=  "CGSEXT" *)      real CGSEXT;
    (* desc=  "CGDEXT" *)      real CGDEXT;
    (* desc=  "CGBOV" *)       real CGBOV;
    (* desc=  "CJST" *)        real CJST;
    (* desc=  "CJDT" *)        real CJDT;
    (* desc=  "RSGEO" *)       real RSGEO;
    (* desc=  "RDGEO" *)       real RDGEO;
    (* desc=  "CFGEO" *)       real CFGEO;
    (* desc=  "T_TOTAL_K" *)   real T_TOTAL_K;
    (* desc=  "T_TOTAL_C" *)   real T_TOTAL_C;
    (* desc=  "T_DELTA_SH" *)  real T_DELTA_SH;

    `ifdef __DEBUG__
        (* desc=  "IGBACC" *)    real IGBACC;
        (* desc=  "IGBINV" *)    real IGBINV;
        (* desc=  "DIDSDVG" *)   real DIDSDVG;
        (* desc=  "DIDSDVS" *)   real DIDSDVS;
        (* desc=  "DIDSDVD" *)   real DIDSDVD;
        (* desc=  "DIGSDVG" *)   real DIGSDVG;
        (* desc=  "DIGSDVS" *)   real DIGSDVS;
        (* desc=  "DIGSDVD" *)   real DIGSDVD;
        (* desc=  "DIGDDVG" *)   real DIGDDVG;
        (* desc=  "DIGDDVS" *)   real DIGDDVS;
        (* desc=  "DIGDDVD" *)   real DIGDDVD;
        (* desc=  "DIIIDVG" *)   real DIIIDVG;
        (* desc=  "DIIIDVS" *)   real DIIIDVS;
        (* desc=  "DIIIDVD" *)   real DIIIDVD;
        (* desc=  "DIGIDLDVG" *) real DIGIDLDVG;
        (* desc=  "DIGIDLDVS" *) real DIGIDLDVS;
        (* desc=  "DIGIDLDVD" *) real DIGIDLDVD;
        (* desc=  "DIGISLDVG" *) real DIGISLDVG;
        (* desc=  "DIGISLDVS" *) real DIGISLDVS;
        (* desc=  "DIGISLDVD" *) real DIGISLDVD;

        `ifdef __SHMOD__
            (* desc=  "CGT" *)        real CGT;
            (* desc=  "CST" *)        real CST;
            (* desc=  "CDT" *)        real CDT;
            (* desc=  "DIDSDVTH" *)   real DIDSDVTH;
            (* desc=  "DIGSDVTH" *)   real DIGSDVTH;
            (* desc=  "DIGDDVTH" *)   real DIGDDVTH;
            (* desc=  "DIIIDVTH" *)   real DIIIDVTH;
            (* desc=  "DIGIDLDVTH" *) real DIGIDLDVTH;
            (* desc=  "DIGISLDVTH" *) real DIGISLDVTH;
            (* desc=  "DITHDVTH" *)   real DITHDVTH;
        `endif

        (* desc=  "ITH" *)      real ITH;
        (* desc=  "DITHDVG" *)  real DITHDVG;
        (* desc=  "DITHDVS" *)  real DITHDVS;
        (* desc=  "DITHDVD" *)  real DITHDVD;
    `endif
`endif

// Variables Inside the Model
integer devsign;

real NFINtotal;
real DevTemp;
real ids0, ids0_ov_dqi, ids, vgs, vds, vdsx, sigvds, vch, etaiv;
real vgs_noswap, vds_noswap, vgd_noswap;
real qd, qg, qs, qb;
real ni, epssub, epssp, epsratio, Eg, Eg0, Nc;
real Lg, deltaL, deltaL1, deltaLCV, Leff, Leff1, LeffCV, LeffCV_acc, Weff0, WeffCV0;
real cox, cdsc, cbox;
real nbody, phib, deltaPhi;
real T0, T0y, T1, T1y, T2, T2y, T3, T3y, T4, T4a, T5, T6, T7, T8, T9;
real Vtm, Vtm0, nVtm;
real beta, beta0 ;
real wf, wr;

// Temperature Effects
real Tnom, TRatio, dvth_temp, delTemp, ThetaSS;
real K0_t, K0SI_t, K2SI_t, K1_t, K2SAT_t, A1_t, A2_t;
real AIGBINV_t, AIGBACC_t, AIGC_t, AIGS_t, AIGD_t;
real BETA0_t, SII0_t, BGISL_t, BGIDL_t, igtemp, PTWG_t, PTWGR_t;
real ALPHA0_t, ALPHA1_t, ALPHAII0_t, ALPHAII1_t;
real CJS_t, CJSWS_t, CJSWGD_t, CJD_t, CJSWD_t, CJSWGS_t;
real PBS_t, PBSWS_t, PBSWGS_t, PBD_t, PBSWD_t, PBSWGD_t;
real JSS_t, JSWS_t, JSWGS_t, JSD_t, JSWD_t, JSWGD_t;
real JTSS_t, JTSD_t, JTSSWS_t, JTSSWD_t, JTSSWGS_t, JTSSWGD_t;
real NJTS_t, NJTSD_t, NJTSSW_t, NJTSSWD_t, NJTSSWG_t, NJTSSWGD_t;
real K2_t;
real K0SISAT_t, K2SISAT_t;

// Variables for analytical surface potential
real q0;
real T10, T11, T12;
real e0, e1, e2;

// Accumulation Model
real vgsfb, vgsfbeff;

// Short Channel Effect
real ETA0_t, ETA0R_t;
real scl, vbi, phist, dvth_vtroll, dvth_dibl, dvth_rsce, dvth_all;
real tmp, Theta_SCE, Theta_SW, Theta_DIBL, Theta_RSCE, Theta_DITS;

// Lateral Non-uniform Doping Effect
real Mnud;

// Body Effect for BULKMOD=1
real ves, vesx, vesmax, veseff;
real Mob;

// Quantum mechanical correction [units are MKS]
real coxeff, Tcen0, Tcen, dvch_qm, MTcen;
real E0, E0prime, E1, E1prime, mx, mxprime, md, mdprime;
real gprime, gfactor, gam0, gam1, kT;

// Drain Saturation Voltage
real Vdseff, qis, qid, qbs, Dmobs;

// Midpoint Potential and Charge
real qia, qia2, qba, dqi;
real qb0;
real eta_mu, eta_mu_cv, Eeffm, Eeffm_cv, Dmob, Dmob_cv, u0, ueff, u0_a, u0r;
real UA_t, UAR_t, UC_t, UCR_t, UCS_t, UD_t, UDR_t, U0_t, U0R_t, ETAMOB_t, Eeffs, EeffFactor;

real Dr, WeffWRFactor;
real RSourceGeo, RDrainGeo;
real RDSWMIN_i, RDWMIN_i, RSWMIN_i;
real Rdrain, Rsource;

real rdstemp, Rdsi, Rdss;
real RSDR_t, RSDRR_t, RDDR_t, RDDRR_t;

real DIBLfactor, PVAGfactor, diffVds, VaDIBL, Vgst2Vtm, Moc, Mclm;
real MclmCV, inv_MclmCV;

real Dvsat, Vdsat, inv_MEXP, DvsatCV, Nsat;
real VSAT_t, VSAT1_t, VSAT1R_t, VSATCV_t, MEXP_t, MEXPR_t, Esat, EsatL, Esat1, Esat1L, EsatCV, EsatCVL;
real WVCox, Ta, Tb, Tc;

// Asymmetry Model
real VSAT1_a, MEXP_a, PTWG_a, RSDR_a, RDDR_a, PDIBL1_a, VSAT_a;

// Geometry dependent Source/Drain Resistance
real mu_max, mu_rsd, rhorsd, afin, thetarsp;
real Rsp, lt, arsd_total, prsd_total, alpha;
real eta, RrsdTML, Rrsdside, Rrsd;
real Rdsgeo, Arsd, Prsd;

// Geometry dependent fringing capacitance
real Hg, Wg, Trsd, Hrsd, Cgg_top, Cgg_side, Cfr_geo, Acorner, Ccorner;

// Gate Electrode Resistance
`ifdef __RGATEMOD__
    real ggeltd, Rgeltd;
`endif

// Gate Current
real Vaux_Igbinv, igbinv, igsd_mult, igsd_mult0, igbs, igbd;
real Voxacc, Vaux_Igbacc, vfbzb, igbacc;
real igcs, igcd, igc0, Vdseffx, T1_exp;
real igisl, igidl, vfbsd, igs, igd, vgs_eff, vgd_eff;
real Aechvb, Bechvb, Toxratio, Toxratioedge;

// Impact Ionization current
real Iii, Vdiff, Vdsatii, VgsStep, Ratio, ALPHAII;

// Accumulation Capacitance
real cox_acc;
real qg_acc, qb_acc;
real vge;

// Parasitic Capacitance
real qgs_ov, qgd_ov, qgs_fr, qgd_fr, qds_fr;
real qgs_parasitic, qgd_parasitic, Qes, Qed, Qeg;
real vgs_overlap, vgd_overlap, vge_overlap;
real cgsp, cgdp, csbox, cdbox, cgbox, vfbsdcv;

// Junction Current and Capacitance
real Ies, Ied, ves_jct, ved_jct, vec;
real Czbs, Czbssw, Czbsswg, Czbd, Czbdsw, Czbdswg;
real pb2, arg, sarg, Qec;
real Qesj, Qesj1, Qesj2, Qesj3, Qedj, Qedj1, Qedj2, Qedj3;
real Isbs, Isbd, Nvtms, Nvtmd;
real SslpRev, IVjsmRev, VjsmRev, SslpFwd, IVjsmFwd, VjsmFwd, XExpBVS;
real DslpRev, IVjdmRev, VjdmRev, DslpFwd, IVjdmFwd, VjdmFwd, XExpBVD;
real igentemp, idsgen, LINTIGEN_i;

// NQS Gate Resistance
`ifdef __NQSMOD1__
    real gcrg, XRCRG1_i, XRCRG2_i;
    real IdovVds;
`endif

// NQS Charge Deficit Model
`ifdef __NQSMOD2__
    real xdpart, gtau, gcrg, XRCRG1_i, XRCRG2_i;
    real IdovVds;
`endif

// Flicker Noise
real LINTNOI_i;
real litl, Esatnoi, Leffnoi, Leffnoisq, DelClm;
real N0, Nl, Nstar, Ssi, Swi, FNPowerAt1Hz;

// Thermal Noise
real NTNOI_i, qinv;
real Gtnoi, sid;
real gspr, gdpr;

// Variables Controlled by Correlated Thermal Noise Switch
`ifdef  __TNOIMOD1__
    real Abulk, etaa, gamma, delta, epsilon, gche;
    real npart_beta, npart_theta, ctnoi, npart_c;
    real noiGd0, GammaGd0, C0, sf;
`endif

// Self Heating
`ifdef __SHMOD__
    real gth, cth;
`endif

// Binning
real Inv_L, Inv_NFIN, Inv_LNFIN;
real NBODY_i, PHIG_i, CFD_i, CFS_i, COVS_i, COVD_i, CGSO_i, CGDO_i;
real CGSL_i, CGDL_i, CGBL_i, CKAPPAS_i, CKAPPAD_i, CKAPPAB_i;
real QMFACTOR_i, QMTCENCV_i, QMTCENCVA_i, KSATIV_i, KSATIVR_i, KSATIV_a;
real CDSC_i, CDSCD_i, CDSCD_a, CDSCDR_i, CIT_i, DVT0_i, CITR_i, CIT_a;
real DVT1_i, DVT1SS_i, PHIN_i, ETA0_i, ETA0_a, ETA0R_i, DSUB_i, VSAT_i, VSATR_i, VSATR_t;
real DVTP0_i, DVTP1_i ;
real K0_i, K01_i, K0SI_i, K0SI1_i, K2SI_i, K2SI1_i, PHIBE_i, K1_i, K11_i, K2SAT_i, K2SAT1_i;
real DELTAVSAT_i, PSAT_i, DELTAVSATCV_i, PSATCV_i, VSAT1_i, VSAT1R_i, PTWG_i, PTWGR_i, VSATCV_i;
real UP_i, U0_i, U0R_i, ETAMOB_i, NGATE_i, RDSW_i, UPR_i;
real PRWGS_i, PRWGD_i, WR_i, PDIBL1_i, PDIBL1R_i, PDIBL2_i,PDIBL2R_i, PDIBL2_a ;
real DROUT_i, PVAG_i;
real AIGBINV_i, AIGBINV1_i, BIGBINV_i, CIGBINV_i, EIGBINV_i, NIGBINV_i;
real AIGBACC_i, AIGBACC1_i, BIGBACC_i, CIGBACC_i, NIGBACC_i;
real AIGC_i, AIGC1_i, BIGC_i, CIGC_i, PIGCD_i;
real AIGS_i, AIGS1_i, BIGS_i, CIGS_i, NTOX_i, POXEDGE_i;
real AIGD_i, AIGD1_i, BIGD_i, CIGD_i;
real AGIDL_i, BGIDL_i, CGIDL_i, EGIDL_i, PGIDL_i;
real AGISL_i, BGISL_i, CGISL_i, EGISL_i, PGISL_i;
real ALPHA0_i, ALPHA1_i, ALPHAII0_i, ALPHAII1_i, BETA0_i;
real BETAII0_i, BETAII1_i, BETAII2_i, ESATII_i;
real LII_i, SII0_i, SII1_i, SII2_i, SIID_i, TII_i;
real MEXP_i, MEXPR_i;
real PCLM_i, PCLMG_i, PCLMCV_i, PCLM_a, PCLMR_i;
real A1_i, A2_i, A11_i, A21_i;
real K1RSCE_i, LPE0_i, DVTSHIFT_i, DVTSHIFT_a, DVTSHIFTR_i ;
real UA_i, UC_i, EU_i, UD_i, UCS_i, UAR_i, EUR_i, UCR_i, UDR_i, UA_a, UD_a, UC_a, EU_a;
real UA1_i, UA1R_i, UC1_i, UD1_i, UCSTE_i, UTE_i, UTL_i, EMOBT_i, UC1R_i, UD1R_i, UTER_i, UTLR_i;
real PTWGT_i;
real AT_i, ATCV_i, ATR_i;
real RDW_i, RSW_i;
real PRT_i, KT1_i, TSS_i, IIT_i, IGT_i, TGIDL_i;
real NTGEN_i, AIGEN_i, BIGEN_i;
real K0SISAT_i, K0SISAT1_i;
real K2SISAT_i, K2SISAT1_i;
real K2_i, K21_i;

// Variables of Unified Finfet Compact Model
real Cins, Ach, Weff_UFCM, qdep,rc, vth_fixed_factor_Sub, vth_fixed_factor_SI, qm, Qdep_ov_Cins, qi_acc_for_QM;
real fieldnormalizationfactor, auxQMfact, QMFACTORCVfinal;
real psipclamp, sqrtpsip, nq, F0;

real LSP_i;

`Cfringe_2d_vars();

//===================================================
//               analog block begins
//===================================================
analog begin

    // ************************************************
    // *      Geometry dependent calculations         *
    // ************************************************
    begin : CMGBiasIndepCalc

        // Variable Initialization to Prevent Hidden States
        qid         =  0.0;
        qis         =  0.0;
        qba         =  0.0;
        T11         =  0.0;
        T12         =  0.0;
        ids         =  0.0;
        sigvds      =  0.0;
        Iii         =  0.0;
        qd          =  0.0;
        qg          =  0.0;
        qs          =  0.0;
        qb          =  0.0;
        Weff0       =  0.0;
        WeffCV0     =  0.0;
        CJS_t       =  0.0;
        CJSWS_t     =  0.0;
        CJSWGS_t    =  0.0;
        CJD_t       =  0.0;
        CJSWD_t     =  0.0;
        CJSWGD_t    =  0.0;
        PBS_t       =  0.0;
        PBSWS_t     =  0.0;
        PBSWGS_t    =  0.0;
        PBD_t       =  0.0;
        PBSWD_t     =  0.0;
        PBSWGD_t    =  0.0;
        JSS_t       =  0.0;
        JSWS_t      =  0.0;
        JSWGS_t     =  0.0;
        JSD_t       =  0.0;
        JSWD_t      =  0.0;
        JSWGD_t     =  0.0;
        JTSS_t      =  0.0;
        JTSSWS_t    =  0.0;
        JTSSWGS_t   =  0.0;
        JTSD_t      =  0.0;
        JTSSWD_t    =  0.0;
        JTSSWGD_t   =  0.0;
        NJTS_t      =  0.0;
        NJTSSW_t    =  0.0;
        NJTSSWG_t   =  0.0;
        NJTSD_t     =  0.0;
        NJTSSWD_t   =  0.0;
        NJTSSWGD_t  =  0.0;
        Ies         =  0.0;
        Ied         =  0.0;
        Czbs        =  0.0;
        Czbssw      =  0.0;
        Czbsswg     =  0.0;
        Czbd        =  0.0;
        Czbdsw      =  0.0;
        Czbdswg     =  0.0;
        Qes         =  0.0;
        Qed         =  0.0;
        Qeg         =  0.0;
        Isbs        =  0.0;
        Isbd        =  0.0;
        Nvtms       =  0.0;
        Nvtmd       =  0.0;
        SslpRev     =  0.0;
        IVjsmRev    =  0.0;
        VjsmRev     =  0.0;
        SslpFwd     =  0.0;
        IVjsmFwd    =  0.0;
        VjsmFwd     =  0.0;
        DslpRev     =  0.0;
        IVjdmRev    =  0.0;
        VjdmRev     =  0.0;
        DslpFwd     =  0.0;
        IVjdmFwd    =  0.0;
        VjdmFwd     =  0.0;
        XExpBVS     =  0.0;
        XExpBVD     =  0.0;
        idsgen      =  0.0;
        q0          =  0.0;
        Tcen        =  0.0;
        MTcen       =  0.0;
        Rdrain      =  0.0;
        Rsource     =  0.0;
        Cfr_geo     =  0.0;
        igbinv      =  0.0;
        igbs        =  0.0;
        igbd        =  0.0;
        igbacc      =  0.0;
        igcs        =  0.0;
        igcd        =  0.0;
        igidl       =  0.0;
        igisl       =  0.0;
        igs         =  0.0;
        igd         =  0.0;
        cox_acc     =  0.0;
        CGSO_i      =  0.0;
        CGDO_i      =  0.0;
        qb_acc      =  0.0;
        qg_acc      =  0.0;
        qgs_fr      =  0.0;
        qgd_fr      =  0.0;
        qds_fr      =  0.0;
        qgs_parasitic  =  0.0;
        qgd_parasitic  =  0.0;
        FNPowerAt1Hz   =  0.0;
        Gtnoi       =  0.0;
        gspr        =  0.0;
        gdpr        =  0.0;
        Dr          =  1.0;
        CDSCDR_i    =  0.0;
        ETA0R_i     =  0.0;
        VSAT1R_i    =  0.0;
        VSAT1R_t    =  0.0;
        MEXPR_i     =  0.0;
        MEXPR_t     =  0.0;
        PTWGR_i     =  0.0;
        PTWGR_t     =  0.0;
        PDIBL1R_i   =  0.0;
        PDIBL2R_i   =  0.0;
        PHIBE_i     =  0.0;
        K1_i        =  0.0;
        K11_i       =  0.0;
        K2SAT_i     =  0.0;
        K2SAT1_i    =  0.0;
        KSATIVR_i   =  0.0;
        K2_i        =  0.0;
        K21_i       =  0.0;
        UC_i        =  0.0;
        UC1_i       =  0.0;
        UC_t        =  0.0;
        U0R_i       =  0.0;
        UPR_i       =  0.0;
        EUR_i       =  0.0;
        ATR_i       =  0.0;
        CITR_i      =  0.0;
        ETA0R_i     =  0.0;
        DVTP0_i     =  0.0;
        DVTP1_i     =  0.0;
        PDIBL2R_i   =  0.0;
        PCLMR_i     =  0.0;
        LeffCV_acc  =  0.0;
        RDDRR_t     =  0.0;
        RSDRR_t     =  0.0;
        Rdsi        =  0.0;
        T3y         =  0.0;
        Tcen0       =  0.0;
        veseff      =  0.0;
        U0R_t       =  0.0;
        UAR_t       =  0.0;
        UCR_t       =  0.0;
        UDR_t       =  0.0;
        VSAT_a      =  0.0;
        DVTSHIFTR_i =  0.0;
        UA1R_i      =  0.0;
        UAR_i       =  0.0;
        UC1R_i      =  0.0;
        UCR_i       =  0.0;
        UD1R_i      =  0.0;
        UDR_i       =  0.0;
        UTER_i      =  0.0;
        UTLR_i      =  0.0;
        VSATR_i     =  0.0;
        VSATR_t     =  0.0;
        u0r         =  0.0;

        // Thermal Noise
        sid         =  0.0;

        `ifdef __TNOIMOD1__
            ctnoi  =  0.0;
            sf     =  0.0;
            C0     =  0.0;
            gamma  =  0.0;
            delta  =  0.0;
        `endif

        `ifdef __RGATEMOD__
            ggeltd  =  0.0;
        `endif

        // Unified FinFET Model
        qm    =  1.0;
        Cins  =  1.0;
        Ach   =  1.0;
        Weff_UFCM  =  1.0;
        qdep  =  -1.0;
        rc    =  1.0;
        vth_fixed_factor_Sub  =  1.0;
        vth_fixed_factor_SI   =  1.0;
        qi_acc_for_QM         =  0.0;
        fieldnormalizationfactor  =  0.0;
        auxQMfact  =  0.0;
        QMFACTORCVfinal  =  0.0;
        psipclamp  =  1.0;
        sqrtpsip   =  1.0;
        nq  =  1.0;
        F0  =  0.0;
        e0  =  0.0;
        e1  =  0.0;
        e2  =  0.0;
        Qdep_ov_Cins  =  0.0;

        // Constants
        if ( TYPE == `ntype ) begin
            devsign  =  1;
        end else begin
            devsign  =  -1;
        end

        epssub    =  EPSRSUB * `EPS0;
        epssp     =  EPSRSP * `EPS0;
        cbox      =  EPSROX * `EPS0 / EOTBOX;
        epsratio  =  EPSRSUB / EPSROX;

        if ($port_connected(t) == 1) begin
            `ifdef __SHMOD__
                if (SHMOD == 0) begin
                    if (SH_WARN == 1) begin
                        $strobe("The optional 5th terminal is present but not active because SHMOD=0.");
                    end
                end
            `else
                Temp(t) <+ 0.0;
                if (SH_WARN == 1) begin
                    $strobe("The optional 5th terminal is present but not active because the model was not compiled with self-heating enabled (__SHMOD__ was not activated).");
                end
            `endif
        end

        // Constants for Quantum Mechanical Effects
        mx       =  0.916 * `MEL;
        mxprime  =  0.190 * `MEL;
        md       =  0.190 * `MEL;
        mdprime  =  0.417 * `MEL;
        gprime   =  4.0;
        gfactor  =  2.0;

        // Effective Channel Length for I-V / C-V
        Lg       =  L + XL;
        deltaL   =  LINT + LL * pow(Lg, -LLN);
        deltaL1  =  LINT + LL * pow(Lg+DLBIN, -LLN);
        deltaLCV =  DLC + LLC * pow(Lg, -LLN);
        Leff     =  Lg - 2.0 * deltaL;
        Leff1    =  Lg + DLBIN - 2.0 * deltaL1;        //Used in the binning equations only
        LeffCV   =  Lg - 2.0 * deltaLCV;
        if (BULKMOD != 0) LeffCV_acc  =  LeffCV - DLCACC;

        // Total Fins
        NFINtotal  =  NFIN * NF;

        // Range Checking on Leff and Leff1
        if (Leff <= 0.0) begin
            $strobe("Fatal: Leff = %e is not positive.", Leff);
            $finish(0);
        end else if (Leff <= 1.0e-9) begin
            $strobe("Warning: Leff = %e <= 1.0e-9.", Leff);
        end

        if (Leff1 <= 0.0) begin
            $strobe("Fatal: Leff1 = %e is not positive.", Leff1);
            $finish(0);
        end else if (Leff1 <= 1.0e-9) begin
            $strobe("Warning: Leff1 = %e <= 1.0e-9.", Leff1);
        end

        // Binning
        Inv_L      =  1.0e-6 / (Leff1);
        Inv_NFIN   =  1.0 / NFIN;
        Inv_LNFIN  =  1.0e-6 / (Leff1 * NFIN);

        // Nbody Binning Equation for UFCM Parameters
        NBODY_i  =  NBODY + Inv_L * LNBODY + Inv_NFIN * NNBODY + Inv_LNFIN * PNBODY;

        if (NBODYN1 != 0.0) begin
            NBODY_i  =  NBODY_i * (1.0 + NBODYN1/NFIN * lln(1.0 + NFIN/NBODYN2));
        end

        // Model Parameters for Unified FinFET Compact Model by Juan Duarte 10/2013
        case (GEOMOD)
            0: begin // Double Gate
                if (!$param_given(TFIN_TOP) || !$param_given(TFIN_BASE)) begin
                    Weff_UFCM  =  2.0 * HFIN;
                    Cins  =  Weff_UFCM * EPSROX * `EPS0 / EOT;
                    Ach   =  HFIN * TFIN;
                    rc    =  (2.0 * Cins / (Weff_UFCM * Weff_UFCM * epssub / Ach));
                    Qdep_ov_Cins  =  -`q * NBODY_i * Ach / Cins;
                end else begin
                    Weff_UFCM  =  2.0 * sqrt( HFIN * HFIN + (TFIN_TOP - TFIN_BASE) * (TFIN_TOP - TFIN_BASE) / 4.0);
                    Cins  =  Weff_UFCM * EPSROX * `EPS0 / EOT;
                    Ach   =  HFIN * (TFIN_TOP + TFIN_BASE) / 2.0;
                    rc    =  (2.0 * Cins / (Weff_UFCM * Weff_UFCM * epssub / Ach));
                    Qdep_ov_Cins  =  -`q * NBODY_i * Ach / Cins;
                end
            end
            1: begin // Triple Gate
                if (!$param_given(TFIN_TOP) || !$param_given(TFIN_BASE)) begin
                    Weff_UFCM  =  2.0 * HFIN + TFIN;
                    Cins  =  Weff_UFCM * EPSROX * `EPS0 / EOT;
                    Ach   =  HFIN * TFIN;
                    rc    =  (2.0 * Cins / (Weff_UFCM * Weff_UFCM * epssub / Ach));
                    Qdep_ov_Cins  =  -`q * NBODY_i * Ach / Cins;
                end else begin
                    Weff_UFCM  =  2.0 * sqrt(HFIN * HFIN + (TFIN_TOP - TFIN_BASE) * (TFIN_TOP - TFIN_BASE) / 4.0) + TFIN_TOP;
                    Cins  =  Weff_UFCM * EPSROX * `EPS0 / EOT;
                    Ach   =  HFIN * (TFIN_TOP + TFIN_BASE) / 2.0;
                    rc    =  (2.0 * Cins /(Weff_UFCM * Weff_UFCM * epssub / Ach));
                    Qdep_ov_Cins  =  -`q * NBODY_i * Ach / Cins;
                end
            end
            2: begin // Quadruple Gate
                if (!$param_given(TFIN_TOP) || !$param_given(TFIN_BASE)) begin
                    Weff_UFCM  =  2.0 * HFIN + 2.0 * TFIN;
                    Cins  =  Weff_UFCM * EPSROX * `EPS0 / EOT;
                    Ach   =  HFIN * TFIN;
                    rc    =  (2.0 * Cins / (Weff_UFCM * Weff_UFCM * epssub / Ach));
                    Qdep_ov_Cins  =  -`q * NBODY_i * Ach / Cins;
                end else begin
                    Weff_UFCM  =  2.0 * sqrt(HFIN * HFIN + (TFIN_TOP - TFIN_BASE) * (TFIN_TOP - TFIN_BASE) / 4.0) + TFIN_TOP + TFIN_BASE;
                    Cins  =  Weff_UFCM * EPSROX * `EPS0 / EOT;
                    Ach   =  HFIN * (TFIN_TOP + TFIN_BASE) / 2.0;
                    rc    =  (2.0 * Cins / (Weff_UFCM * Weff_UFCM * epssub / Ach));
                    Qdep_ov_Cins  =  -`q * NBODY_i * Ach / Cins;
                end
            end
            3:  begin // Cylindrical Gate
                Weff_UFCM  =  `M_PI * D;
                Cins  =  2.0 * `M_PI * EPSROX * `EPS0 / ln(1.0 + 2.0 * EOT / D);
                Ach   =  `M_PI * D * D / 4.0;
                rc    =  (2.0 * Cins / (Weff_UFCM * Weff_UFCM * epssub / Ach));
                Qdep_ov_Cins  =  -`q * NBODY_i * Ach / Cins;
            end
            4:  begin // Unified Model
                Weff_UFCM  =  W_UFCM;
                Cins  =  CINS_UFCM;
                Ach   =  ACH_UFCM;
                rc    =  (2.0 * Cins / (Weff_UFCM * Weff_UFCM * epssub / Ach));
                Qdep_ov_Cins  =  -`q * NBODY_i * Ach / Cins;
            end
        endcase

        // Cox Definition
        cox  =  Cins / Weff_UFCM;
        if (BULKMOD != 0) begin
            cox_acc =  cox * EOT / EOTACC;
        end

        // Effective Width Calculation
        Weff0    =  Weff_UFCM - DELTAW;
        WeffCV0  =  Weff_UFCM - DELTAWCV;

        // SCE Scaling Length
        scl  =  sqrt((epssub * Ach / Cins) * (1.0 + Ach * Cins / (2.0 * epssub * Weff_UFCM * Weff_UFCM)));

        // Binning Equations (Process Parameters)
        PHIG_i   =  PHIG + Inv_L * LPHIG + Inv_NFIN * NPHIG + Inv_LNFIN * PPHIG;
        NGATE_i  =  NGATE + Inv_L * LNGATE + Inv_NFIN * NNGATE + Inv_LNFIN * PNGATE;

        // Binning Equations (Model Parameters)
        CIT_i       =  CIT + Inv_L * LCIT + Inv_NFIN * NCIT + Inv_LNFIN * PCIT;
        CDSC_i      =  CDSC + Inv_L * LCDSC + Inv_NFIN * NCDSC + Inv_LNFIN * PCDSC;
        CDSCD_i     =  CDSCD + Inv_L * LCDSCD + Inv_NFIN * NCDSCD + Inv_LNFIN * PCDSCD;
        DVT0_i      =  DVT0 + Inv_L * LDVT0 + Inv_NFIN * NDVT0 + Inv_LNFIN * PDVT0;
        DVT1_i      =  DVT1 + Inv_L * LDVT1 + Inv_NFIN * NDVT1 + Inv_LNFIN * PDVT1;
        DVT1SS_i    =  DVT1SS + Inv_L * LDVT1SS + Inv_NFIN * NDVT1SS + Inv_LNFIN * PDVT1SS;
        PHIN_i      =  PHIN + Inv_L * LPHIN + Inv_NFIN * NPHIN + Inv_LNFIN * PPHIN;
        ETA0_i      =  ETA0 + Inv_L * LETA0 + Inv_NFIN * NETA0 + Inv_LNFIN * PETA0;
        DSUB_i      =  DSUB + Inv_L * LDSUB + Inv_NFIN * NDSUB + Inv_LNFIN * PDSUB;
        K1RSCE_i    =  K1RSCE + Inv_L * LK1RSCE + Inv_NFIN * NK1RSCE + Inv_LNFIN * PK1RSCE;
        LPE0_i      =  LPE0 + Inv_L * LLPE0 + Inv_NFIN * NLPE0 + Inv_LNFIN * PLPE0;
        DVTSHIFT_i  =  DVTSHIFT + Inv_L * LDVTSHIFT + Inv_NFIN * NDVTSHIFT + Inv_LNFIN * PDVTSHIFT;
        K0_i        =  K0 + Inv_L * LK0 + Inv_NFIN * NK0 + Inv_LNFIN * PK0;
        K01_i       =  K01 + Inv_L * LK01 + Inv_NFIN * NK01 + Inv_LNFIN * PK01;
        K0SI_i      =  K0SI + Inv_L * LK0SI + Inv_NFIN * NK0SI + Inv_LNFIN * PK0SI;
        K0SI1_i     =  K0SI1 + Inv_L * LK0SI1 + Inv_NFIN * NK0SI1 + Inv_LNFIN * PK0SI1;
        K2SI_i      =  K2SI + Inv_L * LK2SI + Inv_NFIN * NK2SI + Inv_LNFIN * PK2SI;
        K2SI1_i     =  K2SI1 + Inv_L * LK2SI1 + Inv_NFIN * NK2SI1 + Inv_LNFIN * PK2SI1;
        K0SISAT_i   =  K0SISAT + Inv_L * LK0SISAT + Inv_NFIN * NK0SISAT + Inv_LNFIN * PK0SISAT;
        K0SISAT1_i  =  K0SISAT1 + Inv_L * LK0SISAT1 + Inv_NFIN * NK0SISAT1 + Inv_LNFIN * PK0SISAT1;
        K2SISAT_i   =  K2SISAT + Inv_L * LK2SISAT + Inv_NFIN * NK2SISAT + Inv_LNFIN * PK2SISAT;
        K2SISAT1_i  =  K2SISAT1 + Inv_L * LK2SISAT1 + Inv_NFIN * NK2SISAT1 + Inv_LNFIN * PK2SISAT1;

        if (BULKMOD != 0) begin
            if (BULKMOD == 2) begin
                K2_i      =  K2 + Inv_L * LK2 + Inv_NFIN * NK2 + Inv_LNFIN * PK2;
                K21_i     =  K21 + Inv_L * LK21 + Inv_NFIN * NK21 + Inv_LNFIN * PK21;
                K2SAT_i   =  K2SAT + Inv_L * LK2SAT + Inv_NFIN * NK2SAT + Inv_LNFIN * PK2SAT;
                K2SAT1_i  =  K2SAT1 + Inv_L * LK2SAT1 + Inv_NFIN * NK2SAT1 + Inv_LNFIN * PK2SAT1;
            end
            PHIBE_i    =  PHIBE + Inv_L * LPHIBE + Inv_NFIN * NPHIBE + Inv_LNFIN * PPHIBE;
            K1_i       =  K1 + Inv_L * LK1 + Inv_NFIN * NK1 + Inv_LNFIN * PK1;
            K11_i      =  K11 + Inv_L * LK11 + Inv_NFIN * NK11 + Inv_LNFIN * PK11;
        end
        QMFACTOR_i    =  QMFACTOR + Inv_L * LQMFACTOR + Inv_NFIN * NQMFACTOR + Inv_LNFIN * PQMFACTOR;
        QMTCENCV_i    =  QMTCENCV + Inv_L * LQMTCENCV + Inv_NFIN * NQMTCENCV + Inv_LNFIN * PQMTCENCV;
        QMTCENCVA_i   =  QMTCENCVA + Inv_L * LQMTCENCVA + Inv_NFIN * NQMTCENCVA + Inv_LNFIN * PQMTCENCVA;
        VSAT_i        =  VSAT + Inv_L * LVSAT + Inv_NFIN * NVSAT + Inv_LNFIN * PVSAT;
        VSAT1_i       =  VSAT1 + Inv_L * LVSAT1 + Inv_NFIN * NVSAT1 + Inv_LNFIN * PVSAT1;
        VSATCV_i      =  VSATCV + Inv_L * LVSATCV + Inv_NFIN * NVSATCV + Inv_LNFIN * PVSATCV;
        DELTAVSAT_i   =  DELTAVSAT + Inv_L * LDELTAVSAT + Inv_NFIN * NDELTAVSAT + Inv_LNFIN * PDELTAVSAT;
        PSAT_i        =  PSAT + Inv_L * LPSAT + Inv_NFIN * NPSAT + Inv_LNFIN * PPSAT;
        DELTAVSATCV_i =  DELTAVSATCV + Inv_L * LDELTAVSATCV + Inv_NFIN * NDELTAVSATCV + Inv_LNFIN * PDELTAVSATCV;
        PSATCV_i      =  PSATCV + Inv_L * LPSATCV + Inv_NFIN * NPSATCV + Inv_LNFIN * PPSATCV;
        KSATIV_i      =  KSATIV + Inv_L * LKSATIV + Inv_NFIN * NKSATIV + Inv_LNFIN * PKSATIV;
        MEXP_i        =  MEXP + Inv_L * LMEXP + Inv_NFIN * NMEXP + Inv_LNFIN * PMEXP;
        PTWG_i        =  PTWG + Inv_L * LPTWG + Inv_NFIN * NPTWG + Inv_LNFIN * PPTWG;
        U0_i          =  U0 + Inv_L * LU0 + Inv_NFIN * NU0 + Inv_LNFIN * PU0;
        ETAMOB_i      =  ETAMOB + Inv_L * LETAMOB + Inv_NFIN * NETAMOB + Inv_LNFIN * PETAMOB;
        UP_i          =  UP + Inv_L * LUP + Inv_NFIN * NUP + Inv_LNFIN * PUP;
        UA_i          =  UA + Inv_L * LUA + Inv_NFIN * NUA + Inv_LNFIN * PUA;
        if (BULKMOD != 0) begin
            UC_i       =  UC + Inv_L * LUC + Inv_NFIN * NUC + Inv_LNFIN * PUC;
        end
        EU_i          =  EU + Inv_L * LEU + Inv_NFIN * NEU + Inv_LNFIN * PEU;
        UD_i          =  UD + Inv_L * LUD + Inv_NFIN * NUD + Inv_LNFIN * PUD;
        UCS_i         =  UCS + Inv_L * LUCS + Inv_NFIN * NUCS + Inv_LNFIN * PUCS;
        PCLM_i        =  PCLM + Inv_L * LPCLM + Inv_NFIN * NPCLM + Inv_LNFIN * PPCLM;
        PCLMG_i       =  PCLMG + Inv_L * LPCLMG + Inv_NFIN * NPCLMG + Inv_LNFIN * PPCLMG;
        PCLMCV_i      =  PCLMCV + Inv_L * LPCLMCV + Inv_NFIN * NPCLMCV + Inv_LNFIN * PPCLMCV;
        A1_i          =  A1 + Inv_L * LA1 + Inv_NFIN * NA1 + Inv_LNFIN * PA1;
        A11_i         =  A11 + Inv_L * LA11 + Inv_NFIN * NA11 + Inv_LNFIN * PA11;
        A2_i          =  A2 + Inv_L * LA2 + Inv_NFIN * NA2 + Inv_LNFIN * PA2;
        A21_i         =  A21 + Inv_L * LA21 + Inv_NFIN * NA21 + Inv_LNFIN * PA21;
        RDSW_i        =  RDSW + Inv_L * LRDSW + Inv_NFIN * NRDSW + Inv_LNFIN * PRDSW;
        RSW_i         =  RSW + Inv_L * LRSW + Inv_NFIN * NRSW + Inv_LNFIN * PRSW;
        RDW_i         =  RDW + Inv_L * LRDW + Inv_NFIN * NRDW + Inv_LNFIN * PRDW;
        PRWGD_i       =  PRWGD + Inv_L * LPRWGD + Inv_NFIN * NPRWGD + Inv_LNFIN * PPRWGD;
        PRWGS_i       =  PRWGS + Inv_L * LPRWGS + Inv_NFIN * NPRWGS + Inv_LNFIN * PPRWGS;
        WR_i          =  WR + Inv_L * LWR + Inv_NFIN * NWR + Inv_LNFIN * PWR;
        PDIBL1_i      =  PDIBL1 + Inv_L * LPDIBL1 + Inv_NFIN * NPDIBL1 + Inv_LNFIN * PPDIBL1;
        PDIBL2_i      =  PDIBL2 + Inv_L * LPDIBL2 + Inv_NFIN * NPDIBL2 + Inv_LNFIN * PPDIBL2;
        DROUT_i       =  DROUT + Inv_L * LDROUT + Inv_NFIN * NDROUT + Inv_LNFIN * PDROUT;
        PVAG_i        =  PVAG + Inv_L * LPVAG + Inv_NFIN * NPVAG + Inv_LNFIN * PPVAG;
        AIGBINV_i     =  AIGBINV + Inv_L * LAIGBINV  + Inv_NFIN * NAIGBINV + Inv_LNFIN * PAIGBINV;
        AIGBINV1_i    =  AIGBINV1 + Inv_L * LAIGBINV1 + Inv_NFIN * NAIGBINV1 + Inv_LNFIN * PAIGBINV1;
        BIGBINV_i     =  BIGBINV + Inv_L * LBIGBINV + Inv_NFIN * NBIGBINV + Inv_LNFIN * PBIGBINV;
        CIGBINV_i     =  CIGBINV + Inv_L * LCIGBINV + Inv_NFIN * NCIGBINV + Inv_LNFIN * PCIGBINV;
        EIGBINV_i     =  EIGBINV + Inv_L * LEIGBINV + Inv_NFIN * NEIGBINV + Inv_LNFIN * PEIGBINV;
        NIGBINV_i     =  NIGBINV + Inv_L * LNIGBINV + Inv_NFIN * NNIGBINV + Inv_LNFIN * PNIGBINV;
        AIGBACC_i     =  AIGBACC + Inv_L * LAIGBACC + Inv_NFIN * NAIGBACC + Inv_LNFIN * PAIGBACC;
        AIGBACC1_i    =  AIGBACC1 + Inv_L * LAIGBACC1 + Inv_NFIN * NAIGBACC1 + Inv_LNFIN * PAIGBACC1;
        BIGBACC_i     =  BIGBACC + Inv_L * LBIGBACC + Inv_NFIN * NBIGBACC + Inv_LNFIN * PBIGBACC;
        CIGBACC_i     =  CIGBACC + Inv_L * LCIGBACC + Inv_NFIN * NCIGBACC + Inv_LNFIN * PCIGBACC;
        NIGBACC_i     =  NIGBACC + Inv_L * LNIGBACC + Inv_NFIN * NNIGBACC + Inv_LNFIN * PNIGBACC;
        AIGC_i        =  AIGC + Inv_L * LAIGC + Inv_NFIN * NAIGC + Inv_LNFIN * PAIGC;
        AIGC1_i       =  AIGC1 + Inv_L * LAIGC1 + Inv_NFIN * NAIGC1 + Inv_LNFIN * PAIGC1;
        BIGC_i        =  BIGC + Inv_L * LBIGC + Inv_NFIN * NBIGC + Inv_LNFIN * PBIGC;
        CIGC_i        =  CIGC + Inv_L * LCIGC + Inv_NFIN * NCIGC + Inv_LNFIN * PCIGC;
        PIGCD_i       =  PIGCD + Inv_L * LPIGCD + Inv_NFIN * NPIGCD + Inv_LNFIN * PPIGCD;
        AIGS_i        =  AIGS + Inv_L * LAIGS + Inv_NFIN * NAIGS + Inv_LNFIN * PAIGS;
        AIGS1_i       =  AIGS1 + Inv_L * LAIGS1 + Inv_NFIN * NAIGS1 + Inv_LNFIN * PAIGS1;
        BIGS_i        =  BIGS + Inv_L * LBIGS + Inv_NFIN * NBIGS + Inv_LNFIN * PBIGS;
        CIGS_i        =  CIGS + Inv_L * LCIGS + Inv_NFIN * NCIGS + Inv_LNFIN * PCIGS;
        AIGD_i        =  AIGD + Inv_L * LAIGD + Inv_NFIN * NAIGD + Inv_LNFIN * PAIGD;
        AIGD1_i       =  AIGD1 + Inv_L * LAIGD1 + Inv_NFIN * NAIGD1 + Inv_LNFIN * PAIGD1;
        BIGD_i        =  BIGD + Inv_L * LBIGD + Inv_NFIN * NBIGD + Inv_LNFIN * PBIGD;
        CIGD_i        =  CIGD + Inv_L * LCIGD + Inv_NFIN * NCIGD + Inv_LNFIN * PCIGD;
        NTOX_i        =  NTOX + Inv_L * LNTOX + Inv_NFIN * NNTOX + Inv_LNFIN * PNTOX;
        POXEDGE_i     =  POXEDGE + Inv_L * LPOXEDGE + Inv_NFIN * NPOXEDGE + Inv_LNFIN * PPOXEDGE;
        AGIDL_i       =  AGIDL + Inv_L * LAGIDL + Inv_NFIN * NAGIDL + Inv_LNFIN * PAGIDL;
        BGIDL_i       =  BGIDL + Inv_L * LBGIDL + Inv_NFIN * NBGIDL + Inv_LNFIN * PBGIDL;
        CGIDL_i       =  CGIDL + Inv_L * LCGIDL + Inv_NFIN * NCGIDL + Inv_LNFIN * PCGIDL;
        EGIDL_i       =  EGIDL + Inv_L * LEGIDL + Inv_NFIN * NEGIDL + Inv_LNFIN * PEGIDL;
        PGIDL_i       =  PGIDL + Inv_L * LPGIDL + Inv_NFIN * NPGIDL + Inv_LNFIN * PPGIDL;
        AGISL_i       =  AGISL + Inv_L * LAGISL + Inv_NFIN * NAGISL + Inv_LNFIN * PAGISL;
        BGISL_i       =  BGISL + Inv_L * LBGISL + Inv_NFIN * NBGISL + Inv_LNFIN * PBGISL;
        CGISL_i       =  CGISL + Inv_L * LCGISL + Inv_NFIN * NCGISL + Inv_LNFIN * PCGISL;
        EGISL_i       =  EGISL + Inv_L * LEGISL + Inv_NFIN * NEGISL + Inv_LNFIN * PEGISL;
        PGISL_i       =  PGISL + Inv_L * LPGISL + Inv_NFIN * NPGISL + Inv_LNFIN * PPGISL;
        ALPHA0_i      =  ALPHA0 + Inv_L * LALPHA0 + Inv_NFIN * NALPHA0 + Inv_LNFIN * PALPHA0;
        ALPHA1_i      =  ALPHA1 + Inv_L * LALPHA1 + Inv_NFIN * NALPHA1 + Inv_LNFIN * PALPHA1;
        ALPHAII0_i    =  ALPHAII0 + Inv_L * LALPHAII0 + Inv_NFIN * NALPHAII0 + Inv_LNFIN * PALPHAII0;
        ALPHAII1_i    =  ALPHAII1 + Inv_L * LALPHAII1 + Inv_NFIN * NALPHAII1 + Inv_LNFIN * PALPHAII1;
        BETA0_i       =  BETA0 + Inv_L * LBETA0 + Inv_NFIN * NBETA0 + Inv_LNFIN * PBETA0;
        BETAII0_i     =  BETAII0 + Inv_L * LBETAII0 + Inv_NFIN * NBETAII0 + Inv_LNFIN * PBETAII0;
        BETAII1_i     =  BETAII1 + Inv_L * LBETAII1 + Inv_NFIN * NBETAII1 + Inv_LNFIN * PBETAII1;
        BETAII2_i     =  BETAII2 + Inv_L * LBETAII2 + Inv_NFIN * NBETAII2 + Inv_LNFIN * PBETAII2;
        ESATII_i      =  ESATII + Inv_L * LESATII + Inv_NFIN * NESATII + Inv_LNFIN * PESATII;
        LII_i         =  LII + Inv_L * LLII + Inv_NFIN * NLII + Inv_LNFIN * PLII;
        SII0_i        =  SII0 + Inv_L * LSII0 + Inv_NFIN * NSII0 + Inv_LNFIN * PSII0;
        SII1_i        =  SII1 + Inv_L * LSII1 + Inv_NFIN * NSII1 + Inv_LNFIN * PSII1;
        SII2_i        =  SII2 + Inv_L * LSII2 + Inv_NFIN * NSII2 + Inv_LNFIN * PSII2;
        SIID_i        =  SIID + Inv_L * LSIID + Inv_NFIN * NSIID + Inv_LNFIN * PSIID;
        TII_i         =  TII + Inv_L * LTII + Inv_NFIN * NTII + Inv_LNFIN * PTII;
        CFS_i         =  CFS + Inv_L * LCFS + Inv_NFIN * NCFS + Inv_LNFIN * PCFS;
        CFD_i         =  CFD + Inv_L * LCFD + Inv_NFIN * NCFD + Inv_LNFIN * PCFD;
        COVS_i        =  COVS + Inv_L * LCOVS + Inv_NFIN * NCOVS + Inv_LNFIN * PCOVS;
        COVD_i        =  COVD + Inv_L * LCOVD + Inv_NFIN * NCOVD + Inv_LNFIN * PCOVD;
        CGSL_i        =  CGSL + Inv_L * LCGSL + Inv_NFIN * NCGSL + Inv_LNFIN * PCGSL;
        CGDL_i        =  CGDL + Inv_L * LCGDL + Inv_NFIN * NCGDL + Inv_LNFIN * PCGDL;
        CGBL_i        =  CGBL + Inv_L * LCGBL + Inv_NFIN * NCGBL + Inv_LNFIN * PCGBL;
        CKAPPAS_i     =  CKAPPAS + Inv_L * LCKAPPAS + Inv_NFIN * NCKAPPAS + Inv_LNFIN * PCKAPPAS;
        CKAPPAD_i     =  CKAPPAD + Inv_L * LCKAPPAD + Inv_NFIN * NCKAPPAD + Inv_LNFIN * PCKAPPAD;
        CKAPPAB_i     =  CKAPPAB + Inv_L * LCKAPPAB + Inv_NFIN * NCKAPPAB + Inv_LNFIN * PCKAPPAB;
        NTGEN_i       =  NTGEN + Inv_L * LNTGEN + Inv_NFIN * NNTGEN + Inv_LNFIN * PNTGEN;
        AIGEN_i       =  AIGEN + Inv_L * LAIGEN + Inv_NFIN * NAIGEN + Inv_LNFIN * PAIGEN;
        BIGEN_i       =  BIGEN + Inv_L * LBIGEN + Inv_NFIN * NBIGEN + Inv_LNFIN * PBIGEN;

        if (ASYMMOD != 0) begin
            CDSCDR_i     =  CDSCDR + Inv_L * LCDSCDR + Inv_NFIN * NCDSCDR + Inv_LNFIN * PCDSCDR;
            CITR_i       =  CITR + Inv_L * LCITR + Inv_NFIN * NCITR + Inv_LNFIN * PCITR;
            ETA0R_i      =  ETA0R + Inv_L * LETA0R + Inv_NFIN * NETA0R + Inv_LNFIN * PETA0R;
            VSAT1R_i     =  VSAT1R + Inv_L * LVSAT1R + Inv_NFIN * NVSAT1R + Inv_LNFIN * PVSAT1R;
            MEXPR_i      =  MEXPR + Inv_L * LMEXPR + Inv_NFIN * NMEXPR + Inv_LNFIN * PMEXPR;
            PTWGR_i      =  PTWGR + Inv_L * LPTWGR + Inv_NFIN * NPTWGR + Inv_LNFIN * PPTWGR;
            PDIBL1R_i    =  PDIBL1R + Inv_L * LPDIBL1R + Inv_NFIN * NPDIBL1R + Inv_LNFIN * PPDIBL1R;
            PDIBL2R_i    =  PDIBL2R + Inv_L * LPDIBL2R + Inv_NFIN * NPDIBL2R + Inv_LNFIN * PPDIBL2R;
            PCLMR_i      =  PCLMR + Inv_L * LPCLMR + Inv_NFIN * NPCLMR + Inv_LNFIN * PPCLMR;
            DVTSHIFTR_i  =  DVTSHIFTR + Inv_L * LDVTSHIFTR + Inv_NFIN * NDVTSHIFTR + Inv_LNFIN * PDVTSHIFTR;
            VSATR_i      =  VSATR + Inv_L * LVSATR + Inv_NFIN * NVSATR + Inv_LNFIN * PVSATR;
            KSATIVR_i    =  KSATIVR + Inv_L * LKSATIVR + Inv_NFIN * NKSATIVR + Inv_LNFIN * PKSATIVR;
            U0R_i        =  U0R + Inv_L * LU0R + Inv_NFIN * NU0R + Inv_LNFIN * PU0R;
            UAR_i        =  UAR + Inv_L * LUAR + Inv_NFIN * NUAR + Inv_LNFIN * PUAR;
            UPR_i        =  UPR + Inv_L * LUPR + Inv_NFIN * NUPR + Inv_LNFIN * PUPR;
            if (BULKMOD != 0) begin
                UCR_i    =  UCR + Inv_L * LUCR + Inv_NFIN * NUCR + Inv_LNFIN * PUCR;
            end
            EUR_i        =  EUR + Inv_L * LEUR + Inv_NFIN * NEUR + Inv_LNFIN * PEUR;
            UDR_i        =  UDR + Inv_L * LUDR + Inv_NFIN * NUDR + Inv_LNFIN * PUDR;
        end

        `ifdef __NQSMOD1__
            if (NQSMOD == 1 && XRCRG1 != 0.0) begin
                XRCRG1_i  =  XRCRG1 + Inv_L * LXRCRG1 + Inv_NFIN * NXRCRG1 + Inv_LNFIN * PXRCRG1;
                XRCRG2_i  =  XRCRG2 + Inv_L * LXRCRG2 + Inv_NFIN * NXRCRG2 + Inv_LNFIN * PXRCRG2;
            end else begin
                XRCRG1_i  =  0.0;
                XRCRG2_i  =  0.0;
            end
        `else
            if (NQSMOD == 1) begin
                $strobe(" Although the model selector NQSMOD is set to 1, the NQS gate resistance model is not activated in the Verilog-A code. Please uncomment \"`define __NQSMOD1__\" in bsimcmg.va to activate it.");
            end
        `endif

        `ifdef __NQSMOD2__
            if (NQSMOD == 2 && XRCRG1 != 0.0) begin
                XRCRG1_i  =  XRCRG1 + Inv_L * LXRCRG1 + Inv_NFIN * NXRCRG1 + Inv_LNFIN * PXRCRG1;
                XRCRG2_i  =  XRCRG2 + Inv_L * LXRCRG2 + Inv_NFIN * NXRCRG2 + Inv_LNFIN * PXRCRG2;
            end else begin
                XRCRG1_i  =  0.0;
                XRCRG2_i  =  0.0;
            end
        `else
            if (NQSMOD == 2) begin
                $strobe(" Although the model selector NQSMOD is set to 2, the NQS gate resistance model is not activated in the Verilog-A code. Please uncomment \"`define __NQSMOD2__\" in bsimcmg.va to activate it.");
            end
        `endif

        UTE_i    =  UTE + Inv_L * LUTE + Inv_NFIN * NUTE + Inv_LNFIN * PUTE;
        UTL_i    =  UTL + Inv_L * LUTL + Inv_NFIN * NUTL + Inv_LNFIN * PUTL;
        EMOBT_i  =  EMOBT + Inv_L * LEMOBT + Inv_NFIN * NEMOBT + Inv_LNFIN * PEMOBT;
        UA1_i    =  UA1 + Inv_L * LUA1 + Inv_NFIN * NUA1 + Inv_LNFIN * PUA1;

        if (BULKMOD != 0) begin
            UC1_i  =  UC1 + Inv_L * LUC1 + Inv_NFIN * NUC1 + Inv_LNFIN * PUC1;
        end
        UD1_i    =  UD1 + Inv_L * LUD1 + Inv_NFIN * NUD1 + Inv_LNFIN * PUD1;
        UCSTE_i  =  UCSTE + Inv_L * LUCSTE + Inv_NFIN * NUCSTE + Inv_LNFIN * PUCSTE;
        PTWGT_i  =  PTWGT + Inv_L * LPTWGT + Inv_NFIN * NPTWGT + Inv_LNFIN * PPTWGT;
        AT_i     =  AT + Inv_L * LAT + Inv_NFIN * NAT + Inv_LNFIN * PAT;
        ATCV_i   =  ATCV + Inv_L * LATCV  + Inv_NFIN * NATCV  + Inv_LNFIN * PATCV;
        PRT_i    =  PRT + Inv_L * LPRT + Inv_NFIN * NPRT + Inv_LNFIN * PPRT;
        KT1_i    =  KT1 + Inv_L * LKT1 + Inv_NFIN * NKT1 + Inv_LNFIN * PKT1;
        TSS_i    =  TSS + Inv_L * LTSS + Inv_NFIN * NTSS + Inv_LNFIN * PTSS;
        IIT_i    =  IIT + Inv_L * LIIT + Inv_NFIN * NIIT + Inv_LNFIN * PIIT;
        TGIDL_i  =  TGIDL + Inv_L * LTGIDL + Inv_NFIN * NTGIDL + Inv_LNFIN * PTGIDL;
        IGT_i    =  IGT + Inv_L * LIGT + Inv_NFIN * NIGT + Inv_LNFIN * PIGT;

        if (ASYMMOD != 0) begin
            UTER_i  =  UTER + Inv_L * LUTER + Inv_NFIN * NUTER + Inv_LNFIN * PUTER;
            UTLR_i  =  UTLR + Inv_L * LUTLR + Inv_NFIN * NUTLR + Inv_LNFIN * PUTLR;
            UA1R_i  =  UA1R + Inv_L * LUA1R + Inv_NFIN * NUA1R + Inv_LNFIN * PUA1R;
            UD1R_i  =  UD1R + Inv_L * LUD1R + Inv_NFIN * NUD1R + Inv_LNFIN * PUD1R;
            ATR_i   =  ATR + Inv_L * LATR + Inv_NFIN * NATR + Inv_LNFIN * PATR;
            if (BULKMOD != 0) begin
                UC1R_i  =  UC1R + Inv_L * LUC1R + Inv_NFIN * NUC1R + Inv_LNFIN * PUC1R;
            end
        end

        // Geometrical Scaling
        // NFIN Scaling
        if (PHIGN1 != 0.0) begin
            PHIG_i  =  PHIG_i * (1.0 + PHIGN1 / NFIN * lln(1.0 + NFIN / PHIGN2));
        end

        if (ETA0N1 != 0.0) begin
            ETA0_i  =  ETA0_i * (1.0 + ETA0N1 / NFIN * lln(1.0 + NFIN / ETA0N2));
        end

        if (CDSCN1 != 0.0) begin
            CDSC_i  =  CDSC_i * (1.0 + CDSCN1 / NFIN * lln(1.0 + NFIN / CDSCN2));
        end

        if (CDSCDN1 != 0.0) begin
            CDSCD_i  =  CDSCD_i * (1.0 + CDSCDN1 / NFIN * lln(1.0 + NFIN / CDSCDN2));
        end

        if (CDSCDRN1 != 0.0) begin
            CDSCDR_i  =  CDSCDR_i * (1.0 + CDSCDRN1 / NFIN * lln(1.0 + NFIN / CDSCDRN2));
        end

        if (VSATN1 != 0.0) begin
            VSAT_i  =  VSAT_i * (1.0 + VSATN1 / NFIN * lln(1.0 + NFIN / VSATN2));
        end

        if (VSAT1N1 != 0.0) begin
            VSAT1_i  =  VSAT1_i * (1.0 + VSAT1N1 / NFIN * lln(1.0 + NFIN / VSAT1N2));
        end

        if (VSAT1RN1 != 0.0) begin
            VSAT1R_i  =  VSAT1R_i * (1.0 + VSAT1RN1 / NFIN * lln(1.0 + NFIN / VSAT1RN2));
        end

        if (U0N1 != 0.0) begin
            U0_i  =  U0_i * (1.0 + U0N1 / NFIN * lln(1.0 + NFIN / U0N2));
        end

        if ($param_given(NFINNOM)) begin
            PHIG_i  =  PHIG_i * (1.0 + (NFIN - NFINNOM) * PHIGLT * Leff) ;
            ETA0_i  =  ETA0_i * (1.0 + (NFIN - NFINNOM) * ETA0LT * Leff);
            U0_i    =  U0_i * (1.0 + (NFIN - NFINNOM) * U0LT * Leff);
        end

        if (U0N1R != 0.0) begin
            U0R_i  =  U0R_i * (1.0 + U0N1R / NFIN * lln(1.0 + NFIN / U0N2R));
        end

        // Length Scaling
        PHIG_i  =  PHIG_i + PHIGL * Leff;
        if (LPA > 0.0) begin
            U0_i  =  U0_i * (1.0 - UP_i * pow(Leff, -LPA));
        end else begin
            U0_i  =  U0_i * (1.0 - UP_i);
        end
        UA_i  =  UA_i + AUA * lexp(-Leff / BUA);
        UD_i  =  UD_i + AUD * lexp(-Leff / BUD);
        EU_i  =  EU_i + AEU * lexp(-Leff / BEU);

        if (ASYMMOD != 0) begin
            if (LPAR > 0.0) begin
                U0R_i  =  U0R_i * (1.0 - UPR_i * pow(Leff, -LPAR));
            end else begin
                U0R_i  =  U0R_i * (1.0 - UPR_i);
            end
            UAR_i  =  UAR_i + AUAR * lexp(-Leff / BUAR);
            UDR_i  =  UDR_i + AUDR * lexp(-Leff / BUDR);
            EUR_i  =  EUR_i + AEUR * lexp(-Leff / BEUR);
        end

        if (RDSMOD == 1) begin
            RSW_i  =  RSW_i + ARSW * lexp(-Leff / BRSW);
            RDW_i  =  RDW_i + ARDW * lexp(-Leff / BRDW);
        end else begin
            RDSW_i  =  RDSW_i + ARDSW * lexp(-Leff / BRDSW);
        end

        PCLM_i  =  PCLM_i + APCLM * lexp(-Leff / BPCLM);
        if (ASYMMOD != 0) begin
            PCLMR_i  =  PCLMR_i + APCLMR * pow(Leff, -BPCLMR);
        end

        MEXP_i  =  MEXP_i  + AMEXP  * pow(Leff, -BMEXP);
        if (ASYMMOD != 0) begin
            MEXPR_i  =  MEXPR_i + AMEXPR * pow(Leff, -BMEXPR);
        end

        PTWG_i  =  PTWG_i + APTWG * lexp(-Leff / BPTWG);
        if (ASYMMOD != 0) begin
            PTWGR_i  =  PTWGR_i + APTWG * lexp(-Leff / BPTWG);
        end

        VSAT_i    =  VSAT_i + AVSAT * lexp(-Leff / BVSAT);
        VSAT1_i   =  VSAT1_i + AVSAT1 * lexp(-Leff / BVSAT1);
        if (ASYMMOD != 0) begin
            VSAT1R_i =  VSAT1R_i + AVSAT1 * lexp(-Leff / BVSAT1);
        end

        PSAT_i    =  PSAT_i + APSAT * lexp(-Leff / BPSAT);
        PSATCV_i  =  PSATCV_i + APSATCV * lexp(-Leff / BPSATCV);
        VSATCV_i  =  VSATCV_i + AVSATCV * lexp(-LeffCV / BVSATCV);

        // Scaling for DITS Parameters
        DVTP0_i   =  DVTP0 + ADVTP0 * lexp(-Leff / BDVTP0);
        DVTP1_i   =  DVTP1 + ADVTP1 * lexp(-Leff / BDVTP1);

        // Geometrical Scaling for Toxeff / Charge Centroid Tcen
        if (QMTCENCV_i > 0.0 || QMTCENCVA_i > 0.0) begin
            MTcen  =  1.0 + AQMTCEN * lexp(- (2.0 * Ach / Weff_UFCM) / BQMTCEN);
            Tcen0  =  (2.0 * Ach / Weff_UFCM ) * MTcen;
        end

        // **************************************
        // *        Parameter Checking          *
        // **************************************

        if (LeffCV <= 1.0e-9) begin
            $strobe("Warning: LeffCV = %e <= 1.0e-9.", LeffCV);
        end

        if (BULKMOD != 0) begin
            if (LeffCV_acc <= 1.0e-9) begin
                $strobe("Warning: LeffCV_acc = %e <= 1.0e-9.", LeffCV_acc);
            end
        end

        if (Weff0 <= 1.0e-9) begin
            $strobe("Warning: Weff0 = %e <= 1.0e-9.", Weff0);
        end

        if (WeffCV0 <= 1.0e-9) begin
            $strobe("Warning: WeffCV0 = %e <= 1.0e-9.", WeffCV0);
        end

        if (NBODY_i <= 0.0) begin
            $strobe("Fatal: NBODY_i = %e is not positive.", NBODY_i);
            $finish(0);
        end else if (NBODY_i <= 1.0e18) begin
            $strobe("Warning: NBODY_i = %e m^-3 may be too small.", NBODY_i);
        end

        if (NGATE_i < 0.0) begin
            $strobe("Fatal: NGATE_i = %e is negative.", NGATE_i);
            $finish(0);
        end else if (NGATE_i != 0.0 && NGATE_i <= 1.0e24) begin
            $strobe("Warning: NGATE_i = %e may be too small.", NGATE_i);
        end else if (NGATE_i > 1.0e31) begin
            $strobe("Fatal: NGATE_i = %e is too high.", NGATE_i);
            $finish(0);
        end

        if (DVT0_i < 0.0) begin
            $strobe("Warning: DVT0_i = %e is negative.", DVT0_i);
        end

        if (PHIG_i <= 0.0) begin
            $strobe("Fatal: PHIG_i = %e is not positive.", PHIG_i);
            $finish(0);
        end

        if (VSAT_i <= 0.0) begin
            $strobe("Fatal: VSAT_i = %e is not positive.", VSAT_i);
            $finish(0);
        end

        if (VSAT1_i <= 0.0) begin
            $strobe("Fatal: VSAT1_i = %e is not positive.", VSAT1_i);
            $finish(0);
        end

        if (ASYMMOD != 0 && VSAT1R_i <= 0.0) begin
            $strobe("Fatal: VSAT1R_i = %e is not positive.", VSAT1R_i);
            $finish(0);
        end

        if (DVT1_i <= 0.0) begin
            $strobe("Fatal: DVT1_i = %e is not positive.", DVT1_i);
            $finish(0);
        end

        if (DVT1SS_i <= 0.0) begin
            $strobe("Fatal: DVT1SS_i = %e is not positive.", DVT1SS_i);
            $finish(0);
        end

        if (CDSC_i < 0.0) begin
            $strobe("Warning: CDSC_i = %e is negative.", CDSC_i);
        end

        if (CDSCD_i < 0.0) begin
            $strobe("Warning: CDSCD_i = %e is negative.", CDSCD_i);
        end

        if (ASYMMOD != 0 && CDSCDR_i < 0.0) begin
            $strobe("Warning: CDSCDR_i = %e is negative.", CDSCDR_i);
        end

        if (DSUB_i <= 0.0) begin
            $strobe("Fatal: DSUB_i = %e is not positive.", DSUB_i);
            $finish(0);
        end

        if (ETA0_i < 0.0) begin
            $strobe("Warning: ETA0_i = %e is negative, setting it to 0", ETA0_i);
            ETA0_i  =  0.0;
        end

        if (ETA0R_i < 0.0) begin
            $strobe("Warning: ETA0R_i = %e is negative, setting it to 0", ETA0R_i);
            ETA0R_i  =  0.0;
        end

        if (LPE0_i < -Leff) begin
            $strobe("Warning: LPE0_i = %e is less than -Leff. Clipping LPE0_i to 0", LPE0_i);
            LPE0_i  =  0.0;
        end

        if (K0SI_i <= 0.0) begin
            $strobe("Warning: K0SI_i = %e is not positive, setting it to 0.", K0SI_i);
            K0SI_i  =  0.0;
        end

        if (K2SI_i <= 0.0) begin
            $strobe("Warning: K2SI_i = %e is not positive, setting it to 0.", K2SI_i);
            K2SI_i  =  0.0;
        end

        if (PHIBE_i < 0.2 && BULKMOD != 0) begin
            $strobe("Warning: PHIBE_i = %e is less than 0.2, setting it to 0.2.", PHIBE_i);
            PHIBE_i  =  0.2;
        end

        if (PHIBE_i > 1.2 && BULKMOD != 0) begin
            $strobe("Warning: PHIBE_i = %e is larger than 1.2, setting it to 1.2.", PHIBE_i);
            PHIBE_i  =  1.2;
        end

        if (PSAT_i < 2.0) begin
            $strobe("Warning: PSAT_i = %e is less than 2.0, setting it to 2.0.", PSAT_i);
            PSAT_i  =  2.0;
        end

        if (PSATCV_i < 2.0) begin
            $strobe("Warning: PSATCV_i = %e is less than 2.0, setting it to 2.0.", PSATCV_i);
            PSATCV_i  =  2.0;
        end

        if (U0_i < 0.0) begin
            $strobe("Warning: U0_i = %e is negative, setting it to the default value.", U0_i);
            U0_i  =  0.03;
        end

        if (UA_i < 0.0) begin
            $strobe("Warning: UA_i = %e is negative, setting it to 0.", UA_i);
            UA_i  =  0.0;
        end

        if (EU_i < 0.0) begin
            $strobe("Warning: EU_i = %e is negative, setting it to 0.", EU_i);
            EU_i  =  0.0;
        end

        if (UD_i < 0.0) begin
            $strobe("Warning: UD_i = %e is negative, setting it to 0.", UD_i);
            UD_i  =  0.0;
        end

        if (UCS_i < 0.0) begin
            $strobe("Warning: UCS_i = %e is negative, setting it to 0.", UCS_i);
            UCS_i  =  0.0;
        end

        if (ETAMOB_i < 0.0) begin
            $strobe("Warning: ETAMOB_i = %e is negative, setting it to 0", ETAMOB_i);
            ETAMOB_i  = 0.0;
        end

        RDSWMIN_i  =  RDSWMIN;
        if (RDSWMIN_i < 0.0) begin
            $strobe("Warning: RDSWMIN = %e is negative.  Set to zero", RDSWMIN_i);
            RDSWMIN_i  =  0.0;
        end

        if (RDSW_i < 0.0) begin
            $strobe("Warning: RDSW_i = %e is negative.  Set to zero", RDSW_i);
            RDSW_i  =  0.0;
        end

        RSWMIN_i  =  RSWMIN;
        if (RSWMIN_i < 0.0) begin
            $strobe("Warning: RSWMIN = %e is negative.  Set to zero", RSWMIN_i);
            RSWMIN_i  =  0.0;
        end

        if (RSW_i < 0.0) begin
            $strobe("Warning: RSW_i = %e is negative.  Set to zero", RSW_i);
            RSW_i  =  0.0;
        end

        RDWMIN_i  =  RDWMIN;
        if (RDWMIN_i < 0.0) begin
            $strobe("Warning: RDWMIN = %e is negative.  Set to zero", RDWMIN_i);
            RDWMIN_i  =  0.0;
        end

        if (RDW_i < 0) begin
            $strobe("Warning: RDW_i = %e is negative.  Set to zero", RDW_i);
            RDW_i  =  0.0;
        end

        if (PRWGD_i < 0.0) begin
            $strobe("Warning: PRWGD_i = %e is negative.  Set to zero", PRWGD_i);
            PRWGD_i  =  0.0;
        end

        if (PRWGS_i < 0.0) begin
            $strobe("Warning: PRWGS_i = %e is negative.  Set to zero", PRWGS_i);
            PRWGS_i  =  0.0;
        end

        if (PCLM_i < 0) begin
            $strobe("Warning: PCLM_i = %e is negative.", PCLM_i);
        end

        if (PDIBL1_i < 0.0) begin
            $strobe("Warning: PDIBL1_i = %e is negative.", PDIBL1_i);
        end

        if (ASYMMOD != 0) begin
            if (PDIBL1R_i < 0.0) begin
                $strobe("Warning: PDIBL1R_i = %e is negative.", PDIBL1R_i);
            end
            if (PDIBL2R_i < 0.0) begin
                $strobe("Warning: PDIBL2R_i = %e is negative.", PDIBL2R_i);
            end
            if (U0R_i < 0) begin
                $strobe("Warning: U0R_i = %e is negative, setting it to 0.", U0R_i);
                U0R_i  =  0.0;
            end
            if (UAR_i < 0.0) begin
                $strobe("Warning: UAR_i = %e is negative, setting it to 0.", UAR_i);
                UAR_i  =  0.0;
            end
            if (EUR_i < 0.0) begin
                $strobe("Warning: EUR_i = %e is negative, setting it to 0.", EUR_i);
                EUR_i  =  0.0;
            end
            if (UDR_i < 0.0) begin
                $strobe("Warning: UDR_i = %e is negative, setting it to 0.", UDR_i);
                UDR_i  =  0.0;
            end
        end

        if (PDIBL2_i < 0.0) begin
            $strobe("Warning: PDIBL2_i = %e is negative.", PDIBL2_i);
        end

        if (DROUT_i <= 0.0) begin
            $strobe("Fatal: DROUT_i = %e is non-positive.", DROUT_i);
            $finish(0);
        end

        if (MEXP_i < 2.0) begin
            $strobe("Warning: MEXP_i = %e < 2. Setting MEXP_i = 2.", MEXP_i);
            MEXP_i  =  2.0;
        end

        if (ASYMMOD != 0) begin
            if (MEXPR_i < 2.0) begin
                $strobe("Warning: MEXPR_i = %e < 2.  Setting MEXPR_i = 2.", MEXPR_i);
                MEXPR_i  =  2.0;
            end
        end

        if (PTWG_i < 0) begin
            $strobe("Warning: PTWG_i = %e is negative, setting it to 0.", PTWG_i );
            PTWG_i  =  0.0;
        end

        if (QMTCENCV_i > 0.0) begin
            if (QM0 <= 0.0) begin
                $strobe("Fatal: QM0 = %e is non-positive.", QM0);
                $finish(0);
            end
        end

        if (BULKMOD != 0 && QMTCENCVA_i > 0.0) begin
            if (QM0ACC <= 0.0) begin
                $strobe("Fatal: QM0ACC = %e is non-positive.", QM0ACC);
                $finish(0);
            end
        end

        if (CGIDL_i < 0.0) begin
            $strobe("Warning: CGIDL_i = %e < 0.  Setting CGIDL_i = 0.", CGIDL_i);
            CGIDL_i  =  0.0;
        end

        if (CGISL_i < 0.0) begin
            $strobe("Warning: CGISL_i = %e < 0.  Setting CGISL_i = 0.", CGISL_i);
            CGISL_i  =  0.0;
        end

        if (IGBMOD != 0) begin
            if (NIGBINV_i <= 0.0) begin
                $strobe("Fatal: NIGBINV_i = %e is non-positive.", NIGBINV_i);
                $finish(0);
            end
            if (NIGBACC_i <= 0.0) begin
                $strobe("Fatal: NIGBACC_i = %e is non-positive.", NIGBACC_i);
                $finish(0);
            end
        end

        if (IGCMOD != 0) begin
            if (POXEDGE_i <= 0.0) begin
                $strobe("Fatal: POXEDGE_i = %e is non-positive.", POXEDGE_i);
                $finish(0);
            end
            if (PIGCD_i <= 0.0) begin
                $strobe("Fatal: PIGCD_i = %e is non-positive.", PIGCD_i);
                $finish(0);
            end
        end

        if (IGCMOD != 0 || IGBMOD != 0) begin
            if (TOXREF <= 0) begin
                $strobe("Fatal: TOXREF = %e is non-positive.", TOXREF);
                $finish(0);
            end
        end

        if (LINTIGEN >= (Leff / 2.0)) begin
            $strobe("Warning: LINTIGEN = %e is too large - Leff for r/g current is negative.  Re-setting LINTIGEN = 0.", LINTIGEN);
            LINTIGEN_i  =  0.0;
        end else begin
            LINTIGEN_i  =  LINTIGEN;
        end

        if (NTGEN_i <= 0.0) begin
            $strobe("Fatal: NTGEN_i = %e is non-positive.", NTGEN_i);
            $finish(0);
        end

        `ifdef __NQSMOD1__
            if (NQSMOD == 1 && XRCRG1_i != 0.0 && XRCRG1_i < 1.0e-3) begin
                $strobe("Warning: XRCRG1_i = %e.  Gate resistance may be too large. Disabling NQS Gate Resistance.", XRCRG1_i);
                XRCRG1_i  =  0.0;
            end
        `endif

        if (IIMOD == 2) begin
            if (BETAII0_i < 0.0) begin
                $strobe("Warning:  BETAII0_i = %e is negative.",BETAII0_i);
            end
            if (BETAII1_i < 0.0) begin
                $strobe("Warning:  BETAII1_i = %e is negative.", BETAII1_i);
            end
            if (BETAII2_i < 0.0) begin
                $strobe("Warning:  BETAII2_i = %e is negative.", BETAII2_i);
            end
            if (ESATII_i < 0.0) begin
                $strobe("Warning:  ESATII_i = %e is negative.", ESATII_i);
            end
            if (LII_i < 0.0) begin
                $strobe("Warning:  LII_i = %e is negative.", LII_i);
            end
            if (SII1_i < 0.0) begin
                $strobe("Warning:  SII1_i = %e is negative.", SII1);
            end
            if (SII2_i < 0.0) begin
                $strobe("Warning:  SII2_i = %e is negative.", SII2_i);
            end
            if (SIID_i < 0.0) begin
                $strobe("Warning:  SIID_i = %e is negative.", SIID_i);
            end
        end

        if (EF <= 0.0) begin
            $strobe("Fatal: EF = %e is non-positive.", EF);
            $finish(0);
        end else if (EF > 2.0) begin
            $strobe("Fatal: EF = %e > 2.0.", EF);
            $finish(0);
        end

        if (LINTNOI >= (Leff / 2.0)) begin
            $strobe("Warning: LINTNOI = %e is too large - Leff for noise is negative.  Re-setting LINTNOI = 0.", LINTNOI);
            LINTNOI_i  =  0.0;
        end else begin
            LINTNOI_i  =  LINTNOI;
        end

        if (NTNOI < 0) begin
            $strobe("Warning: NTNOI = %e is negative. Set to zero.", NTNOI);
            NTNOI_i  =  0.0;
        end else begin
            NTNOI_i  =  NTNOI;
        end

        // Self-Heating
        `ifdef __SHMOD__
            if (SHMOD != 0 && RTH0 > 0.0) begin
                T1y  =  WTH0;
                if (BSHEXP != 0.0)  begin
                    T1y  =  WTH0 * pow(NF, BSHEXP);
                end
                T2y  =  FPITCH;
                if (ASHEXP != 0.0) begin
                    T2y = FPITCH * pow(NFINtotal, ASHEXP);
                end
                gth  =  (T1y + T2y)/ RTH0;
                cth  =  CTH0 * (T1y + T2y);
            end else begin
                gth  =  1.0;
                cth  =  0.0;
            end
        `else
            if (SHMOD != 0) begin
                $strobe("Although the model selector SHMOD is set to 1, the self heating model is not activated in the Verilog-A code.  Please uncomment \"`define __SHMOD__\" in bsimcmg.va to activate it.");
            end
        `endif

        // Gate Electrode Resistance
        `ifdef __RGATEMOD__
            if (RGATEMOD != 0) begin
                Rgeltd  =  (RGEXT / NGCON + (RGFIN * NFIN) / (NGCON == 2 ? 12.0 : 3.0)) / NF;
                ggeltd  =  1.0 / max(1.0e-3, Rgeltd);
            end
        `else
            if (RGATEMOD != 0)
                $strobe("Although the model selector RGATEMOD is set to 1, the gate electrode resistance model is not activated in the Verilog-A code.  Please uncomment \"`define __RGATEMOD__\" in bsimcmg.va to activate it.");
        `endif

        // Geometry-Dependent Source/Drain Resistance
        if (RGEOMOD == 0) begin
            RSourceGeo  =  RSHS * NRS;
            RDrainGeo   =  RSHD * NRD;
        end else begin
            // Area and perimeter calculation
            if (HEPI > 0.0) begin
                Arsd  =  FPITCH * HFIN + (TFIN + (FPITCH - TFIN) * CRATIO) * HEPI;
            end else begin
                Arsd  =  FPITCH * max(1.0e-9, HFIN + HEPI);
            end
            Prsd  =  FPITCH + DELTAPRSD;

            // Resistivity Calculation
            if ($param_given(RHORSD)) begin
                rhorsd  =  RHORSD;
            end else begin
                mu_max  =  (TYPE == `ntype) ? 1417.0 : 470.5;
                if (TYPE == `ntype) begin
                    mu_rsd  =  (52.2 + (mu_max - 52.2) / (1.0 + pow(NSD / 9.68e22, 0.680)) - 43.4 / (1.0 + pow(3.43e26 / NSD, 2.0))) * 1.0e-4;
                end else begin
                    mu_rsd  =  (44.9 + (mu_max - 44.9) / (1.0 + pow(NSD / 2.23e22, 0.719)) - 29.0 / (1.0 + pow(6.10e26 / NSD, 2.0))) * 1.0e-4;
                end
                rhorsd  =  1.0 / (`q * NSD * mu_rsd);
            end

            // Component: Spreading Resistance (extension -> hdd)
            thetarsp  =  55.0 * `M_PI / 180.0;
            afin      =  min(Arsd, max(1.0e-18, TFIN * (HFIN + min(0.0, HEPI))));
            T1y        =  `COT(thetarsp);
            Rsp       =  rhorsd * T1y / (sqrt(`M_PI) * NFIN) * (1.0 / sqrt(afin) - 2.0 / sqrt(Arsd) + sqrt(afin / (Arsd*Arsd)));

            // Component: Contact Resistance
            arsd_total  =  Arsd * NFIN + ARSDEND;
            prsd_total  =  Prsd * NFIN + PRSDEND;
            lt          =  sqrt(RHOC * arsd_total / (rhorsd * prsd_total));
            alpha       =  LRSD / lt;
            T0y          =  lexp(alpha + alpha);

            if (SDTERM == 1.0) begin
                eta  =  rhorsd * lt / RHOC;
                T1y   =  T0y * (1.0 + eta);
                T2y  =  T1y + 1.0 - eta;
                T3y  =  T1y - 1.0 + eta;
            end else begin
                T2y  =  T0y + 1.0;
                T3y  =  T0y - 1.0;
            end
            RrsdTML  =  rhorsd * lt * T2y / (arsd_total * T3y);

            if (HEPI < -1.0e-10) begin
                Rrsdside  =  RHOC / (-HEPI * TFIN * NFIN);
                Rrsd      =  (RrsdTML + Rsp) * Rrsdside / ((RrsdTML + Rsp) + Rrsdside);
            end else begin
                Rrsd      =  RrsdTML + Rsp;
            end

            Rdsgeo      =  Rrsd / NF * max(0.0, RGEOA + RGEOB * TFIN + RGEOC * FPITCH + RGEOD * LRSD + RGEOE * HEPI);
            RSourceGeo  =  Rdsgeo;
            RDrainGeo   =  Rdsgeo;
        end

        // Clamping of Source/Drain Resistances
        if (RSourceGeo <= 1.0e-3) begin
            RSourceGeo  =  1.0e-3;
        end

        if (RDrainGeo <= 1.0e-3) begin
            RDrainGeo  =  1.0e-3;
        end

        if (RDSMOD == 1) begin
            if (RSWMIN_i <= 0.0) begin
                RSWMIN_i  =  0.0;
            end
            if (RDWMIN_i <= 0.0) begin
                RDWMIN_i  =  0.0;
            end
            if (RSW_i <= 0.0) begin
                RSW_i  =  0.0;
            end
            if (RDW_i <= 0.0) begin
                RDW_i  =  0.0;
            end
        end else begin
            if (RDSWMIN_i <= 0.0) begin
                RDSWMIN_i  =  0.0;
            end
            if (RDSW_i <= 0.0) begin
                RDSW_i  =  0.0;
            end
        end  // End of Clamping of Source/Drain Resistances

        if (CGEOMOD != 1) begin
            if ($param_given(CGSO)) begin
                CGSO_i  =  CGSO;
            end else begin
                if ($param_given(DLC) && DLC > 0.0) begin
                    CGSO_i  =  max(0.0, DLC * cox - CGSL_i);
                end else begin
                    CGSO_i  =  0.3 * TFIN * cox;
                end
            end
            if ($param_given(CGDO)) begin
                CGDO_i  =  CGDO;
            end else begin
                if ($param_given(DLC) && DLC > 0.0) begin
                    CGDO_i  =  max(0.0, DLC * cox - CGDL_i);
                end else begin
                    CGDO_i  =  0.3 * TFIN * cox;
                end
            end
        end

        // Parasitic Source/Drain to Gate Fringe Capacitance Model
        if (CGEOMOD == 2) begin
            if ($param_given(LSP))
                LSP_i = LSP;
            else
                LSP_i = 0.2*(L + XL);
            Hg    =  TGATE + TMASK;
            Trsd  =  0.5 * (FPITCH - TFIN);
            Wg    =  max(0.0, Trsd - TOXP);
            Hrsd  =  max(0.0, HEPI + TSILI);

            // Top Component
            if (TMASK > 0.0) begin
                // Capacitance Model by Chung-Hsun Lin (IBM)
                T0y       =  3.467e-11 * lln(1.0e-7 * EPSRSP / (3.9 * LSP_i));
                T1y       =  0.942 * Hrsd * epssp / LSP_i;
                Cgg_top  =  (T0y + T1y) * (TFIN + (FPITCH - TFIN) * CRATIO);
            end else begin
                `Cfringe_2d(cfr_top_trigate, Hg, Hrsd, LSP_i, TFIN, LRSD, Lg, TOXP, 0.85, Cgg_top)
            end

            // Side Component
            if (TMASK > 0) begin
                `Cfringe_2d(cfr_side_dblgate, Wg, Trsd, LSP_i, HFIN, LRSD, Lg, TOXP, 0.70, Cgg_side)
            end else begin
                `Cfringe_2d(cfr_side_trigate, Wg, Trsd, LSP_i, HFIN, LRSD, Lg, TOXP, 0.85, Cgg_side)
            end

            // Corner Component
            if (TMASK > 0.0) begin
                Acorner  =  0.0;
            end else begin
                if (HEPI > 0.0) begin
                    Acorner  =  (FPITCH - TFIN) * (HEPI * CRATIO + TSILI);
                end else begin
                    Acorner  =  (FPITCH - TFIN) * Hrsd;
                end
            end
            Ccorner  =  (NFIN * Acorner + ARSDEND + ASILIEND) * epssp / LSP_i;
            Cfr_geo  =  (Ccorner + Cgg_top * NFIN + CGEOE * Cgg_side * NFIN * 2.0) * NF;
            Cfr_geo  =  Cfr_geo * max(0.0, CGEOA + CGEOB * TFIN + CGEOC * FPITCH + CGEOD * LRSD);
        end

        // Source/Gate/Drain-to-Substrate Parasitic Capacitances
        T0y     =  CSDESW * lln(1.0 + HFIN / EOTBOX);
        csbox  =  cbox * ASEO + T0y * max(0.0, PSEO - FPITCH * NFINtotal);
        cdbox  =  cbox * ADEO + T0y * max(0.0, PDEO - FPITCH * NFINtotal);
        cgbox  =  (CGBO * NF * NGCON + CGBN * NFINtotal) * Lg;

        // Mobility Degradation
        EeffFactor    =  1.0e-8 / (epsratio * (EOT));
        WeffWRFactor  =  1.0 / (pow((Weff0) * 1.0e6, WR_i) * NFINtotal);
        litl          =  sqrt(epsratio * EOT * 0.5 * TFIN);

        if (!$param_given(THETASCE)) begin
            tmp  =  DVT1_i * Leff / scl + 1.0e-6;
            if (tmp < 40.0) begin
                Theta_SCE  =  0.5 / (cosh(tmp) - 1.0);
            end else begin
                Theta_SCE  =  exp(-tmp);
            end
        end else begin
            Theta_SCE  =  THETASCE;
        end

        if (!$param_given(THETASW)) begin
            tmp  =  DVT1SS_i * Leff / scl + 1.0e-6;
            if (tmp < 40.0) begin
                Theta_SW  =  0.5 / (cosh(tmp) - 1.0);
            end else begin
                Theta_SW  =  exp(-tmp);
            end
        end else begin
            Theta_SW  =  THETASW;
        end

        if (!$param_given(THETADIBL)) begin
            tmp  =  DSUB_i * Leff / scl + 1.0e-6;
            if (tmp < 40.0) begin
                Theta_DIBL  =  0.5 / (cosh(tmp) - 1.0);
            end else begin
                Theta_DIBL  =  exp(-tmp);
            end
        end else begin
            Theta_DIBL  =  THETADIBL;
        end

        Theta_RSCE  =  sqrt(1.0 + LPE0_i / Leff) - 1.0;

        tmp  =  DSUB_i * Leff / scl + 1.0e-6;
        if (tmp < 40.0) begin
            T0y  =  1.0 / max((1.0 + DVTP2 * (cosh(tmp) - 2.0)), 1.0e-6);
        end else begin
            T0y  =  exp(-tmp) / max((exp(-tmp) + DVTP2), 1.0e-6);
        end

        Theta_DITS  =  T0y;
        nbody       =  NBODY_i;
        qbs         =  `q * nbody * Ach / Cins;

        // Gate Current
        if (TYPE == `ntype) begin
            Aechvb  =  4.97232e-7;  // NMOS
            Bechvb  =  7.45669e11;  // NMOS
        end else begin
            Aechvb  =  3.42537e-7;  // PMOS
            Bechvb  =  1.16645e12;  // PMOS
        end

        T0y  =  TOXG * TOXG;
        T1y  =  TOXG * POXEDGE_i;
        T2y =  T1y * T1y;
        Toxratio      =  lexp(NTOX_i * lln(TOXREF / TOXG)) / T0y;
        Toxratioedge  =  lexp(NTOX_i * lln(TOXREF / T1y)) / T2y;
        igsd_mult0    =  Weff0 * Aechvb * Toxratioedge;

        if (TNOM < -`P_CELSIUS0) begin
            $strobe("Warning: (TNOM=%e) < -`P_CELSIUS0. Set to 27 C.", TNOM);
            Tnom  =  `REFTEMP;
        end else begin
            Tnom  =  TNOM + `CONSTCtoK;
        end
    end // initial_step

    // ************************************************
    // *     Temperature Dependence Calculations      *
    // ************************************************
    `ifdef __SHMOD__
        if (SHMOD != 0 && RTH0 > 0.0) begin
            DevTemp  =  $temperature + Temp(rth_branch) + DTEMP;
        end else begin
        `endif
        DevTemp  =  $temperature + DTEMP;
        `ifdef __SHMOD__
        end
    `endif

    begin : CMGTempDepCalc
        TRatio   =  DevTemp / Tnom;
        delTemp  =  DevTemp - Tnom;
        Vtm      =  `KboQ * DevTemp;
        Vtm0     =  `KboQ * Tnom;
        Eg       =  BG0SUB - TBGASUB * DevTemp * DevTemp / (DevTemp + TBGBSUB);
        Eg0      =  BG0SUB - TBGASUB * Tnom * Tnom / (Tnom + TBGBSUB);
        T1       =  (DevTemp / 300.15) * sqrt(DevTemp / 300.15);
        ni       =  NI0SUB * T1 * lexp(BG0SUB / (2.0 * `KboQ * 300.15) - Eg / (2.0 * Vtm));
        Nc       =  NC0SUB * T1;
        ThetaSS  =  hypsmooth(1.0 + TSS_i * delTemp - 1.0e-6, 1.0e-3);

        // Quantum Mechanical Vth Correction (Ref: Trivedi et al., EDL 2005)
        kT       =  Vtm * `q;
        T0y       =  `HBAR * `M_PI / (2*Ach/Weff_UFCM );
        E0       =  T0y * T0y / (2.0 * mx);
        E0prime  =  T0y * T0y / (2.0 * mxprime);
        E1       =  4.0 * E0;
        E1prime  =  4.0 * E0prime;
        T1       =  gprime * mdprime / (gfactor * md);
        gam0     =  1.0 + T1 * lexp((E0 - E0prime) / kT);
        gam1     =  gam0 + lexp((E0 - E1) / kT) + T1 * lexp((E0 - E1prime) / kT);
        T2       =  -Vtm * lln(gfactor * md / (`M_PI * `HBAR * `HBAR * Nc) * kT / (2.0 * Ach / Weff_UFCM) * gam1);
        dvch_qm  =  QMFACTOR_i * (E0 / `q + T2);

        // Temperature Dependence
        ETA0_t   =  Tempdep(ETA0_i, TETA0, delTemp, TEMPMOD);
        ETA0R_t  =  Tempdep(ETA0R_i, TETA0R, delTemp, TEMPMOD);
        T1       =  U0_i * pow(TRatio, UTE_i);
        U0_t     =  T1 + hypmax(UTL_i * delTemp, -0.9 * T1, 1.0e-4);
        u0       =  U0_t;
        if (ASYMMOD == 1) begin
            T1     =  U0R_i * pow(TRatio, UTER_i);
            U0R_t  =  T1 + hypmax(UTLR_i * delTemp, -0.9 * T1, 1.0e-4);
            u0r    =  U0R_t;
        end

        ETAMOB_t  =  Tempdep(ETAMOB_i, EMOBT_i, delTemp, TEMPMOD);
        UA_t      =  UA_i + hypmax(UA1_i*delTemp, -UA_i, 1.0e-6);
        if (ASYMMOD != 0) begin
            UAR_t  =  UAR_i + hypmax(UA1R_i * delTemp, -UAR_i, 1.0e-6);
        end

        if (BULKMOD != 0) begin
            if (TEMPMOD == 0) begin
                UC_t  =  Tempdep(UC_i, UC1_i, delTemp, 0);
                if (ASYMMOD != 0) begin
                    UCR_t  =  Tempdep(UCR_i, UC1R_i, delTemp, 0);
                end
            end else begin
                UC_t  =  UC_i + UC1_i * delTemp;
                if (ASYMMOD != 0) begin
                    UCR_t  =  UCR_i + UC1R_i * delTemp;
                end
            end
        end

        UD_t  =  UD_i * pow(TRatio, UD1_i);
        if (ASYMMOD != 0) begin
            UDR_t  =  UDR_i * pow(TRatio, UD1R_i);
        end

        UCS_t  =  UCS_i * pow(TRatio, UCSTE_i);

        rdstemp  =  hypsmooth(1.0 + PRT_i * delTemp - 1.0e-6, 1.0e-3);
        RSDR_t   =  Tempdep(RSDR, TRSDR, delTemp, TEMPMOD);
        if (ASYMMOD != 0) begin
            RSDRR_t  =  Tempdep(RSDRR, TRSDR, delTemp, TEMPMOD);
        end

        RDDR_t  = Tempdep(RDDR, TRDDR, delTemp, TEMPMOD);
        if (ASYMMOD != 0) begin
            RDDRR_t  =  Tempdep(RDDRR, TRDDR, delTemp, TEMPMOD);
        end

        VSAT_t  =  Tempdep(VSAT_i, -AT_i, delTemp, TEMPMOD);
        if (VSAT_t < 1000) begin
            $strobe("Warning: VSAT(%f) = %e is less than 1K, setting it to 1K.", DevTemp, VSAT_t);
            VSAT_t =  1000;
        end

        if (ASYMMOD != 0) begin
            VSATR_t   = Tempdep(VSATR_i, -ATR_i, delTemp, TEMPMOD);
            if (VSATR_t < 1000) begin
                $strobe("Warning: VSATR(%f) = %e is less than 1K, setting it to 1K.", DevTemp, VSATR_t);
                VSATR_t =  1000;
            end
        end

        VSAT1_t  =  Tempdep(VSAT1_i, -AT_i, delTemp, TEMPMOD);
        if (VSAT1_t < 1000) begin
            $strobe("Warning: VSAT1(%f) = %e is less than 1K, setting it to 1K.", DevTemp, VSAT1_t);
            VSAT1_t  =  1000;
        end

        if (ASYMMOD != 0) begin
            VSAT1R_t   = Tempdep(VSAT1R_i, -AT_i, delTemp, TEMPMOD);
            if (VSAT1R_t < 1000) begin
                $strobe("Warning: VSAT1R(%f) = %e is less than 1K, setting it to 1K.", DevTemp, VSAT1R_t);
                VSAT1R_t  =  1000;
            end
        end

        VSATCV_t  =  Tempdep(VSATCV_i, -ATCV_i, delTemp, TEMPMOD);
        if (VSATCV_t < 1000) begin
            $strobe("Warning: VSATCV(%f) = %e is less than 1K, setting it to 1K.", DevTemp, VSATCV_t);
            VSATCV_t  =  1000;
        end

        MEXP_t  =  hypsmooth(MEXP_i * (1.0 + TMEXP * delTemp) - 2.0, 1.0e-3) + 2.0;
        if (ASYMMOD != 0) begin
            MEXPR_t  =  hypsmooth(MEXPR_i * (1.0 + TMEXPR * delTemp) - 2.0, 1.0e-3) + 2.0;
        end

        PTWG_t  =  Tempdep(PTWG_i, -PTWGT_i, delTemp, TEMPMOD);
        if (ASYMMOD != 0) begin
            PTWGR_t  =  Tempdep(PTWGR_i, -PTWGT_i, delTemp, TEMPMOD);
        end

        dvth_temp   =  (KT1_i + KT1L / Leff) * (TRatio - 1.0);
        BETA0_t     =  BETA0_i * pow(TRatio, IIT_i);
        SII0_t      =  SII0_i * (hypsmooth(1.0 + TII_i * (TRatio - 1.0) - 0.01, 1.0e-3) + 0.01);

        K0_t        =  K0_i + K01_i * delTemp;
        K0SI_t      =  K0SI_i + hypmax(K0SI1_i * delTemp, -K0SI_i, 1.0e-6);
        K2SI_t      =  K2SI_i + hypmax(K2SI1_i * delTemp, -K2SI_i, 1.0e-6);
        K1_t        =  K1_i + hypmax(K11_i * delTemp, -K1_i, 1.0e-6);
        K2SAT_t     =  K2SAT_i + K2SAT1_i * delTemp;
        A1_t        =  A1_i + A11_i * delTemp;
        A2_t        =  A2_i + A21_i * delTemp;
        K2_t        =  K2_i + hypmax(K21_i * delTemp, -K2_i, 1.0e-6);
        K0SISAT_t   =  K0SISAT_i + K0SISAT1_i * delTemp;
        K2SISAT_t   =  K2SISAT_i + K2SISAT1_i * delTemp;
        AIGBINV_t   =  AIGBINV_i + hypmax(AIGBINV1_i * delTemp, -AIGBINV_i, 1.0e-6);
        AIGBACC_t   =  AIGBACC_i + hypmax(AIGBACC1_i * delTemp, -AIGBACC_i, 1.0e-6);
        AIGC_t      =  AIGC_i + hypmax(AIGC1_i * delTemp, -AIGC_i, 1.0e-6);
        AIGS_t      =  AIGS_i + hypmax(AIGS1_i * delTemp, -AIGS_i, 1.0e-6);
        AIGD_t      =  AIGD_i + hypmax(AIGD1_i * delTemp, -AIGD_i, 1.0e-6);
        BGIDL_t     =  BGIDL_i * hypsmooth(1.0 + TGIDL_i * delTemp - 1.0e-6, 1.0e-3);
        BGISL_t     =  BGISL_i * hypsmooth(1.0 + TGIDL_i * delTemp - 1.0e-6, 1.0e-3);
        ALPHA0_t    =  ALPHA0_i + hypmax(ALPHA01 * delTemp, -ALPHA0_i, 1.0e-6);
        ALPHA1_t    =  ALPHA1_i + hypmax(ALPHA11 * delTemp, -ALPHA1_i, 1.0e-6);
        ALPHAII0_t  =  ALPHAII0_i + hypmax(ALPHAII01 * delTemp, -ALPHAII0_i, 1.0e-25);
        ALPHAII1_t  =  ALPHAII1_i + hypmax(ALPHAII11 * delTemp, -ALPHAII1_i, 1.0e-20);
        igtemp      =  lexp(IGT_i * lln(TRatio));
        igsd_mult   =  igsd_mult0 * igtemp;

        if (BULKMOD != 0) begin
            CJS_t     =  Tempdep(CJS, TCJ, delTemp, TEMPMOD);
            CJD_t     =  Tempdep(CJD, TCJ, delTemp, TEMPMOD);
            CJSWS_t   =  Tempdep(CJSWS, TCJSW, delTemp, TEMPMOD);
            CJSWD_t   =  Tempdep(CJSWD, TCJSW, delTemp, TEMPMOD);
            CJSWGS_t  =  Tempdep(CJSWGS, TCJSWG, delTemp, TEMPMOD);
            CJSWGD_t  =  Tempdep(CJSWGD, TCJSWG, delTemp, TEMPMOD);

            PBS_t     =  hypsmooth(PBS - TPB * delTemp - 0.01, 1.0e-3) + 0.01;
            PBD_t     =  hypsmooth(PBD - TPB * delTemp - 0.01, 1.0e-3) + 0.01;
            PBSWS_t   =  hypsmooth(PBSWS - TPBSW * delTemp - 0.01, 1.0e-3) + 0.01;
            PBSWD_t   =  hypsmooth(PBSWD - TPBSW * delTemp - 0.01, 1.0e-3) + 0.01;
            PBSWGS_t  =  hypsmooth(PBSWGS - TPBSWG * delTemp - 0.01, 1.0e-3) + 0.01;
            PBSWGD_t  =  hypsmooth(PBSWGD - TPBSWG * delTemp - 0.01, 1.0e-3) + 0.01;

            T0        =  Eg0 / Vtm0 - Eg / Vtm;
            T1        =  lln(TRatio);
            T3        =  lexp((T0 + XTIS * T1) / NJS);
            JSS_t     =  JSS * T3;
            JSWS_t    =  JSWS * T3;
            JSWGS_t   =  JSWGS * T3;

            T3        =  lexp((T0 + XTID * T1) / NJD);
            JSD_t     =  JSD * T3;
            JSWD_t    =  JSWD * T3;
            JSWGD_t   =  JSWGD * T3;

            JTSS_t     =  JTSS * lexp(Eg0 * XTSS * (TRatio - 1.0) / Vtm);
            JTSD_t     =  JTSD * lexp(Eg0 * XTSD * (TRatio - 1.0) / Vtm);
            JTSSWS_t   =  JTSSWS * lexp(Eg0 * XTSSWS * (TRatio - 1.0) / Vtm);
            JTSSWD_t   =  JTSSWD * lexp(Eg0 * XTSSWD * (TRatio - 1.0) / Vtm);
            JTSSWGS_t  =  JTSSWGS * (sqrt(JTWEFF / Weff0) + 1.0) * lexp(Eg0 * XTSSWGS * (TRatio - 1.0) / Vtm);
            JTSSWGD_t  =  JTSSWGD * (sqrt(JTWEFF / Weff0) + 1.0) * lexp(Eg0 * XTSSWGD * (TRatio - 1.0) / Vtm);

            // All NJT's Smoothed to 0.01 to Prevent Divide-by-zero / Negative Values
            NJTS_t      =  hypsmooth(NJTS * (1.0 + TNJTS * (TRatio-1.0)) - 0.01, 1.0e-3) + 0.01;
            NJTSD_t     =  hypsmooth(NJTSD * (1.0 + TNJTSD * (TRatio-1.0)) - 0.01, 1.0e-3) + 0.01;
            NJTSSW_t    =  hypsmooth(NJTSSW * (1.0 + TNJTSSW * (TRatio-1.0)) - 0.01, 1.0e-3) + 0.01;
            NJTSSWD_t   =  hypsmooth(NJTSSWD * (1.0 + TNJTSSWD * (TRatio-1.0)) - 0.01, 1.0e-3) + 0.01;
            NJTSSWG_t   =  hypsmooth(NJTSSWG * (1.0 + TNJTSSWG * (TRatio-1.0)) - 0.01, 1.0e-3) + 0.01;
            NJTSSWGD_t  =  hypsmooth(NJTSSWGD * (1.0 + TNJTSSWGD * (TRatio-1.0)) - 0.01, 1.0e-3) + 0.01;
        end

        if (!$param_given(VFBSD)) begin
            if (NGATE > 0.0) begin
                vfbsd  =  devsign * (hypsmooth(0.5 * Eg - Vtm * lln(NGATE / ni), 1.0e-4) - (0.5 * Eg - devsign * (0.5 * Eg - hypsmooth(0.5 * Eg - Vtm * lln(NSD / ni), 1.0e-4))));
            end else begin
                vfbsd  =  devsign * (PHIG_i - (EASUB + 0.5 * Eg - devsign * (0.5 * Eg - hypsmooth(0.5 * Eg - Vtm * lln(NSD / ni), 1.0e-4))));
            end
        end else begin
            vfbsd  =  VFBSD;
        end

        if (!$param_given(VFBSDCV)) begin
            vfbsdcv  =  vfbsd;
        end else begin
            vfbsdcv  =  VFBSDCV;
        end

        `ifdef __SHMOD__
            if (SHMOD != 0 && RTH0 > 0.0) begin
                T0    =  Vtm * lln(nbody / ni);
                phib  =  sqrt(T0 * T0 + 1.0e-6);
            end else begin
                phib  =  Vtm * lln(nbody / ni);
            end
        `else
            phib  =  Vtm * lln(nbody/ni);
        `endif

        `ifdef __SHMOD__
            if (SHMOD != 0 && RTH0 > 0.0) begin
                T0   =  Vtm * lln(nbody * NSD / (ni * ni));
                vbi  =  sqrt(T0 * T0 + 1.0e-6);
            end else begin
                vbi  =  Vtm * lln(nbody * NSD / (ni * ni));
            end
        `else
            vbi  =  Vtm * lln(nbody * NSD / (ni * ni));
        `endif

        // deltaPhi definition and Polysilicon Depletion
        // deltaPhi: workfunction difference between the gate and the n+ source.
        deltaPhi  =  devsign*(PHIG_i - (EASUB + (TYPE == `ntype ? 0 : Eg)));

        // Mobility Degradation
        eta_mu     =  0.5 * ETAMOB_t;
        eta_mu_cv  =  0.5;
        if ( TYPE != `ntype ) begin
            eta_mu     =  1.0 / 3.0 * ETAMOB_t;
            eta_mu_cv  =  1.0 / 3.0;
        end

        // Junction Current and Capacitance
        if (BULKMOD != 0) begin
            // Source-Side Junction Current
            Isbs  =  ASEJ * JSS_t + PSEJ * JSWS_t + TFIN * NFINtotal * JSWGS_t;
            if (Isbs > 0.0) begin
                Nvtms     =  Vtm * NJS;
                XExpBVS   =  lexp(-BVS / Nvtms) * XJBVS;
                T2        =  max(IJTHSFWD / Isbs, 10.0);
                Tb        =  1.0 + T2 - XExpBVS;
                VjsmFwd   =  Nvtms * lln(0.5 * (Tb + sqrt(Tb * Tb + 4.0 * XExpBVS)));
                T0        =  lexp(VjsmFwd / Nvtms);
                IVjsmFwd  =  Isbs * (T0 - XExpBVS / T0 + XExpBVS - 1.0);
                SslpFwd   =  Isbs * (T0 + XExpBVS / T0) / Nvtms;
                T2        =  hypsmooth(IJTHSREV / Isbs - 10.0, 1.0e-3) + 10.0;
                VjsmRev   =  -BVS - Nvtms * lln((T2 - 1.0) / XJBVS);
                T1        =  XJBVS * lexp(-(BVS + VjsmRev) / Nvtms);
                IVjsmRev  =  Isbs * (1.0 + T1);
                SslpRev   =  -Isbs * T1 / Nvtms;
            end

            // Drain-Side Junction Current
            Isbd  =  ADEJ * JSD_t + PDEJ * JSWD_t + TFIN * NFINtotal * JSWGD_t;
            if (Isbd > 0.0) begin
                Nvtmd     =  Vtm * NJD;
                XExpBVD   =  lexp(-BVD / Nvtmd) * XJBVD;
                T2        =  max(IJTHDFWD / Isbd, 10.0);
                Tb        =  1.0 + T2 - XExpBVD;
                VjdmFwd   =  Nvtmd * lln(0.5 * (Tb + sqrt(Tb * Tb + 4.0 * XExpBVD)));
                T0        =  lexp(VjdmFwd / Nvtmd);
                IVjdmFwd  =  Isbd * (T0 - XExpBVD / T0 + XExpBVD - 1.0);
                DslpFwd   =  Isbd * (T0 + XExpBVD / T0) / Nvtmd;
                T2        =  hypsmooth(IJTHDREV / Isbd - 10.0, 1.0e-3) + 10.0;
                VjdmRev   =  -BVD - Nvtmd * lln((T2 - 1.0) / XJBVD);
                T1        =  XJBVD * lexp(-(BVD + VjdmRev) / Nvtmd);
                IVjdmRev  =  Isbd * (1.0 + T1);
                DslpRev   =  -Isbd * T1 / Nvtmd;
            end

            // Junction Capacitance
            Czbs     =  CJS_t * ASEJ;
            Czbssw   =  CJSWS_t * PSEJ;
            Czbsswg  =  CJSWGS_t * Weff0 * NFINtotal;
            Czbd     =  CJD_t * ADEJ;
            Czbdsw   =  CJSWD_t * PDEJ;
            Czbdswg  =  CJSWGD_t * Weff0 * NFINtotal;
        end

        // Generation-Recombination Current
        T0        =  Eg / Vtm * (TRatio - 1.0);
        T1        =  T0 / NTGEN_i;
        igentemp  =  lexp(T1);

    end  //  End of temperature dependent calculations

    // ************************************************
    // *     Bias dependent calculations follow       *
    // ************************************************

    // Load Terminal Voltages
    vgs_noswap  =  devsign * V(`IntrinsicGate, si);
    vds_noswap  =  devsign * V(di, si);
    vgd_noswap  =  devsign * V(`IntrinsicGate, di);
    ves_jct     =  devsign * V(e, si);
    ved_jct     =  devsign * V(e, di);
    vge         =  devsign * V(`IntrinsicGate, e);

    // Source-Drain Interchange
    sigvds  =  1.0;
    if (vds_noswap < 0.0) begin
        sigvds  =  -1.0;
        vgs     =  vgs_noswap - vds_noswap;
        vds     =  -1.0 * vds_noswap;
        ves     =  ved_jct;
    end else begin
        vgs     =  vgs_noswap;
        vds     =  vds_noswap;
        ves     =  ves_jct;
    end
    vgsfb  =  vgs - deltaPhi;

    // Initialize Certain Variables to Zero to Prevent Unnecessary Update
    etaiv  =  0.0;
    Qes    =  0.0;
    Qesj   =  0.0;
    Qeg    =  0.0;
    Qed    =  0.0;
    Qedj   =  0.0;

    // Vds Smoothing
    vdsx  =  sqrt (vds * vds + 0.01) - 0.1;

    // Ves Smoothing
    if (BULKMOD != 0) begin
        vesx    =  ves - 0.5 * (vds - vdsx);
        vesmax  =  0.95 * PHIBE_i;
        T2      =  vesmax - vesx - 1.0e-3;
        veseff  =  vesmax - 0.5 * (T2 + sqrt(T2 * T2 + 0.004 * vesmax));
    end

    // Asymmetry Model
    T0  =  tanh(0.6 * vds_noswap / Vtm);
    wf  =  0.5 + 0.5 * T0;
    wr  =  1.0 - wf;
    if (ASYMMOD != 0) begin
        CDSCD_a     =  CDSCDR_i * wr + CDSCD_i * wf;
        ETA0_a      =  ETA0R_t * wr + ETA0_t * wf;
        PDIBL1_a    =  PDIBL1R_i * wr + PDIBL1_i * wf;
        PDIBL2_a    =  PDIBL2R_i * wr + PDIBL2_i * wf;
        MEXP_a      =  MEXPR_t * wr + MEXP_t * wf;
        PTWG_a      =  PTWGR_t * wr + PTWG_t * wf;
        VSAT1_a     =  VSAT1R_t * wr + VSAT1_t * wf;
        RSDR_a      =  RSDRR_t * wr + RSDR_t * wf;
        RDDR_a      =  RDDRR_t * wr + RDDR_t * wf;
        PCLM_a      =  PCLMR_i * wr + PCLM_i * wf;
        VSAT_a      =  VSATR_t * wr + VSAT_t * wf;
        KSATIV_a    =  KSATIVR_i * wr + KSATIV_i * wf;
        DVTSHIFT_a  =  DVTSHIFTR_i * wr + DVTSHIFT_i * wf;
        CIT_a       =  CITR_i * wr + CIT_i * wf;
        u0_a        =  u0r * wr + u0 * wf;
        UA_a        =  UAR_t*wr + UA_t * wf;
        UD_a        =  UDR_t * wr + UD_t * wf;
        UC_a        =  UCR_t * wr + UC_t * wf;
        EU_a        =  EUR_i * wr + EU_i * wf;
    end else begin
        CDSCD_a     =  CDSCD_i;
        ETA0_a      =  ETA0_t;
        PDIBL1_a    =  PDIBL1_i;
        PDIBL2_a    =  PDIBL2_i;
        MEXP_a      =  MEXP_t;
        PTWG_a      =  PTWG_t;
        VSAT1_a     =  VSAT1_t;
        RSDR_a      =  RSDR_t;
        RDDR_a      =  RDDR_t;
        PCLM_a      =  PCLM_i;
        VSAT_a      =  VSAT_t;
        KSATIV_a    =  KSATIV_i;
        DVTSHIFT_a  =  DVTSHIFT_i;
        CIT_a       =  CIT_i;
        u0_a        =  u0;
        UA_a        =  UA_t;
        UD_a        =  UD_t;
        UC_a        =  UC_t;
        EU_a        =  EU_i;
    end

    // Drain Saturation Voltage
    inv_MEXP  =  1.0 / MEXP_a;

    // SCE, DIBL, SS Degradation Effects (Ref: BSIM4 Model)
    phist  =  0.4 + phib + PHIN_i;
    T1     =  2.0 * (Cins / Weff_UFCM) / (rc + 2.0);
    cdsc   =  Theta_SW * (CDSC_i + CDSCD_a * vdsx);

    if (!$param_given(NVTM))
        nVtm  =  Vtm * ThetaSS * (1.0 + (CIT_a + cdsc) / T1);
    else nVtm  =  NVTM;

    // temp deped UFCM
    qdep                  =  Qdep_ov_Cins / nVtm;
    vth_fixed_factor_SI   =  ln(Cins * nVtm/(`q * Nc * 2.0 * Ach));
    vth_fixed_factor_Sub  =  ln((qdep * rc) * (qdep * rc) / ((exp(qdep * rc) - qdep * rc - 1.0))) + vth_fixed_factor_SI;
    q0                    =  10.0 * nVtm / rc + 2.0 * qbs;

    // New QM parameter calculation: fieldnormalizationfactor, auxQMfact, QMFACTORCVfinal
    fieldnormalizationfactor  =  Vtm * Cins / (Weff_UFCM * epssub);
    auxQMfact                 =  pow(((3.0 / 4.0) * 3.0 * `HBAR * 2.0 * `M_PI * `q / (4.0 * sqrt(2.0 * mx))), 2.0 / 3.0);
    QMFACTORCVfinal           =  QMFACTORCV * auxQMfact * pow(fieldnormalizationfactor, 2.0 / 3.0) * (1/(`q * Vtm));

    dvth_vtroll  =  -DVT0_i * Theta_SCE * (vbi - phist);
    dvth_dibl    =  -ETA0_a * Theta_DIBL * vdsx + (DVTP0_i * Theta_DITS * pow(vdsx, DVTP1_i));
    dvth_rsce    =  K1RSCE_i * Theta_RSCE * sqrt(phist);
    dvth_all     =  dvth_vtroll + dvth_dibl + dvth_rsce + dvth_temp + DVTSHIFT_a;
    vgsfb        =  vgsfb - dvth_all;

    // Vgs Clamping for Inversion Region Calculation in Accumulation
    beta0     =  u0_a * cox * Weff0 / Leff;
    T0        =  -(dvch_qm + nVtm * lln(2.0 * cox * Imin / (beta0 * nVtm * `q * Nc * TFIN)));
    T1        =  vgsfb + T0 + DELVTRAND;
    vgsfbeff  =  hypsmooth(T1 , 1.0e-4) - T0;

    // Core Model Calculation at Source Side
    vch  =  0.0 + dvch_qm;

    if (BULKMOD != 0) begin
        T1  =  hypsmooth(2.0 * phib + vch - ves, 0.1);
        T3  =  (-K1_t / (2.0 * nVtm)) * (sqrt(T1) - sqrt(2.0 * phib));
        T0  =  -qdep - T3 + vth_fixed_factor_Sub + QMFACTORCVfinal * pow(-qdep, 2.0/3.0);
        T1  =  -qdep - T3 + vth_fixed_factor_SI;
    end else begin
        T0  =  -qdep + vth_fixed_factor_Sub + QMFACTORCVfinal * pow(-qdep, 2.0/3.0);
        T1  =  -qdep + vth_fixed_factor_SI;
    end
    T2  =  (vgsfbeff - vch) / nVtm;
    F0  =  -T2 + T1;
    T3  =  0.5 * (T2 - T0);
    qm  =  exp(T3);
    if (qm > 1.0e-7) begin
        T7  =  ln(1.0 + qm);
        qm  =  2.0 * (1.0 - sqrt(1.0 + T7 * T7));
        T8  =  (qm * ALPHA_UFCM + qdep) * rc;
        T4  =  T8 / (exp(T8) - T8 - 1.0);
        T5  =  T8 * T4;
        e0  =  F0 - qm + ln(-qm) + ln(T5) + QMFACTORCVfinal * pow(-(qm + qdep), 2.0 / 3.0);
        e1  =  -1.0 + 1.0 / qm + (2.0 / T8 - T4 - 1.0) * rc - (2.0 / 3.0) * QMFACTORCVfinal * pow(-(qm + qdep), -1.0 / 3.0);
        e2  =  -1.0 / (qm * qm) - (2.0 / 9.0) * QMFACTORCVfinal * pow(-(qm + qdep), -4.0/3.0);
        qm  =  qm - (e0 / e1) * (1.0 + (e0 * e2) / (2.0 * e1 * e1));
        T8  =  (qm * ALPHA_UFCM + qdep) * rc;
        T4  =  T8 / (exp(T8) - T8 - 1.0);
        T5  =  T8 * T4;
        e0  =  F0 - qm + ln(-qm) + ln(T5) + QMFACTORCVfinal * pow(-(qm + qdep), 2.0 / 3.0);
        e1  =  -1.0 + 1.0 / qm + (2.0 / T8 - T4 - 1.0) * rc - (2.0 / 3.0) * QMFACTORCVfinal * pow(-(qm + qdep), -1.0/3.0);
        e2  =  -1.0 / (qm * qm) - (2.0 / 9.0) * QMFACTORCVfinal * pow(-(qm + qdep), -4.0/3.0);
        qm  =  qm - (e0 / e1) * (1.0 + (e0 * e2) / (2.0 * e1 * e1));
    end else begin
        qm  =  -qm * qm;
    end
    qis  =  -qm * nVtm;

    // Drain Saturation Voltage
    Eeffs  =  EeffFactor * (qbs + eta_mu * qis);
    qb0    =  1.0e-2 / cox;
    T2     =  pow(0.5 * (1.0 + abs((qis) / qb0)), UCS_t);
    if (BULKMOD != 0) begin
        T3  =  (UA_a + UC_a * veseff) * pow(abs(Eeffs), EU_a) + UD_a / T2;
    end else begin
        T3  =  UA_a * pow(abs(Eeffs), EU_a) + UD_a / T2;
    end
    Dmobs  =  1.0 + T3;
    Dmobs  =  Dmobs / U0MULT;

    if (RDSMOD == 1) begin
        Rdss  =  0.0;
    end else if (RDSMOD == 0) begin
        T4    =  1.0 + PRWGS_i * qis;
        T1    =  1.0 / T4;
        T0    =  0.5 * (T1 + sqrt(T1 * T1 + 0.01));
        Rdss  =  (RDSWMIN_i + RDSW_i * T0) * WeffWRFactor * NFINtotal * rdstemp;
    end else begin
        T4    =  1.0 + PRWGS_i * qis;
        T1    =  1.0 / T4;
        T0    =  0.5 * (T1 + sqrt(T1 * T1 + 0.01));
        Rdss  =  (RSourceGeo + RDrainGeo + RDSWMIN_i + RDSW_i * T0) * WeffWRFactor * NFINtotal * rdstemp;
    end

    Esat   =  2.0 * VSAT_a / u0_a * Dmobs;
    EsatL  =  Esat * Leff;
    T6     =  KSATIV_a * (qis +  2 * Vtm);

    if (Rdss == 0.0) begin
        Vdsat  =  EsatL * T6 / (EsatL + T6);
    end else begin
        WVCox  =  Weff0 * VSAT_a * cox;
        T0     =  WVCox * Rdss;
        Ta     =  2.0 * T0;
        Tb     =  T6 + EsatL + 3.0 * T6 * T0;
        Tc     =  T6 * (EsatL + 2.0 * T6 * T0);
        Vdsat  =  (Tb - sqrt(Tb * Tb - 2.0 * Ta * Tc)) / Ta;
    end
    Vdsat   =  hypsmooth(Vdsat - 1.0e-3, 1.0e-5) + 1.0e-3;
    T7      =  pow(vds / Vdsat , MEXP_a);
    T8      =  pow(1.0 + T7, inv_MEXP);
    Vdseff  =  vds / T8;

    if (Vdseff > vds) begin
        Vdseff  =  vds;
    end

    // Core Model Calculation at Drain Side
    vch  =  Vdseff + dvch_qm;

    if (BULKMOD != 0) begin
        T1  = hypsmooth(2.0 * phib + vch - ves, 0.1);
        T3  = (-K1_t / (2.0 * nVtm)) * (sqrt(T1) - sqrt(2.0 * phib));
        T0  = -qdep - T3 + vth_fixed_factor_Sub + QMFACTORCVfinal * pow(-qdep, 2.0 / 3.0);
        T1  = -qdep - T3 + vth_fixed_factor_SI;
    end else begin
        T0  = -qdep + vth_fixed_factor_Sub + QMFACTORCVfinal * pow(-qdep, 2.0 / 3.0);
        T1  = -qdep + vth_fixed_factor_SI;
    end
    T2  =  (vgsfbeff - vch) / nVtm;
    F0  =  -T2 + T1;
    T3  =  (T2 - T0) * 0.5;
    qm  =  exp(T3);
    if (qm > 1.0e-7) begin
        T7  =  ln(1.0 + qm);
        qm  =  2.0 * (1.0 - sqrt(1.0 + T7 * T7));
        T8  =  (qm * ALPHA_UFCM + qdep) * rc;
        T4  =  T8 / (exp(T8) - T8 - 1.0);
        T5  =  T8 * T4;
        e0  =  F0 - qm + ln(-qm) + ln(T5) + QMFACTORCVfinal * pow(-(qm + qdep), 2.0 / 3.0);
        e1  =  -1.0 + (1.0 / qm) + (2.0 / T8 - T4 - 1.0) * rc - (2.0 / 3.0) * QMFACTORCVfinal * pow(-(qm + qdep), -1.0 / 3.0);
        e2  =  -1.0 / (qm * qm) - (2.0 / 9.0) * QMFACTORCVfinal * pow(-(qm + qdep), -4.0 / 3.0);
        qm  =  qm - (e0 / e1)*(1.0 + (e0 * e2) / (2.0 * e1 * e1));
        T8  =  (qm * ALPHA_UFCM + qdep) * rc;
        T4  =  T8 / (exp(T8) - T8 - 1.0);
        T5  =  T8 * T4;
        e0  =  F0 - qm + ln(-qm) + ln(T5) + QMFACTORCVfinal * pow(-(qm + qdep), 2.0/3.0);
        e1  =  -1.0 + (1.0 / qm) + (2.0 / T8 - T4 - 1.0) * rc - (2.0 / 3.0) * QMFACTORCVfinal * pow(-(qm + qdep), -1.0 / 3.0);
        e2  =  -1.0 / (qm * qm) - (2.0 / 9.0) * QMFACTORCVfinal * pow(-(qm + qdep), -4.0 / 3.0);
        qm  =  qm - (e0 / e1) * (1.0 + (e0 * e2) / (2.0 * e1 * e1));
    end else begin
        qm  =  -qm * qm;
    end
    qid  =  -qm * nVtm;

    if (BULKMOD != 0) begin
        T9  =  (K1_t / (2.0 * nVtm)) * sqrt(Vtm);
        T0  =  T9 / 2.0;
        T2  =  (vge - (deltaPhi - Eg - Vtm * ln(NBODY / Nc) + DELVFBACC)) / Vtm;
        if ((T2 * Vtm) > phib + T9 * sqrt(phib * Vtm)) begin
            T1   =  sqrt(T2 - 1.0 + T0 * T0) - T0;
            T10  =  1.0 + T1 * T1;
        end else begin
            T3   =  T2 * 0.5 - 3.0 * (1.0 + T9 / sqrt(2.0));
            T10  =  T3 + sqrt(T3 * T3 + 6.0 * T2);
            if (T2 < 0.0) begin
                T4   =  (T2 - T10) / T9;
                T10  = -ln(1.0 - T10 + T4 * T4 );
            end else begin
                T11  = exp(-T10);
                T4   = sqrt(T2 - 1.0 + T11 + T0 * T0) - T0;
                T10  = 1.0 - T11 + T4 * T4;
            end
        end
        T6  =  exp(-T10) - 1.0;
        T7  =  sqrt(T6 + T10);
        if (T10 > 1.0e-15) begin
            e0   =  -(T2 - T10) + T9 * T7;
            e1   =  1.0 - T9 * 0.5 * T6 / T7;
            T8   =  T10 - (e0 / e1);
            T11  =  exp(-T8) - 1.0;
            T12  =  sqrt(T11 + T8);
            qba  =  -T9 * T12 * Vtm;
        end else begin
            if (T10 < -1.0e-15) begin
                e0   =  -(T2 - T10) - T9 * T7;
                e1   =  1.0 + T9 * 0.5 * T6 / T7;
                T8   =  T10 - e0 / e1;
                T12  =  T9 * sqrt(exp(-T8) + T8 - 1.0);
            end else begin
                T12  =  0.0;
                T8   =  0.0;
            end
            qba  =  T12 * Vtm;
        end
        qi_acc_for_QM  =  T9 * exp(-T8 / 2.0) * Vtm;

        psipclamp  =  0.5 * (T8 + 1.0 + sqrt((T8 - 1.0) * (T8 - 1.0) + 0.25 * 2.0 * 2.0));
        sqrtpsip   =  sqrt(psipclamp);
        nq         =  1.0 + T9 / sqrtpsip;
    end

    // Drain Side and Average Potential / Charge
    qia  =  0.5 * (qis + qid);
    dqi  =  qis - qid;

    T0   = pow(Vdseff, 2.0) / 6.25e-4;   // pow(Vdseff,2.0) / pow(25e-3, 2.0)
    if (CHARGEWF != 0.0)
        qia2 = 0.5 * (qis + qid) + CHARGEWF * (1.0 - lexp(-T0)) * 0.5 * dqi;
    else
        qia2 = 0.5 * (qis + qid);

    `ifdef __DEBUG__
        if (qis < 0.0) $strobe("Warning: negative source-side inversion carrier density.  Vgs=%f  Vds=%f  Vbs=%f qis=%e", V(g, s), V(d, s), V(e, s), qis);
        if (qid < 0.0) $strobe("Warning: negative drain-side inversion carrier density.  Vgs=%f  Vds=%f  Vbs=%f qid=%e", V(g, s), V(d, s), V(e, s), qid);
    `endif

    // Toxeff model for quantum mechanical effects
    // Normal operation (Vgs > Vfb)
    if (QMTCENCV_i > 0.0) begin
        T4      =  qia / QM0;
        T5      =  1.0 + pow(T4, PQM);
        Tcen    =  Tcen0 / T5;
        coxeff  =  1.0 / (1.0 / (cox * EOT / TOXP) + Tcen * QMTCENCV_i / epssub);
    end else begin
        coxeff  =  cox;
    end

    // Quantum Mechanical Effect Correction for Accumulation Side Cap (Vgs < Vfb)
    if (BULKMOD != 0 && QMTCENCVA_i != 0.0) begin
        T6       =  1.0 + pow(qi_acc_for_QM / QM0ACC, PQMACC);
        Tcen     =  Tcen0 / T6;
        cox_acc  =  1.0 / (1.0 / cox_acc + Tcen * QMTCENCVA_i / epssub);
    end

    // Multiplication Factor for I-V
    beta  =  u0_a * cox * Weff0 / Leff;

    // Mobility Degradation
    Eeffm  =  EeffFactor * (qba + eta_mu * qia2);
    T2     =  pow(0.5 * (1.0 + abs((qia2) / qb0)), UCS_t);
    if (BULKMOD != 0) begin
        T3  =  (UA_a + UC_a * veseff) * pow(abs(Eeffm), EU_a) + UD_a / T2;
    end else begin
        T3  =  UA_a * pow(abs(Eeffm), EU_a) + UD_a / T2;
    end
    Dmob  =  1.0 + T3;
    Dmob  =  Dmob / U0MULT;
    ueff  =  u0_a / Dmob;

    // Mobility Degradation for C-V
    Eeffm_cv  =  EeffFactor * (qba + eta_mu_cv * qia2);
    T3        =  UA_a * pow(abs(Eeffm_cv), EU_a) + UD_a / T2;
    Dmob_cv   =  1.0 + T3;
    Dmob_cv   =  Dmob_cv / U0MULT;

    // Calculate current and capacitance enhancement factors due to CLM and DIBL
    tmp  =  DROUT_i * Leff / scl + 1.0e-6;

    if (tmp < 40.0) begin
        DIBLfactor  =  0.5 * PDIBL1_a / (cosh(tmp) - 1.0) + PDIBL2_a;
    end else begin
        DIBLfactor  =  PDIBL1_a * exp(-tmp) + PDIBL2_a;
    end

    if (PVAG_i > 0.0) begin
        PVAGfactor  =  1.0 + PVAG_i * qia / EsatL;
    end else begin
        PVAGfactor  =  1.0 / (1.0 - PVAG_i * qia / EsatL);
    end

    if (Vdseff > vds) begin
        Vdseff  =  vds;
    end
    diffVds   =  vds - Vdseff;
    Vgst2Vtm  =  qia + 2.0 * Vtm;
    if (DIBLfactor > 0) begin
        T1      =  Vgst2Vtm;
        T3      =  T1 / (Vdsat + T1);
        VaDIBL  =  T1 / DIBLfactor * T3 * PVAGfactor;
        Moc     =  1.0 + diffVds / VaDIBL;
    end else begin
        Moc     =  1.0;
    end

    if (PCLM_a > 0.0) begin
        if (PCLMG_i < 0.0) begin
            T1  =  1.0 / (1.0 / PCLM_a - PCLMG_i * qia);
        end else begin
            T1  =  PCLM_a + PCLMG_i * qia;
        end
        Mclm  =  1.0 + T1 * lln(1.0 + (vds - Vdseff) / T1 / (Vdsat + EsatL));
    end else begin
        Mclm  =  1.0;
    end

    Moc  =  Moc * Mclm;

    // Current Degradation Factor Due to Velocity Saturation
    Esat1   =  2.0 * VSAT1_a / ueff;
    Esat1L  =  Esat1 * Leff;
    T0      =  lexp(PSAT_i * lln(dqi / Esat1L));
    Ta      =  (1.0 + lexp(1.0 / PSAT_i * lln(DELTAVSAT_i)));
    Dvsat   =  (1.0 + lexp(1.0 / PSAT_i * lln(DELTAVSAT_i + T0))) / Ta;
    Dvsat   =  Dvsat + 0.5 * PTWG_a * qia * dqi * dqi;

    // Non-Saturation Effect
    T0     =  A1_t + A2_t / (qia + 2.0 * nVtm);
    T1     =  T0 * dqi * dqi;
    T2     =  T1 + 1.0 - 0.001;
    T3     =  -1.0 + 0.5 * (T2 + sqrt(T2 * T2 + 0.004));   // max(T1, -1.0)
    Nsat   =  0.5 * (1.0 + sqrt(1.0 + T3));
    Dvsat  =  Dvsat * Nsat;

    // Lateral Non-uniform doping effect (IV-CV Vth shift) factor
    if (K0_t != 0) begin
        T1    =  K0_t / (max(0, K0SI_t + K0SISAT_t * dqi * dqi) * qia + 2.0 * nVtm);
        Mnud  =  lexp(-T1);
    end else begin
        Mnud  =  1.0;
    end

    // Body-Effect Factor for BULKMOD = 2
    if (BULKMOD == 2) begin
        T0      =  hypsmooth((K2_t + K2SAT_t * vdsx), 1.0e-6);
        T1      =  T0 / (max(0, K2SI_t + K2SISAT_t * dqi * dqi) * qia + 2.0 * nVtm);
        T3      =  sqrt(PHIBE_i - veseff) - sqrt(PHIBE_i);
        Mob     =  lexp(- T1 * T3);
    end else
        Mob  =  1.0;

    // Velocity Saturation Factor for C-V
    EsatCV   =  2.0 * VSATCV_t * Dmob_cv / u0_a;
    EsatCVL  =  EsatCV * LeffCV;
    T0       =  lexp(PSATCV_i * lln(dqi / EsatCVL));
    Ta       =  (1.0 + lexp(1.0 / PSATCV_i * lln(DELTAVSATCV_i)));
    DvsatCV  =  (1.0 + lexp(1.0 / PSATCV_i * lln(DELTAVSATCV_i + T0))) / Ta;

    // Channel Length Modulation factor for C-V
    if (PCLMCV_i != 0) begin
        MclmCV  =  1.0 + PCLMCV_i * lln(1.0 + (vds - Vdseff) / PCLMCV_i / (Vdsat + EsatCVL));
    end else begin
        MclmCV  =  1.0;
    end

    // Calculating fixed body charge qb with sign (Here to avoid multiple calculation in NQSMOD=3 case)
    qb  =  -`q * nbody * Ach * LeffCV;

    // ************************************************
    // *     Current and Charges Calculations         *
    // ************************************************
    // Quasi Static I-V Model
    T1     =  qia;
    etaiv  =  q0 / (q0 + qia);
    T2     =  (2.0 - etaiv) * nVtm;
    ids0_ov_dqi  =  T1 + T2;
    ids0   =  ids0_ov_dqi * dqi;

    // S/D Series Resistance
    `include "bsimcmg_rdsmod.include"

    ids  =  NFINtotal * beta * ids0 * Moc * Mnud * Mob / (Dmob * Dvsat * Dr);
    ids  =  ids * IDS0MULT;

    // Quasi Static C-V Model
    `include "bsimcmg_quasi_static_cv.include"

    // Parasitic Capacitances
    // Bias-dependent overlap capacitances (CGEOMOD = 0 and 2)
    qgs_ov  =  0.0;
    qgd_ov  =  0.0;
    if (CGEOMOD != 1) begin
        T1           =  NFINtotal * WeffCV0 * devsign;
        T2           =  devsign * V(`GateEdgeNode, si);
        T0           =  T2 - vfbsdcv + `DELTA_1;
        vgs_overlap  =  0.5 * (T0 - sqrt(T0 * T0 + 4.0 * `DELTA_1));
        qgs_ov       =  T1 * (CGSL_i * (T2 - vfbsdcv - vgs_overlap - 0.5 * CKAPPAS_i * (sqrt(1.0 - 4.0 * vgs_overlap / CKAPPAS_i) - 1.0)) + CGSO_i * T2);
        T2           =  devsign * V(`GateEdgeNode, di);
        T0           =  T2 - vfbsdcv + `DELTA_1;
        vgd_overlap  =  0.5 * (T0 - sqrt(T0 * T0 + 4.0 * `DELTA_1));
        qgd_ov       =  T1 * (CGDL_i * (T2 - vfbsdcv - vgd_overlap - 0.5 * CKAPPAD_i * (sqrt(1.0 - 4.0 * vgd_overlap / CKAPPAD_i) - 1.0)) + CGDO_i * T2);
    end

    if (CGEOMOD == 0) begin
        T1      =  NFINtotal * WeffCV0;      // Fringe caps dont see QM effects
        qgs_fr  =  T1 * CFS_i * V(`GateEdgeNode, si);
        qgd_fr  =  T1 * CFD_i * V(`GateEdgeNode, di);
        qgs_parasitic  =  qgs_ov + qgs_fr;
        qgd_parasitic  =  qgd_ov + qgd_fr;
    end else if (CGEOMOD == 1) begin        // CGEO1SW=1 enables parameters to be in F per fin, per gate-finger, per unit channel width
        if (CGEO1SW == 1) begin
            T0      =  NFINtotal * WeffCV0;
            COVS_i  =  T0 * COVS_i;
            COVD_i  =  T0 * COVD_i;
            cgsp    =  T0 * CGSP;
            cgdp    =  T0 * CGDP;
        end else begin
            cgsp    =  CGSP;
            cgdp    =  CGDP;
        end
        qgs_ov         =  COVS_i * V(`GateEdgeNode, si);
        qgd_ov         =  COVD_i * V(`GateEdgeNode, di);
        qgs_parasitic  =  qgs_ov;
        qgd_parasitic  =  qgd_ov;
        qgs_fr         =  cgsp * V(`GateEdgeNode, s);
        qgd_fr         =  cgdp * V(`GateEdgeNode, d);
    end else begin
        qgs_fr         =  Cfr_geo * V(`GateEdgeNode, si);
        qgd_fr         =  Cfr_geo * V(`GateEdgeNode, di);
        qgs_parasitic  =  qgs_ov + qgs_fr;
        qgd_parasitic  =  qgd_ov + qgd_fr;
    end

    // Drain-to-Source Fringe Capacitance Available for all CGEOMOD
    qds_fr  =  CDSP * V(d, s);

    // Impact Ionization Current (Ref: IIMOD = 1 from BSIM4 Model, IIMOD = 2 from BSIMSOI Model)
    Iii  =  0.0;
    if (IIMOD == 1) begin
        T0 =  (ALPHA0_t + ALPHA1_t * Leff) / Leff;
        if ((T0 <= 0.0) || (BETA0_t <= 0.0))
            Iii  =  0.0;
        else begin
            T1   =  -BETA0_t / (diffVds + 1.0e-30);
            Iii  =  T0 * diffVds * ids * lexp(T1);
        end
    end else if (IIMOD == 2) begin    //End of IIMOD=1
        ALPHAII  =  (ALPHAII0_t + ALPHAII1_t * Leff) / Leff;
        if (ALPHAII <= 0.0) begin
            Iii  =  0.0;
        end else begin
            T0       =  ESATII_i * Leff;
            T1       =  SII0_t * T0 / (1.0 + T0);
            T0       =  1.0 / (1.0 + hypsmooth(SII1_i * vgsfbeff, IIMOD2CLAMP1)); // T0 = 1 / (1 + SII1_i * vgsfbeff)
            T3       =  T0 + SII2_i;
            T2       =  hypsmooth(vgsfbeff * T3, IIMOD2CLAMP2);                   // T2 = vgsfbeff * T3
            T3       =  1.0 / (1.0 + SIID_i * vds);
            VgsStep  =  T1 * T2 * T3;
            Vdsatii  =  VgsStep * (1.0 - LII_i / Leff);
            Vdiff    =  vds - Vdsatii;
            T0       =  BETAII2_i + BETAII1_i * Vdiff + BETAII0_i * Vdiff * Vdiff;
            T1       =  sqrt(T0 * T0 + 1.0e-10);
            Ratio    =  -hypmax( -ALPHAII * lexp(Vdiff / T1), -10.0, IIMOD2CLAMP3);
            Iii      =  Ratio * ids;
        end
    end  // End of IIMOD=2

    // Gate Current (Ref: BSIM4 Model)
    igbinv  =  0.0;
    igbacc  =  0.0;
    igcs    =  0.0;
    igcd    =  0.0;
    igs     =  0.0;
    igd     =  0.0;

    // Igb
    if (IGBMOD != 0) begin
        // Igbinv
        T1      =  (qia - EIGBINV_i) / NIGBINV_i / Vtm;
        Vaux_Igbinv  =  NIGBINV_i * Vtm * lln(1.0 + lexp(T1));
        T2      =  AIGBINV_t - BIGBINV_i * qia;
        T3      =  1.0 + CIGBINV_i * qia;
        T4      =  -9.82222e11 * TOXG * T2 * T3;
        T5      =  lexp(T4);
        T6      =  3.75956e-7;
        igbinv  =  Weff0 * Leff * T6 * Toxratio * vge * Vaux_Igbinv * T5;
        igbinv  =  igbinv * igtemp;

        // Igbacc
        vfbzb  =  deltaPhi - (Eg / 2.0) - phib;
        T0     =  vfbzb - vge;
        T1     =  T0 / NIGBACC_i / Vtm;
        Vaux_Igbacc  =  NIGBACC_i * Vtm * lln(1.0 + lexp(T1));
        if (BULKMOD != 0) begin
            Voxacc  =  qi_acc_for_QM;
        end else begin
            if (vfbzb <= 0)
                Voxacc  =  0.5 * (T0 - 0.02 + sqrt((T0 - 0.02) * (T0 - 0.02) - 0.08 * vfbzb));
            else
                Voxacc  =  0.5 * (T0 - 0.02 + sqrt((T0 - 0.02) * (T0 - 0.02) + 0.08 * vfbzb));
        end
        T2      =  AIGBACC_t - BIGBACC_i * Voxacc;
        T3      =  1.0 + CIGBACC_i * Voxacc;
        T4      =  -7.45669e11 * TOXG * T2 * T3;
        T5      =  lexp(T4);
        T6      =  4.97232e-7;
        igbacc  =  Weff0 * Leff * T6 * Toxratio * vge * Vaux_Igbacc * T5;
        igbacc  =  igbacc * igtemp;

    end

    if (IGCMOD != 0) begin
        // Igcinv
        T1    =  AIGC_t - BIGC_i * qia;
        T2    =  1.0 + CIGC_i * qia;
        T3    =  -Bechvb * TOXG * T1 * T2;
        T4    =  qia * lexp(T3);
        T5    =  (vge + 0.5 * vdsx + 0.5 * (ves_jct + ved_jct));
        igc0  =  Weff0 * Leff * Aechvb * Toxratio * T4 * T5 * igtemp;

        // Gate-Current Partitioning
        Vdseffx  =  sqrt(Vdseff * Vdseff + 0.01) - 0.1;
        T1       =  PIGCD_i * Vdseffx;
        T1_exp   =  lexp(-T1);
        T3       =  T1 + T1_exp - 1.0 + 1.0e-4;
        T4       =  1.0 - (T1 + 1.0) * T1_exp + 1.0e-4;
        T5       =  T1 * T1 + 2.0e-4;
        igcd     =  igc0 * T4 / T5;
        igcs     =  igc0 * T3 / T5;

        // Igs
        T0       =  vgs_noswap - vfbsd;
        vgs_eff  =  sqrt(T0 * T0 + 1.0e-4);
        if (IGCLAMP == 1) begin
            T1  =  hypsmooth((AIGS_t - BIGS_i * vgs_eff), 1.0e-6);
            if (CIGS_i < 0.01) begin
                CIGS_i = 0.01;
            end
        end else begin
            T1  =  AIGS_t - BIGS_i * vgs_eff;
        end
        T2  =  1.0 + CIGS_i * vgs_eff;
        T3  =  -Bechvb * TOXG * POXEDGE_i * T1 * T2;
        T4  =  lexp(T3);
        if (sigvds > 0.0) begin
            igs  =  igsd_mult * DLCIGS * vgs_noswap * vgs_eff * T4;
        end else begin
            igd  =  igsd_mult * DLCIGS * vgs_noswap * vgs_eff * T4;
        end

        // Igd
        T0       =  vgd_noswap - vfbsd;
        vgd_eff  =  sqrt(T0 * T0 + 1.0e-4);
        if (IGCLAMP == 1) begin
            T1  =  hypsmooth((AIGD_t - BIGD_i * vgd_eff), 1.0e-6);
            if (CIGD_i < 0.01) begin
                CIGD_i = 0.01;
            end
        end else begin
            T1  =  AIGD_t - BIGD_i * vgd_eff;
        end
        T2  =  1.0 + CIGD_i * vgd_eff;
        T3  =  -Bechvb * TOXG * POXEDGE_i * T1 * T2;
        T4  =  lexp(T3);

        if (sigvds > 0.0) begin
            igd  =  igsd_mult * DLCIGD * vgd_noswap * vgd_eff * T4;
        end else begin
            igs  =  igsd_mult * DLCIGD * vgd_noswap * vgd_eff * T4;
        end
    end

    // GIDL/GISL Current (Ref: BSIM4 Model)
    igisl  =  0.0;
    igidl  =  0.0;

    if (GIDLMOD != 0) begin
        T0  =  epsratio * EOT;
        // GIDL
        if ((AGIDL_i <= 0.0) || (BGIDL_t <= 0.0)) begin
            T6  =  0.0;
        end else begin
            T1  =  (-vgd_noswap - EGIDL_i + vfbsd) / T0;
            T1  =  hypsmooth(T1, 1.0e-2);
            T2  =  BGIDL_t / (T1 + 1.0e-3);
            T3  =  lexp(PGIDL_i * lln(T1));
            if (BULKMOD != 0) begin
                T4   =  -ved_jct*ved_jct*ved_jct;
                T4a  =  CGIDL_i + abs(T4) + 1.0e-5;
                T5   =  hypsmooth(T4/T4a, 1.0e-6) - 1.0e-6;
                T6   =  AGIDL_i * Weff0 * T3 * lexp(-T2) * T5;
            end else begin
                T6   =  AGIDL_i * Weff0 * T3 * lexp(-T2) * vds_noswap;
            end
        end

        if (sigvds > 0.0) begin
            igidl  =  T6;
        end else begin
            igisl  =  T6;
        end

        // GISL
        if ((AGISL_i <= 0.0) || (BGISL_t <= 0.0)) begin
            T6  =  0.0;
        end else begin
            T1  =  (-vgs_noswap - EGISL_i + vfbsd) / T0;
            T1  =  hypsmooth(T1, 1.0e-2);
            T2  =  BGISL_t / (T1 + 1.0e-3);
            T3  =  lexp(PGISL_i * lln(T1));
            if (BULKMOD != 0) begin
                T4   =  -ves_jct * ves_jct * ves_jct;
                T4a  =  CGISL_i + abs(T4) + 1.0e-5;
                T5   =  hypsmooth(T4/T4a, 1.0e-6) - 1.0e-6;
                T6   =  AGISL_i * Weff0 * T3 * lexp(-T2) * T5;
            end else
                T6   =  AGISL_i * Weff0 * T3 * lexp(-T2) * (-vds_noswap);
        end

        if (sigvds > 0.0) begin
            igisl  =  T6;
        end else begin
            igidl  =  T6;
        end

    end  // End of GIDLMOD

    // Junction Current
    if (BULKMOD != 0) begin
        // Source-Side Junction Current
        if (Isbs > 0.0) begin
            if (ves_jct < VjsmRev) begin
                T0   =  ves_jct / Nvtms;
                T1   =  lexp(T0) - 1.0;
                T2   =  IVjsmRev + SslpRev * (ves_jct - VjsmRev);
                Ies  =  T1 * T2;
            end else if (ves_jct <= VjsmFwd) begin
                T0   =  ves_jct / Nvtms;
                T1   =  (BVS + ves_jct) / Nvtms;
                T2   =  lexp(-T1);
                Ies  =  Isbs * (lexp(T0) + XExpBVS - 1.0 - XJBVS * T2);
            end else begin
                Ies  =  IVjsmFwd + SslpFwd * (ves_jct - VjsmFwd);
            end
        end else begin
            Ies  =  0.0;
        end

        // Source-Side Junction Tunneling Current
        if (JTSS_t > 0.0) begin
            if ((VTSS - ves_jct) < (VTSS * 1.0e-3)) begin
                T0   =  -ves_jct / Vtm0 / NJTS_t;
                T1   =  lexp(T0 * 1.0e3) - 1.0;
                Ies  =  Ies - ASEJ * JTSS_t * T1;
            end else begin
                T0   =  -ves_jct / Vtm0 / NJTS_t;
                T1   =  lexp(T0 * VTSS / (VTSS - ves_jct)) - 1.0;
                Ies  =  Ies - ASEJ * JTSS_t * T1;
            end
        end

        if (JTSSWS_t > 0.0) begin
            if ((VTSSWS - ves_jct) < (VTSSWS * 1.0e-3)) begin
                T0   =  -ves_jct / Vtm0 / NJTSSW_t;
                T1   =  lexp(T0 * 1.0e3) - 1.0;
                Ies  =  Ies - PSEJ * JTSSWS_t * T1;
            end else begin
                T0   =  -ves_jct / Vtm0 / NJTSSW_t;
                T1   =  lexp(T0 * VTSSWS / (VTSSWS - ves_jct)) - 1.0;
                Ies  =  Ies - PSEJ * JTSSWS_t * T1;
            end
        end

        if (JTSSWGS_t > 0.0) begin
            if ((VTSSWGS - ves_jct) < (VTSSWGS * 1.0e-3)) begin
                T0   =  -ves_jct / Vtm0 / NJTSSWG_t;
                T1   =  lexp(T0 * 1.0e3) - 1.0;
                Ies  =  Ies - Weff0 * NFINtotal * JTSSWGS_t * T1;
            end else begin
                T0   =  -ves_jct / Vtm0 / NJTSSWG_t;
                T1   =  lexp(T0 * VTSSWGS / (VTSSWGS - ves_jct)) - 1.0;
                Ies  =  Ies - Weff0 * NFINtotal * JTSSWGS_t * T1;
            end
        end

        // Drain-Side Junction Current
        if (Isbd > 0.0) begin
            if (ved_jct < VjdmRev) begin
                T0   =  ved_jct / Nvtmd;
                T1   =  lexp(T0) - 1.0;
                T2   =  IVjdmRev + DslpRev * (ved_jct - VjdmRev);
                Ied  =  T1 * T2;
            end else if (ved_jct <= VjdmFwd) begin
                T0   =  ved_jct / Nvtmd;
                T1   =  (BVD + ved_jct) / Nvtmd;
                T2   =  lexp(-T1);
                Ied  =  Isbd * (lexp(T0) + XExpBVD - 1.0 - XJBVD * T2);
            end else
                Ied  =  IVjdmFwd + DslpFwd * (ved_jct - VjdmFwd);
        end else
            Ied  =  0.0;

        // Drain-Side Junction Tunneling Current
        if (JTSD_t > 0.0) begin
            if ((VTSD - ved_jct) < (VTSD * 1.0e-3)) begin
                T0 =  -ved_jct / Vtm0 / NJTSD_t;
                T1 =  lexp(T0 * 1.0e3) - 1.0;
                Ied  =  Ied - ADEJ * JTSD_t * T1;
            end else begin
                T0   =  -ved_jct / Vtm0 / NJTSD_t;
                T1   =  lexp(T0 * VTSD/ (VTSD - ved_jct)) - 1.0;
                Ied  =  Ied - ADEJ * JTSD_t * T1;
            end
        end
        if (JTSSWD_t > 0.0) begin
            if ((VTSSWD - ved_jct) < (VTSSWD * 1.0e-3)) begin
                T0   =  -ved_jct / Vtm0 / NJTSSWD_t;
                T1   =  lexp(T0 * 1.0e3) - 1.0;
                Ied  =  Ied - PDEJ * JTSSWD_t * T1;
            end else begin
                T0   =  -ved_jct / Vtm0 / NJTSSWD_t;
                T1   =  lexp(T0 * VTSSWD / (VTSSWD - ved_jct)) - 1.0;
                Ied  =  Ied - PDEJ * JTSSWD_t * T1;
            end
        end
        if (JTSSWGD_t > 0.0) begin
            if ((VTSSWGD - ved_jct) < (VTSSWGD * 1.0e-3)) begin
                T0   =  -ved_jct / Vtm0 / NJTSSWGD_t;
                T1   =  lexp(T0 * 1.0e3) - 1.0;
                Ied  =  Ied - Weff0 * NFINtotal * JTSSWGD_t * T1;
            end else begin
                T0   =  -ved_jct / Vtm0 / NJTSSWGD_t;
                T1   =  lexp(T0 * VTSSWGD / (VTSSWGD - ved_jct)) - 1.0;
                Ied  =  Ied - Weff0 * NFINtotal * JTSSWGD_t * T1;
            end
        end

        // Junction Capacitance (No Swapping)
        // Source-Substrate Junction
        `BSIM6JunctnCap(ves_jct, Czbs, PBS_t, SJS, MJS, MJS2, Qesj1)
        `BSIM6JunctnCap(ves_jct, Czbssw, PBSWS_t, SJSWS, MJSWS, MJSWS2, Qesj2)
        `BSIM6JunctnCap(ves_jct, Czbsswg, PBSWGS_t, SJSWGS, MJSWGS, MJSWGS2, Qesj3)
        Qesj  =  Qesj1 + Qesj2 + Qesj3;

        // Drain-Substrate Junction
        `BSIM6JunctnCap(ved_jct, Czbd, PBD_t, SJD, MJD, MJD2, Qedj1)
        `BSIM6JunctnCap(ved_jct, Czbdsw, PBSWD_t, SJSWD, MJSWD, MJSWD2, Qedj2)
        `BSIM6JunctnCap(ved_jct, Czbdswg, PBSWGD_t, SJSWGD, MJSWGD, MJSWGD2, Qedj3)
        Qedj  =  Qedj1 + Qedj2 + Qedj3;

    end // BULKMOD=0

    Qes  =  Qesj + csbox * ves_jct;
    Qed  =  Qedj + cdbox * ved_jct;

    // Gate-to-Substrate Parasitic Capacitance
    // Bias Independent Component
    Qeg  =  cgbox * devsign * V(e, `GateEdgeNode);
    if (BULKMOD != 0) begin
        // Bias Dependent Component
        T2           =  devsign * V(`GateEdgeNode, e);
        T3           =  T2 - deltaPhi + Eg / 2.0 + phib - DELVFBACC;
        T0           =  T3 + `DELTA_1;
        vge_overlap  =  0.5 * (T0 + sqrt(T0 * T0 + 4.0 * `DELTA_1));
        Qeg          =  Qeg - NFINtotal * LeffCV * (CGBL_i * (T3 - vge_overlap + 0.5 * CKAPPAB_i * ( sqrt(1.0 + 4.0 * vge_overlap / CKAPPAB_i) - 1.0 )));
    end

    // Generation-Recombination Component
    T0      =  vds;
    T1      =  T0 * (AIGEN_i + BIGEN_i * T0 * T0);
    idsgen  =  HFIN * TFIN * (Leff - 2.0 * LINTIGEN_i) * igentemp * T1;

    // NQS Gate Resistance (Ref: BSIM4 Model)
    T0      =  ueff * coxeff * Weff0 / Leff;

    `ifdef __NQSMOD1__
        if (NQSMOD == 1 && XRCRG1_i != 0) begin
            IdovVds  =  beta * ids0_ov_dqi * Moc / (Dmob * Dvsat * Dr);
            gcrg     =  NFINtotal * XRCRG1_i * (IdovVds + XRCRG2_i * Vtm * T0);
        end
    `endif

    `ifdef __NQSMOD2__
        if (NQSMOD == 2) begin
            IdovVds  =  beta * ids0_ov_dqi * Moc / (Dmob * Dvsat * Dr);
            gcrg     =  NFINtotal * XRCRG1_i * (IdovVds + XRCRG2_i * Vtm * T0);
            gtau     =  gcrg / (cox * Weff0 * Leff);
        end
    `endif

    // *** Multiply all current and charge components by NFINtotal ***
    // Note: Do not multiply ids, qg, qs, qd, qb, Ies, Ied, Qbs, Qbd with NFINtotal
    //       since it is already considered.

    igidl   =  NFINtotal * igidl;
    igisl   =  NFINtotal * igisl;
    igcd    =  NFINtotal * igcd;
    igcs    =  NFINtotal * igcs;
    igs     =  NFINtotal * igs;
    igd     =  NFINtotal * igd;
    igbinv  =  NFINtotal * igbinv;
    igbacc  =  NFINtotal * igbacc;
    idsgen  =  NFINtotal * idsgen;

    // Gate to Body Tunneling Current Empirical Partition for BULKMOD = 0
    igbs  =  0.0;
    igbd  =  0.0;
    if (BULKMOD == 0) begin
        igbs  =  (igbinv + igbacc) * wf;
        igbd  =  (igbinv + igbacc) * wr;
    end

    // Noise Models
    Esatnoi  =  2.0 * VSAT_a / ueff;  // Thermal noise and flicker noise

    // Flicker Noise (Ref: BSIM4 Model from K. K. Hung et al. TED 1990)
    if (NOIA > 0.0 || NOIB > 0.0 || NOIC > 0.0) begin
        Leffnoi    =  Leff - 2.0 * LINTNOI_i;
        Leffnoisq  =  Leffnoi * Leffnoi;
        if (EM <= 0.0) begin
            DelClm  =  0.0;
        end else begin
            T0      =  (diffVds / litl + EM) / Esatnoi;
            DelClm  =  litl * lln(T0);
            if (DelClm < 0.0) begin
                DelClm  =  0.0;
            end
        end
        T1     =  `q * `q * `q * Vtm * abs(ids) * ueff;
        T2     =  1.0e10 * coxeff * Leffnoisq;
        N0     =  coxeff * qis / `q;
        Nl     =  coxeff * qid / `q;
        Nstar  =  Vtm / `q * (coxeff + CIT_a);
        T3     =  NOIA * lln((N0 + Nstar) / (Nl + Nstar));
        T4     =  NOIB * (N0 - Nl);
        T5     =  0.5 * NOIC * (N0 * N0 - Nl * Nl);
        T6     =  `q * Vtm * ids * ids;
        T7     =  1.0e10 * Leffnoisq * Weff0 * NFINtotal;
        T8     =  NOIA + NOIB * Nl + NOIC * Nl * Nl;
        T9     =  (Nl + Nstar) * (Nl + Nstar);
        Ssi    =  T1 / T2 * (T3 + T4 + T5) + T6 / T7 * DelClm * T8 / T9;
        T10    =  NOIA * `q * Vtm;
        T11    =  Weff0 * NFINtotal * Leffnoi * 1.0e10 * Nstar * Nstar;
        Swi    =  T10 / T11 * ids * ids;
        T1     =  Swi + Ssi;
        if (T1 > 0.0) begin
            FNPowerAt1Hz  =  (Ssi * Swi) / T1;
        end else begin
            FNPowerAt1Hz  =  0.0;
        end
    end else begin
        FNPowerAt1Hz  =  0.0;
    end

    // Thermal Noise
    case (TNOIMOD)
        0 : begin   // Charge-based model (BSIM4 - TNOIMOD=0)
            T0     =  ueff * qinv;
            T1     =  T0 * Rdsi + Leff * Leff;
            Gtnoi  =  (T0 / T1) * NTNOI_i;
            sid    =  4.0 * Vtm * `q * Gtnoi;
        end
        1: begin   // Correlated Thermal Noise by Navid, November 2013, Reference BSIMSOI4.5.0
            `ifdef __TNOIMOD1__
                Abulk    =  1.0;
                Vgst2Vtm =  KSATIV_a * (qis + 2.0 * Vtm);
                etaa     =  1.0 - Vdseff * Abulk / Vgst2Vtm ;
                T0       =  1.0 - etaa;
                T1       =  1.0 + etaa;
                T2       =  T1 + 2.0 * Abulk * Vtm / (qia + 1.0e-10);
                T3       =  T2 * T2;
                T4       =  T0 * T0;
                T5       =  T3 * T3;
                T6       =  1.0 / (1.0 + Vdseff / EsatL) ;
                gamma    =  T6 * (0.5 * T1 + T0 * T0 / (6.0 * T2));
                delta    =  ((T1 / T3) - (5.0 * T1 + T2) * T4 / (15.0 * T5) + T4 * T4 / (9.0 * T5 * T2)) / (6.0 * T6 * T6 * T6);
                T7       =  T0 / T2;
                epsilon  =  (T7 + T7 * T7 * T7 / 3.0) / (6.0 * T6);
                T8       =  qia / EsatL;
                T8       =  T8 * T8;
                npart_c  =  RNOIC * (1.0 + T8 * TNOIC * Leff);
                T9       =  gamma * delta ;
                if (T9 > 0.0) begin
                    ctnoi  =  epsilon / sqrt( gamma * delta) * (2.5316 * npart_c);
                end else begin
                    ctnoi  =  1.0;
                end
                if (ctnoi > 1) begin
                    ctnoi = 1.0;
                end
                if (ctnoi < 0) begin
                    ctnoi = 0.0;
                end
                npart_beta  =  RNOIA * (1.0 + T8 * TNOIA * Leff);
                npart_theta =  RNOIB * (1.0 + T8 * TNOIB * Leff);
                gamma       =  gamma * (3.0 * npart_beta * npart_beta);
                delta       =  delta * (3.75 * npart_theta * npart_theta);
                T9          =  qia * 0.5 * T1;
                gche        =  beta * T9 * T6;
                noiGd0      =  NFINtotal * beta * qia / (1.0 + gche * Rdsi);
                GammaGd0    =  gamma * noiGd0;
                sid         =  4.0 * Vtm * `q * GammaGd0;
                C0          = NFINtotal * coxeff * WeffCV0 * LeffCV;
                if (gamma > 0.0 && delta > 0.0) begin
                    sf  =  (noiGd0 + 1.0e-15) / sqrt(delta / gamma);
                end else begin
                    sf  =  0.0;
                end
            `else
                $strobe("[Warning!] Although the model selector TNOIMOD is set to 1, the new correlated thermal noise model is not activated. Please uncomment \"`define __TNOIMOD1__\" in the bsimcmg.va.");
            `endif
        end
    endcase

    // Source and Drain Conductance for Thermal Noise Contribution
    if (RDSMOD != 2) begin
        gspr  =  1.0 / Rsource;    // Note: gspr considers all fins
        gdpr  =  1.0 / Rdrain;     // Note: gdpr considers all fins
    end

    // Loading Ids, Gate and Drain charges
    `ifdef __NQSMOD2__
        if (sigvds > 0.0) begin
            I(di, si) <+ devsign * ids;
        end else begin
            I(si, di) <+ devsign * ids;
        end

        if (NQSMOD == 2) begin
            I(`IntrinsicGate, si) <+ devsign * gtau * -V(q);
            I(di, si)             <+ devsign * xdpart * gtau * V(q);
        end else begin  // Quasi-static Stamping (Normal case)
            I(di, si)             <+ devsign * ddt(qd);
            I(`IntrinsicGate, si) <+ devsign * ddt(qg);
        end
    `else
        if (sigvds > 0.0) begin
            I(di, si) <+ devsign * ids;
        end else begin
            I(si, di) <+ devsign * ids;
        end
        I(di, si) <+ devsign * ddt(qd);
        I(`IntrinsicGate, si) <+ devsign * ddt(qg);
    `endif

    // Loading Other Currents
    if (sigvds > 0.0) begin
        I(di, si) <+ devsign * idsgen;
        I(`IntrinsicGate, si) <+ devsign * (igcs + igs);
        I(`IntrinsicGate, di) <+ devsign * (igcd + igd);
        if (BULKMOD != 0) begin
            I(di, e) <+ devsign * (igidl + Iii);
            I(si, e) <+ devsign * igisl;
            I(`IntrinsicGate, e) <+ devsign * (igbinv + igbacc);
        end else begin
            I(di, si) <+ devsign * (igidl + Iii);
            I(si, di) <+ devsign * igisl;
        end
    end else begin
        I(si, di) <+ devsign * idsgen;
        I(`IntrinsicGate, di) <+ devsign * (igcs + igs);
        I(`IntrinsicGate, si) <+ devsign * (igcd + igd);
        if (BULKMOD != 0) begin
            I(si, e) <+ devsign * (igidl + Iii);
            I(di, e) <+ devsign * igisl;
            I(`IntrinsicGate, e) <+ devsign * (igbinv + igbacc);
        end else begin
            I(si, di) <+ devsign * (igidl + Iii);
            I(di, si) <+ devsign * igisl;
        end
    end
    if (BULKMOD == 0) begin
        I(`IntrinsicGate, si) <+ devsign * igbs;
        I(`IntrinsicGate, di) <+ devsign * igbd;
    end

    if (BULKMOD != 0) begin
        I(e, si) <+ devsign * Ies;
        I(e, di) <+ devsign * Ied;
    end
    I(e, si) <+ devsign * ddt(Qes);
    I(e, di) <+ devsign * ddt(Qed);
    I(e, `GateEdgeNode) <+ devsign * ddt(Qeg);

    // Loading other charges
    I(`GateEdgeNode, si) <+ ddt(qgs_parasitic);
    I(`GateEdgeNode, di) <+ ddt(qgd_parasitic);
    I(d, s) <+ ddt(qds_fr);
    if (CGEOMOD == 1) begin
        I(`GateEdgeNode, s) <+ ddt(qgs_fr);
        I(`GateEdgeNode, d) <+ ddt(qgd_fr);
    end

    // Accumulation Charge for Bulk FET
    if (BULKMOD != 0) begin
        I(`IntrinsicGate, si) <+ devsign * ddt(qg_acc);
        I(e, si)              <+ devsign * ddt(qb_acc);
    end

    // External S/D Resistance
    if (RDSMOD == 2) begin
        V(d, di) <+ 0.0;
        V(s, si) <+ 0.0;
    end else begin
        I(d, di) <+ V(d, di) / Rdrain;
        I(s, si) <+ V(s, si) / Rsource;
    end

    // NQSMOD1 Gate Resistance Model
    `ifdef __NQSMOD1__
        if (NQSMOD == 1 && XRCRG1_i != 0)
            I(`GateEdgeNode, gi) <+ V(`GateEdgeNode, gi) * gcrg;
        else
            V(`GateEdgeNode, gi) <+ 0.0;
    `endif

    // NQSMOD2 BSIM4 Charge Deficit Model
    `ifdef __NQSMOD2__
        if (NQSMOD ==2) begin
            I(q) <+ ddt(qg - qb);
            I(q) <+ V(q) * gtau;
            I(q) <+ ddt(V(q));
        end else
            V(q) <+ 0.0;
    `endif

    // Gate Electrode Resistance
    `ifdef __RGATEMOD__
        if (RGATEMOD != 0)
            I(g, ge) <+ V(g, ge) * ggeltd;
        else
            V(g, ge) <+ 0.0;
    `endif

    // Flicker Noise
    I(di,si) <+ flicker_noise(FNPowerAt1Hz, EF, "flicker");

    // Thermal noise for parasitics
    if (RDSMOD != 2) begin
        I(d, di) <+ white_noise(4.0 * Vtm * `q * gdpr, "thermal");
        I(s, si) <+ white_noise(4.0 * Vtm * `q * gspr, "thermal");
    end

    `ifdef __RGATEMOD__
        if (RGATEMOD != 0)
            I(g, ge) <+  white_noise(4.0 * Vtm * `q * ggeltd, "thermal");
    `endif

    // Channel thermal noise and induced gate noise stamping
    // Implementation of correlated noise follows C. C. McAndrew, WCM 2005
    if (TNOIMOD == 0) begin
        I(di, si) <+ white_noise(sid, "thermal");
        `ifdef __TNOIMOD1__
            V(N) <+ 0.0;
        `endif
    end else begin
        `ifdef __TNOIMOD1__
            I(di,si) <+ white_noise(sid * abs(1.0 - ctnoi * ctnoi), "thermal");
            I(di,si) <+ ctnoi * V(N) * sf * SCALEN ;
            if (gamma > 0 && delta > 0) begin
                I(N) <+ V(N) * sf * SCALEN;
                I(N) <+ white_noise(sid/(sf*sf*SCALEN*SCALEN));
            end else begin
                I(N) <+ V(N) ;
            end
            I(`IntrinsicGate,si) <+ ddt(0.5 * C0 * SCALEN * V(N));
            I(`IntrinsicGate,di) <+ ddt(0.5 * C0 * SCALEN * V(N));
        `else
            $strobe("[Warning!] Although the model selector TNOIMOD is set to 1, the new correlated thermal noise model is not activated. Please uncomment \"`define __TNOIMOD1__\" in the bsimcmg.va.");
        `endif
    end

    // Gate Current Shot Noise
    if (IGCMOD != 0) begin
        if (sigvds > 0) begin
            I(`IntrinsicGate, si) <+ white_noise(2.0 * `q * abs(igcs + igs), "shot");
            I(`IntrinsicGate, di) <+ white_noise(2.0 * `q * abs(igcd + igd), "shot");
        end else begin
            I(`IntrinsicGate, di) <+ white_noise(2.0 * `q * abs(igcs + igs), "shot");
            I(`IntrinsicGate, si) <+ white_noise(2.0 * `q * abs(igcd + igd), "shot");
        end
    end

    if (IGBMOD != 0) begin
        if (BULKMOD != 0) begin
            I(`IntrinsicGate, e)  <+ white_noise(2.0 * `q * abs(igbinv + igbacc), "shot");
        end else begin
            I(`IntrinsicGate, si) <+ white_noise(2.0 * `q * abs(igbs), "shot");
            I(`IntrinsicGate, di) <+ white_noise(2.0 * `q * abs(igbd), "shot");
        end
    end

    // Self Heating
    `ifdef __SHMOD__
        if (SHMOD != 0 && RTH0 > 0) begin
            if (RDSMOD != 2) begin
                Pwr(ith_branch) <+ -(devsign * sigvds * V(di,si) * ids + V(d,di) * V(d,di) / Rdrain + V(s,si) * V(s,si) / Rsource);
            end else begin
                Pwr(ith_branch) <+ -(devsign * sigvds * V(di,si) * ids );
            end
        end
        Pwr(rth_branch) <+ Temp(rth_branch) * gth;
        Pwr(rth_branch) <+ ddt(Temp(rth_branch) * cth);
    `endif

    // Operating-Point information
    `ifdef __OPINFO__
        // W & L
        WEFF    =  Weff0;                            // Effective width for IV
        LEFF    =  Leff;                             // Effective length for IV
        WEFFCV  =  WeffCV0;                          // Effective width for CV
        LEFFCV  =  LeffCV;                           // Effective length for CV

        // Currents
        IDS     =  devsign * ids;                    // Intrinsic Drain Current (Electrical)
        if (sigvds > 0) begin                        // Total Source/Drain Currents (Physical)
            if (BULKMOD != 0) begin
                IDEFF  =  IDS + devsign * idsgen - devsign * (igd + igcd) + devsign * (Iii + igidl) - devsign * Ied;
                ISEFF  =  -IDS - devsign * idsgen - devsign * (igs + igcs) + devsign * (igisl) - devsign * Ies;
            end else begin
                IDEFF  =  IDS + devsign * idsgen - devsign * (igd + igcd + igbd) + devsign * (Iii + igidl - igisl);
                ISEFF  =  -IDS - devsign * idsgen - devsign * (igs + igcs + igbs) + devsign* (igisl - igidl);
            end
        end else begin
            if (BULKMOD != 0) begin
                IDEFF  =  -IDS - devsign * idsgen - devsign * (igs + igcs) + devsign * (igisl) - devsign * Ied;
                ISEFF  =  IDS + devsign * idsgen - devsign * (igd + igcd) + devsign * (Iii + igidl) -  devsign * Ies;
            end else begin
                IDEFF  =  -IDS - devsign * idsgen - devsign * (igs + igcs + igbd) + devsign * (igisl - igidl);
                ISEFF  =  IDS + devsign * idsgen - devsign * (igd + igcd + igbs) + devsign * (Iii + igidl - igisl);
            end
        end

        if (BULKMOD == 0) begin                      // Total Gate Current
            IGTOT  =  devsign * (igs + igd + igcs + igcd + igbs + igbd);
        end else begin
            IGTOT  =  devsign * (igs + igd + igcs + igcd + igbacc + igbinv);
        end

        IDSGEN =  sigvds * devsign * idsgen;         // Generation-Recombination Current (Physical)
        III    =  devsign * Iii;                     // Impact Ionization Current
        if (sigvds > 0) begin
            IGIDL  =  devsign * igidl;                // GIDL Current (Physical)
            IGISL  =  devsign * igisl;                // GISL Current (Physical)
        end else begin
            IGIDL  =  devsign * igisl;                // GIDL Current (Physical)
            IGISL  =  devsign * igidl;                // GISL Current (Physical)
        end

        if (BULKMOD != 0) begin
            IJSB  =  -devsign * Ies;                 // Source-Body Junction Current (Physical)
            IJDB  =  -devsign * Ied;                 // Drain-Body Junction Current (Physical)
        end else begin
            IJSB  =  0.0;
            IJDB  =  0.0;
        end

        if (BULKMOD != 0) begin
            ISUB  =  -III - IGIDL - IGISL - IJSB - IJDB - devsign * (igbinv + igbacc);      // Substrate Current
        end else begin
            ISUB  =  0.0;
        end

        // Misc Variables
        BETA    =  beta;                               // Drain Current prefactor per fin per finger
        VDSSAT  =  Vdsat;                              // Drain-Source saturation Voltage
        if (NGATE_i > 0)                               // Flatband Voltage
            VFB  =  -devsign * (phib + Vtm * lln(NGATE_i / ni));
        else
            VFB  =  PHIG_i - (EASUB + 0.5 * Eg + devsign * phib);

        // Threshold Voltage Calculation
        q0   =  10.0 * Vtm / rc + 2.0 * qbs;
        T1   =  Vtm * (Vtm + q0);
        T2   =  cox * cox * T1;
        T3   =  2.0 * `q * ni * epssub * Vtm;
        VTH  =  VFB + devsign * (Vtm * lln(T2 / T3) + dvch_qm + phib + qbs + Vtm + dvth_all - DELVTRAND);

        // Conductances
        GM   =  ddx(devsign * ids,V(`IntrinsicGate));      // Transconductance
        GDS  =  ddx(devsign * ids,V(di));                  // Output Conductance
        if (BULKMOD != 0)
            GMBS  =  ddx(devsign * ids,V(e));               // Body Transconductance
        else
            GMBS  =  0.0;

        // Intrinsic Charges (Physical) (Sriram: Not accurate for NQSMOD= 2 and 3)
        QGI  =  devsign * qg + devsign * qg_acc;
        QDI  =  devsign * qd;
        QSI  =  devsign * qs;
        QBI  =  devsign * (qb + qb_acc);

        // Total Charges (Sriram: Not accurate for NQSMOD= 2 and 3)
        QG  =  devsign * qg + qgs_parasitic + qgd_parasitic + (CGEOMOD == 1 ? qgs_fr + qgd_fr : 0) + devsign * qg_acc - devsign * Qeg;
        QD  =  devsign * qd - qgd_parasitic - (CGEOMOD == 1 ? qgd_fr : 0) - devsign * Qed;
        QS  =  devsign * qs - qgs_parasitic - (CGEOMOD == 1 ? qgs_fr : 0) - devsign * Qes;
        QB  =  devsign * (qb + qb_acc) + devsign * (Qeg + Qes + Qed);

        // Intrinsic Capacitances (Physical)
        CGGI  =  ddx(QGI, V(`IntrinsicGate));
        CGSI  =  ddx(-QGI, V(si));
        CGDI  =  ddx(-QGI, V(di));
        CGEI  =  ddx(-QGI, V(e));

        CSGI  =  ddx(-QSI, V(`IntrinsicGate));
        CSDI  =  ddx(-QSI, V(di));
        CSSI  =  ddx(QSI, V(si));
        CSEI  =  ddx(-QSI, V(e));                    // Should be zero everywhere

        CDGI  =  ddx(-QDI, V(`IntrinsicGate));
        CDDI  =  ddx(QDI, V(di));
        CDSI  =  ddx(-QDI, V(si));
        CDEI  =  ddx(-QDI, V(e));

        CEGI  =  ddx(-QBI, V(`IntrinsicGate));
        CEDI  =  ddx(-QBI, V(di));                   // Should be zero everywhere
        CESI  =  ddx(-QBI, V(si));                   // Should be zero everywhere
        CEEI  =  ddx(QBI, V(e));

        // Total Capacitances
        CGG  =  ddx(QG, V(`IntrinsicGate));
        CGS  =  ddx(-QG, V(si));
        CGD  =  ddx(-QG, V(di));
        CGE  =  ddx(-QG, V(e));

        CSG  =  ddx(-QS, V(`IntrinsicGate));
        CSD  =  ddx(-QS, V(di));
        CSS  =  ddx(QS, V(si));
        CSE  =  ddx(-QS, V(e));

        CDG  =  ddx(-QD, V(`IntrinsicGate));
        CDD  =  ddx(QD, V(di));
        CDS  =  ddx(-QD, V(si));
        CDE  =  ddx(-QD, V(e));

        CEG  =  ddx(-QB, V(`IntrinsicGate));
        CED  =  ddx(-QB, V(di));
        CES  =  ddx(-QB, V(si));
        CEE  =  ddx(QB, V(e));

        // Total extrinsic capacitance
        CGSEXT  =  ddx(-(qgs_parasitic + (CGEOMOD == 1 ? qgs_fr : 0)),V(si));  // Gate-Source Overlap + outer fringing
        CGDEXT  =  ddx(-(qgd_parasitic + (CGEOMOD == 1 ? qgd_fr : 0)),V(di));  // Gate-Drain Overlap + outer fringing
        CGBOV   =  ddx(Qeg,V(e));                                 // Gate-Body Overlap
        CGBOV   = -devsign * CGBOV;

        // Total of Junction Capacitance and Source/Drain-Body Overlap Capacitance
        CJST  =  ddx(Qes, V(si));
        CJST  =  -devsign * CJST;
        CJDT  =  ddx(Qed, V(di));
        CJDT  =  -devsign * CJDT;

        RSGEO  =  RSourceGeo;                                                // External bias independent Source Resistance
        RDGEO  =  RDrainGeo;                                                 // External bias independent Drain Resistance
        CFGEO  =  Cfr_geo;                                                   // Geometric Parasitic Cap for CGEOMOD=1

        // Output for Self-Heating Temperature
        T_TOTAL_K   =  DevTemp;
        T_TOTAL_C   =  DevTemp - `P_CELSIUS0;
        T_DELTA_SH  =  Temp(t);

        `ifdef __DEBUG__
            // Individual Gate Current Components
            IGS   =  devsign * igs;
            IGD   =  devsign * igd;
            IGCS  =  devsign * igcs;
            IGCD  =  devsign * igcd;
            if (BULKMOD == 0) begin
                IGBS  =  devsign * igbs;
                IGBD  =  devsign * igbd;
            end else begin
                IGBINV  =  devsign * igbinv;
                IGBACC  =  devsign * igbacc;
            end

            DIDSDVG  = ddx(ids, V(`IntrinsicGate));
            DIDSDVG  =  devsign * sigvds * DIDSDVG;
            DIDSDVS  = ddx(ids, V(si));
            DIDSDVS  =  devsign * sigvds * DIDSDVS;
            DIDSDVD  = ddx(ids, V(di));
            DIDSDVD  =  devsign * sigvds * DIDSDVD;
            `ifdef __SHMOD__
                DIDSDVTH  = ddx(ids, Temp(t));
                DIDSDVTH  =  devsign * sigvds * DIDSDVTH;
            `endif
            DIGSDVG  = ddx(igs + igcs, V(`IntrinsicGate));
            DIGSDVG  =  devsign * DIGSDVG;
            DIGSDVS  = ddx(igs + igcs, V(si));
            DIGSDVS  =  devsign * DIGSDVS;
            DIGSDVD  = ddx(igs + igcs, V(di));
            DIGSDVD  =  devsign * DIGSDVD;
            `ifdef __SHMOD__
                DIGSDVTH  = ddx(igs + igcs, Temp(t));
                DIGSDVTH  =  devsign * DIGSDVTH;
            `endif
            DIGDDVG  = ddx(igd + igcd, V(`IntrinsicGate));
            DIGDDVG  =  devsign * DIGDDVG;
            DIGDDVS  = ddx(igd + igcd, V(si));
            DIGDDVS  =  devsign * DIGDDVS;
            DIGDDVD  = ddx(igd + igcd, V(di));
            DIGDDVD  =  devsign * DIGDDVD;
            `ifdef __SHMOD__
                DIGDDVTH  = ddx(igd + igcd, Temp(t));
                DIGDDVTH  =  devsign * DIGDDVTH;
            `endif
            DIIIDVG  = ddx(Iii, V(`IntrinsicGate));
            DIIIDVG  =  devsign * DIIIDVG;
            DIIIDVS  = ddx(Iii, V(si));
            DIIIDVS  =  devsign * DIIIDVS;
            DIIIDVD  = ddx(Iii, V(di));
            DIIIDVD  =  devsign * DIIIDVD;
            `ifdef __SHMOD__
                DIIIDVTH  = ddx(Iii, Temp(t));
                DIIIDVTH  =  devsign * DIIIDVTH;
            `endif
            DIGIDLDVG  = ddx(igidl, V(`IntrinsicGate));
            DIGIDLDVG  =  devsign * DIGIDLDVG;
            DIGIDLDVS  = ddx(igidl, V(si));
            DIGIDLDVS  =  devsign * DIGIDLDVS;
            DIGIDLDVD  = ddx(igidl, V(di));
            DIGIDLDVD  =  devsign * DIGIDLDVD;
            `ifdef __SHMOD__
                DIGIDLDVTH  = ddx(igidl, Temp(t));
                DIGIDLDVTH  =  devsign * DIGIDLDVTH;
            `endif
            DIGISLDVG  = ddx(igisl, V(`IntrinsicGate));
            DIGISLDVG  =  devsign * DIGISLDVG;
            DIGISLDVS  = ddx(igisl, V(si));
            DIGISLDVS  =  devsign * DIGISLDVS;
            DIGISLDVD  = ddx(igisl, V(di));
            DIGISLDVD  =  devsign * DIGISLDVD;
            `ifdef __SHMOD__
                DIGISLDVTH  = ddx(igisl, Temp(t));
                DIGISLDVTH  =  devsign * DIGISLDVTH;
            `endif

            `ifdef __SHMOD__
                CGT  =  ddx(QG, Temp(t));
                CST  =  ddx(QS, Temp(t));
                CDT  =  ddx(QD, Temp(t));
            `endif
            ITH  =  ids * vds;
            `ifdef __SHMOD__
                DITHDVTH  =  ddx(ITH, Temp(t));
            `endif
            DITHDVG  =  ddx(ITH, V(`IntrinsicGate));
            DITHDVS  =  ddx(ITH, V(si));
            DITHDVD  =  ddx(ITH, V(di));
        `endif // __DEBUG__
    `endif // __OPINFO__
end //           analog block ends
//================================================
