//======================================================================================
//======================================================================================
// Filename: PSP102_module.include
//======================================================================================
//======================================================================================
//
//  (c) Copyright notice
//
//  Since 2012 until today, PSP has been co-developed by NXP Semiconductors and
//  Delft University of Technology. For this part of the model, each claim undivided
//  ownership and copyrights
//  Until and including 2011, PSP has been co-developed by NXP Semiconductors and
//  Arizona State University. For this part of the model, NXP Semiconductors claims
//  undivided ownership and copyrights.
//
//
//  Version: 102.5.0, December 2013
//
//======================================================================================
//======================================================================================
//
// Further information can be found in the file releasenotesPSP102.txt
//

//  Node definitions
inout      D, G, S, B;
electrical D;
electrical G;
electrical S;
electrical B;

// Extra internal nodes for correlated drain and gate noise
electrical NOI;
electrical NOI2;

// Extra branches for correlated drain and gate noise
branch (NOI) NOII;
branch (NOI) NOIR;
branch (NOI) NOIC;

// Internal nodes for gate and bulk resistors
electrical GP;
electrical BP;
electrical BI;
electrical BS;
electrical BD;

`ifdef NQSmodel

    // Internal nodes for spline collocation
    electrical INT1;
    electrical INT2;
    electrical INT3;
    electrical INT4;
    electrical INT5;
    electrical INT6;
    electrical INT7;
    electrical INT8;
    electrical INT9;

    branch(INT1) SPLINE1;
    branch(INT2) SPLINE2;
    branch(INT3) SPLINE3;
    branch(INT4) SPLINE4;
    branch(INT5) SPLINE5;
    branch(INT6) SPLINE6;
    branch(INT7) SPLINE7;
    branch(INT8) SPLINE8;
    branch(INT9) SPLINE9;

    branch(INT1) RES1;
    branch(INT2) RES2;
    branch(INT3) RES3;
    branch(INT4) RES4;
    branch(INT5) RES5;
    branch(INT6) RES6;
    branch(INT7) RES7;
    branch(INT8) RES8;
    branch(INT9) RES9;

`endif // NQSmodel

//////////////////////////
//
//  Model parameters
//
//////////////////////////

`ifdef LocalModel
    ///////////////////////////////////////////////////
    // PSP local model parameters
    ///////////////////////////////////////////////////

    //  Special model parameters, some are also simulator global variables
    parameter real    LEVEL    =  102                             `P(info="Model level" unit="");

    parameter integer TYPE     =  1        `from(   -1  ,1      ) `P(info="Channel type parameter, +1=NMOS -1=PMOS" unit="");
    parameter integer nmos     =  1        `from(    0,  1      ) `P(desc="MOS channel type");
    parameter integer pmos     =  1        `from(    0,  1      ) `P(desc="MOS channel type");
    parameter real    TR       =  21.0     `from( -273.0,inf    ) `P(info="nominal (reference) temperature" unit="C");

    //  Switch parameters that turn models or effects on or off
    parameter real    SWIGATE  =  0.0      `from(    0.0,1.0    ) `P(info="Flag for gate current, 0=turn off IG" unit="");
    parameter real    SWIMPACT =  0.0      `from(    0.0,1.0    ) `P(info="Flag for impact ionization current, 0=turn off II" unit="");
    parameter real    SWGIDL   =  0.0      `from(    0.0,1.0    ) `P(info="Flag for GIDL current, 0=turn off IGIDL" unit="");
    parameter real    SWJUNCAP =  0.0      `from(    0.0,3.0    ) `P(info="Flag for juncap, 0=turn off juncap" unit="");
    parameter real    SWJUNASYM=  0.0                             `P(info="Flag for asymmetric junctions; 0=symmetric, 1=asymmetric" unit="");
    parameter real    QMC      =  1.0      `from(    0.0,inf    ) `P(info="Quantum-mechanical correction factor" unit="");

    //  Process parameters
    parameter real    VFB      = -1.0                             `P(info="Flatband voltage at TR" unit="V");
    parameter real    STVFB    =  5.0e-4                          `P(info="Temperature dependence of VFB" unit="V/K");
    parameter real    TOX      =  2.0e-09  `from(  1e-10,inf    ) `P(info="Gate oxide thickness" unit="m");
    parameter real    EPSROX   =  3.9      `from(    1.0,inf    ) `P(info="Relative permittivity of gate dielectric" unit="");
    parameter real    NEFF     =  5.0e+23  `from(   1e20,1e26   ) `P(info="Effective substrate doping" unit="m^-3");
    parameter real    VNSUB    =  0.0                             `P(info="Effective doping bias-dependence parameter" unit="V");
    parameter real    NSLP     =  0.05     `from(   1e-3,inf    ) `P(info="Effective doping bias-dependence parameter" unit="V");
    parameter real    DNSUB    =  0.0      `from(    0.0,1.0    ) `P(info="Effective doping bias-dependence parameter" unit="V^-1");
    parameter real    DPHIB    =  0.0                             `P(info="Offset parameter for PHIB" unit="V");
    parameter real    NP       =  1.0e+26  `from(    0.0,inf    ) `P(info="Gate poly-silicon doping" unit="m^-3");
    parameter real    CT       =  0.0      `from(    0.0,inf    ) `P(info="Interface states factor" unit="");
    parameter real    TOXOV    =  2.0e-09  `from(  1e-10,inf    ) `P(info="Overlap oxide thickness" unit="m");
    parameter real    TOXOVD   =  2.0e-09  `from(  1e-10,inf    ) `P(info="Overlap oxide thickness for drain side" unit="m");
    parameter real    NOV      =  5.0e+25  `from(   1e20,1e27   ) `P(info="Effective doping of overlap region" unit="m^-3");
    parameter real    NOVD     =  5.0e+25  `from(   1e20,1e27   ) `P(info="Effective doping of overlap region for drain side" unit="m^-3");

    //  DIBL parameters
    parameter real    CF       =  0.0      `from(    0.0,inf    ) `P(info="DIBL-parameter" unit="");
    parameter real    CFB      =  0.0      `from(    0.0,1.0    ) `P(info="Back bias dependence of CF" unit="V^-1");

    //  Mobility parameters
    parameter real    BETN     =  7e-2     `from(    0.0,inf    ) `P(info="Channel aspect ratio times zero-field mobility" unit="m^2/V/s");
    parameter real    STBET    =  1.0                             `P(info="Temperature dependence of BETN" unit="");
    parameter real    MUE      =  0.5      `from(    0.0,inf    ) `P(info="Mobility reduction coefficient at TR" unit="m/V");
    parameter real    STMUE    =  0.0                             `P(info="Temperature dependence of MUE" unit="");
    parameter real    THEMU    =  1.5      `from(    0.0,inf    ) `P(info="Mobility reduction exponent at TR" unit="");
    parameter real    STTHEMU  =  1.5                             `P(info="Temperature dependence of THEMU" unit="");
    parameter real    CS       =  0.0      `from(    0.0,inf    ) `P(info="Coulomb scattering parameter at TR" unit="");
    parameter real    STCS     =  0.0                             `P(info="Temperature dependence of CS" unit="");
    parameter real    XCOR     =  0.0      `from(    0.0,inf    ) `P(info="Non-universality factor" unit="V^-1");
    parameter real    STXCOR   =  0.0                             `P(info="Temperature dependence of XCOR" unit="");
    parameter real    FETA     =  1.0      `from(    0.0,inf    ) `P(info="Effective field parameter" unit="");

    //  Series-resistance parameters (for resistance modeling as part of intrinsic mobility reduction)
    parameter real    RS       =  30       `from(    0.0,inf    ) `P(info="Series resistance at TR" unit="Ohm");
    parameter real    STRS     =  1.0                             `P(info="Temperature dependence of RS" unit="");
    parameter real    RSB      =  0.0      `from(   -0.5,1.0    ) `P(info="Back-bias dependence of series resistance" unit="V^-1");
    parameter real    RSG      =  0.0      `from(   -0.5,inf    ) `P(info="Gate-bias dependence of series resistance" unit="V^-1");

    //  Velocity saturation parameters
    parameter real    THESAT   =  1.0      `from(    0.0,inf    ) `P(info="Velocity saturation parameter at TR" unit="V^-1");
    parameter real    STTHESAT =  1.0                             `P(info="Temperature dependence of THESAT" unit="");
    parameter real    THESATB  =  0.0      `from(   -0.5,1.0    ) `P(info="Back-bias dependence of velocity saturation" unit="V^-1");
    parameter real    THESATG  =  0.0      `from(   -0.5,inf    ) `P(info="Gate-bias dependence of velocity saturation" unit="V^-1");

    //  Saturation voltage parameters
    parameter real    AX       =  3.0      `from(    2.0,inf    ) `P(info="Linear/saturation transition factor" unit="");

    //  Channel length modulation (CLM) parameters
    parameter real    ALP      =  0.01     `from(    0.0,inf    ) `P(info="CLM pre-factor" unit="");
    parameter real    ALP1     =  0.00     `from(    0.0,inf    ) `P(info="CLM enhancement factor above threshold" unit="V");
    parameter real    ALP2     =  0.00     `from(    0.0,inf    ) `P(info="CLM enhancement factor below threshold" unit="V^-1");
    parameter real    VP       =  0.05     `from(  1e-10,inf    ) `P(info="CLM logarithm dependence factor" unit="V");

    //  Impact ionization (II) parameters
    parameter real    A1       =  1.0      `from(   0.0,inf     ) `P(info="Impact-ionization pre-factor" unit="");
    parameter real    A2       = 10.0      `from(   0.0,inf     ) `P(info="Impact-ionization exponent at TR" unit="V");
    parameter real    STA2     =  0.0                             `P(info="Temperature dependence of A2" unit="V");
    parameter real    A3       =  1.0      `from(   0.0,inf     ) `P(info="Saturation-voltage dependence of impact-ionization" unit="");
    parameter real    A4       =  0.0      `from(   0.0,inf     ) `P(info="Back-bias dependence of impact-ionization" unit="V^-0.5");

    //  Gate current parameters
    parameter real    GCO      =  0.0      `from( -10.0,10.0    ) `P(info="Gate tunnelling energy adjustment" unit="");
    parameter real    IGINV    =  0.0      `from(   0.0,inf     ) `P(info="Gate channel current pre-factor" unit="A");
    parameter real    IGOV     =  0.0      `from(   0.0,inf     ) `P(info="Gate overlap current pre-factor" unit="A");
    parameter real    IGOVD    =  0.0      `from(   0.0,inf     ) `P(info="Gate overlap current pre-factor for drain side" unit="A");
    parameter real    STIG     =  2.0                             `P(info="Temperature dependence of IGINV and IGOV" unit="");
    parameter real    GC2      =  0.375    `from(   0.0,10.0    ) `P(info="Gate current slope factor" unit="");
    parameter real    GC3      =  0.063    `from(  -2.0,2.0     ) `P(info="Gate current curvature factor" unit="");
    parameter real    CHIB     =  3.1      `from(   1.0,inf     ) `P(info="Tunnelling barrier height" unit="V");

    //  Gate Induced Drain/Source Leakage (GIDL) parameters
    parameter real    AGIDL    =  0.0      `from(   0.0,inf     ) `P(info="GIDL pre-factor" unit="A/V^3");
    parameter real    AGIDLD   =  0.0      `from(   0.0,inf     ) `P(info="GIDL pre-factor for drain side" unit="A/V^3");
    parameter real    BGIDL    = 41.0      `from(   0.0,inf     ) `P(info="GIDL probability factor at TR" unit="V");
    parameter real    BGIDLD   = 41.0      `from(   0.0,inf     ) `P(info="GIDL probability factor at TR for drain side" unit="V");
    parameter real    STBGIDL  =  0.0                             `P(info="Temperature dependence of BGIDL" unit="V/K");
    parameter real    STBGIDLD =  0.0                             `P(info="Temperature dependence of BGIDL for drain side" unit="V/K");
    parameter real    CGIDL    =  0.0                             `P(info="Back-bias dependence of GIDL" unit="");
    parameter real    CGIDLD   =  0.0                             `P(info="Back-bias dependence of GIDL for drain side" unit="");

    //  Charge model parameters
    parameter real    COX      =  1.0e-14  `from(    0.0,inf    ) `P(info="Oxide capacitance for intrinsic channel" unit="F");
    parameter real    CGOV     =  1.0e-15  `from(    0.0,inf    ) `P(info="Oxide capacitance for gate-drain/source overlap" unit="F");
    parameter real    CGOVD    =  1.0e-15  `from(    0.0,inf    ) `P(info="Oxide capacitance for gate-drain overlap" unit="F");
    parameter real    CGBOV    =  0.0      `from(    0.0,inf    ) `P(info="Oxide capacitance for gate-bulk overlap" unit="F");
    parameter real    CFR      =  0.0      `from(    0.0,inf    ) `P(info="Outer fringe capacitance" unit="F");
    parameter real    CFRD     =  0.0      `from(    0.0,inf    ) `P(info="Outer fringe capacitance for drain side" unit="F");

    //  Noise parameters
    parameter real    FNT      =  1.0      `from(    0.0,inf    ) `P(info="Thermal noise coefficient" unit="");
    parameter real    FNTEXC   =  0.0      `from(    0.0,inf    ) `P(info="Excess noise coefficient" unit="");
    parameter real    NFA      =  8.0e+22  `from(    0.0,inf    ) `P(info="First coefficient of flicker noise" unit="V^-1/m^4");
    parameter real    NFB      =  3.0e+07  `from(    0.0,inf    ) `P(info="Second coefficient of flicker noise" unit="V^-1/m^2");
    parameter real    NFC      =  0.0      `from(    0.0,inf    ) `P(info="Third coefficient of flicker noise" unit="V^-1");
    parameter real    EF       =  1.0      `from(    0.0,inf    ) `P(info="Flicker noise frequency exponent" unit="");

`ifdef NQSmodel
        //  NQS parameters
        parameter real    SWNQS    =  0.0      `from(    0.0,9.0    ) `P(info="Flag for NQS, 0=off, 1, 2, 3, 5, or 9=number of collocation points" unit="");
        parameter real    MUNQS    =  1.0      `from(    0.0,inf    ) `P(info="Relative mobility for NQS modelling");
`endif // NQSmodel

    // Parasitic resistance parameters
    parameter real    RG       =  0.0   `from( 0.0,inf    ) `P(info="Gate resistance" unit="Ohm");
    parameter real    RBULK    =  0.0   `from( 0.0,inf    ) `P(info="Bulk resistance between node BP and BI" unit="Ohm");
    parameter real    RWELL    =  0.0   `from( 0.0,inf    ) `P(info="Well resistance between node BI and B" unit="Ohm");
    parameter real    RJUNS    =  0.0   `from( 0.0,inf    ) `P(info="Source-side bulk resistance between node BI and BS" unit="Ohm");
    parameter real    RJUND    =  0.0   `from( 0.0,inf    ) `P(info="Drain-side bulk resistance between node BI and BD" unit="Ohm");

    // JUNCAP Parameters
    `include "JUNCAP200_parlist.include"

    //  Other parameters
    parameter real    DTA      =  0.0                             `P(info="Temperature offset w.r.t. ambient temperature" unit="K");

    //  Instance parameters
    parameter real    DELVTO     = 0.0                            `P(type="instance" info="Threshold voltage shift parameter" unit="V");
    parameter real    FACTUO     = 1.0     `from(    0.0,inf    ) `P(type="instance" info="Zero-field mobility pre-factor" unit="");
    parameter real    ABSOURCE   = 1e-12   `from(`AB_cliplow,inf) `P(type="instance" info="Bottom area of source junction" unit="m^2");
    parameter real    LSSOURCE   = 1e-6    `from(`LS_cliplow,inf) `P(type="instance" info="STI-edge length of source junction" unit="m");
    parameter real    LGSOURCE   = 1e-6    `from(`LG_cliplow,inf) `P(type="instance" info="Gate-edge length of source junction" unit="m");
    parameter real    ABDRAIN    = 1e-12   `from(`AB_cliplow,inf) `P(type="instance" info="Bottom area of drain junction" unit="m^2");
    parameter real    LSDRAIN    = 1e-6    `from(`LS_cliplow,inf) `P(type="instance" info="STI-edge length of drain junction" unit="m");
    parameter real    LGDRAIN    = 1e-6    `from(`LG_cliplow,inf) `P(type="instance" info="Gate-edge length of drain junction" unit="m");
    parameter real    AS         = 1E-12   `from(`AB_cliplow,inf) `P(type="instance" info="Bottom area of source junction" unit="m^2");
    parameter real    PS         = 1E-6    `from(`LS_cliplow,inf) `P(type="instance" info="Perimeter of source junction" unit="m");
    parameter real    AD         = 1E-12   `from(`AB_cliplow,inf) `P(type="instance" info="Bottom area of drain junction" unit="m^2");
    parameter real    PD         = 1E-6    `from(`LS_cliplow,inf) `P(type="instance" info="Perimeter of drain junction" unit="m");
    parameter real    JW         = 1E-6    `from(`LG_cliplow,inf) `P(type="instance" info="Gate-edge length of source/drain junction" unit="m");
    parameter real    MULT       = 1.0     `from(    0.0,inf    ) `P(type="instance" info="Number of devices in parallel" unit="");
`else // LocalModel
`ifdef Binning

    `include "PSP102_binpars.include"

`else // Binning
        ///////////////////////////////////////////////////
        // PSP global model parameters
        ///////////////////////////////////////////////////

        //  Special model parameters
        parameter real    LEVEL    =  1020                            `P(info="Model level" unit="");

        parameter integer TYPE     =  1        `from(   -1,  1      ) `P(info="Channel type parameter, +1=NMOS -1=PMOS" unit="");
        parameter integer nmos     =  1        `from(    0,  1      ) `P(desc="MOS channel type");
        parameter integer pmos     =  1        `from(    0,  1      ) `P(desc="MOS channel type");

        // Reference Temperature
        parameter real    TR       =  21.0     `from( -273.0,inf    ) `P(info="nominal (reference) temperature" unit="C");

        //  Switch parameters that turn models or effects on or off
        parameter real    SWIGATE  =  0.0      `from(    0.0,1.0    ) `P(info="Flag for gate current, 0=turn off IG" unit="");
        parameter real    SWIMPACT =  0.0      `from(    0.0,1.0    ) `P(info="Flag for impact ionization current, 0=turn off II" unit="");
        parameter real    SWGIDL   =  0.0      `from(    0.0,1.0    ) `P(info="Flag for GIDL current, 0=turn off IGIDL" unit="");
        parameter real    SWJUNCAP =  0.0      `from(    0.0,3.0    ) `P(info="Flag for juncap, 0=turn off juncap" unit="");
        parameter real    SWJUNASYM=  0.0                             `P(info="Flag for asymmetric junctions; 0=symmetric, 1=asymmetric" unit="");
        parameter real    QMC      =  1.0      `from(    0.0,inf    ) `P(info="Quantum-mechanical correction factor" unit="");

        // Process Parameters
        parameter real    LVARO    =  0.0                             `P(info="Geom. independent difference between actual and programmed gate length" unit="m");
        parameter real    LVARL    =  0.0                             `P(info="Length dependence of LVAR" unit="");
        parameter real    LVARW    =  0.0                             `P(info="Width dependence of LVAR" unit="");
        parameter real    LAP      =  0.0                             `P(info="Effective channel length reduction per side" unit="m");
        parameter real    WVARO    =  0.0                             `P(info="Geom. independent difference between actual and programmed field-oxide opening" unit="m");
        parameter real    WVARL    =  0.0                             `P(info="Length dependence of WVAR" unit="");
        parameter real    WVARW    =  0.0                             `P(info="Width dependence of WVAR" unit="");
        parameter real    WOT      =  0.0                             `P(info="Effective channel width reduction per side" unit="m");
        parameter real    DLQ      =  0.0                             `P(info="Effective channel length reduction for CV" unit="m");
        parameter real    DWQ      =  0.0                             `P(info="Effective channel width reduction for CV" unit="m");
        parameter real    VFBO     = -1.0                             `P(info="Geometry-independent flat-band voltage at TR" unit="V");
        parameter real    VFBL     =  0.0                             `P(info="Length dependence of flat-band voltage" unit="");
        parameter real    VFBW     =  0.0                             `P(info="Width dependence of flat-band voltage" unit="");
        parameter real    VFBLW    =  0.0                             `P(info="Area dependence of flat-band voltage" unit="");
        parameter real    STVFBO   =  5e-4                            `P(info="Geometry-independent temperature dependence of VFB" unit="V/K");
        parameter real    STVFBL   =  0.0                             `P(info="Length dependence of temperature dependence of VFB" unit="");
        parameter real    STVFBW   =  0.0                             `P(info="Width dependence of temperature dependence of VFB" unit="");
        parameter real    STVFBLW  =  0.0                             `P(info="Area dependence of temperature dependence of VFB" unit="");
        parameter real    TOXO     =  2e-9     `from(  1e-10,inf    ) `P(info="Gate oxide thickness" unit="m");
        parameter real    EPSROXO  =  3.9      `from(    1.0,inf    ) `P(info="Relative permittivity of gate dielectric" unit="");
        parameter real    NSUBO    =  3e23     `from(   1e20,inf    ) `P(info="Geometry independent substrate doping" unit="m^-3");
        parameter real    NSUBW    =  0.0                             `P(info="Width dependence of background doping NSUBO due to segregation" unit="");
        parameter real    WSEG     =  1e-8     `from(  1e-10,inf    ) `P(info="Char. length of segregation of background doping NSUBO" unit="m");
        parameter real    NPCK     =  1e24     `from(    0.0,inf    ) `P(info="Pocket doping level" unit="m^-3");
        parameter real    NPCKW    =  0.0                             `P(info="Width dependence of pocket doping NPCK due to segregation" unit="");
        parameter real    WSEGP    =  1e-8     `from(  1e-10,inf    ) `P(info="Char. length of segregation of pocket doping NPCK" unit="m");
        parameter real    LPCK     =  1e-8     `from(  1e-10,inf    ) `P(info="Char. length of lateral doping profile" unit="m");
        parameter real    LPCKW    =  0.0                             `P(info="Width dependence of char. length of lateral doping profile" unit="");
        parameter real    FOL1     =  0.0                             `P(info="First length dependence coefficient for short channel body effect" unit="");
        parameter real    FOL2     =  0.0                             `P(info="Second length dependence coefficient for short channel body effect" unit="");
        parameter real    VNSUBO   =  0.0                             `P(info="Effective doping bias-dependence parameter" unit="V");
        parameter real    NSLPO    =  0.05                            `P(info="Effective doping bias-dependence parameter" unit="V");
        parameter real    DNSUBO   =  0.0                             `P(info="Effective doping bias-dependence parameter" unit="V^-1");
        parameter real    DPHIBO   =  0.0                             `P(info="Geometry independent offset of PHIB" unit="V");
        parameter real    DPHIBL   =  0.0                             `P(info="Length dependence offset of PHIB" unit="V");
        parameter real    DPHIBLEXP=  1.0                             `P(info="Exponent for length dependence of offset of PHIB" unit="");
        parameter real    DPHIBW   =  0.0                             `P(info="Width dependence of offset of PHIB" unit="");
        parameter real    DPHIBLW  =  0.0                             `P(info="Area dependence of offset of PHIB" unit="");
        parameter real    NPO      =  1e26                            `P(info="Geometry-independent gate poly-silicon doping" unit="m^-3");
        parameter real    NPL      =  0.0                             `P(info="Length dependence of gate poly-silicon doping" unit="");
        parameter real    CTO      =  0.0                             `P(info="Geometry-independent interface states factor" unit="");
        parameter real    CTL      =  0.0                             `P(info="Length dependence of interface states factor" unit="");
        parameter real    CTLEXP   =  1.0                             `P(info="Exponent for length dependence of interface states factor" unit="");
        parameter real    CTW      =  0.0                             `P(info="Width dependence of interface states factor" unit="");
        parameter real    CTLW     =  0.0                             `P(info="Area dependence of interface states factor" unit="");
        parameter real    TOXOVO   =  2e-9     `from(  1e-10,inf    ) `P(info="Overlap oxide thickness" unit="m");
        parameter real    TOXOVDO  =  2e-9     `from(  1e-10,inf    ) `P(info="Overlap oxide thickness for drain side" unit="m");
        parameter real    LOV      =  0        `from(    0.0,inf    ) `P(info="Overlap length for gate/drain and gate/source overlap capacitance" unit="m");
        parameter real    LOVD     =  0        `from(    0.0,inf    ) `P(info="Overlap length for gate/drain overlap capacitance" unit="m");
        parameter real    NOVO     =  5e25                            `P(info="Effective doping of overlap region" unit="m^-3");
        parameter real    NOVDO    =  5e25                            `P(info="Effective doping of overlap region for drain side" unit="m^-3");

        // DIBL Parameters
        parameter real    CFL      =  0.0                             `P(info="Length dependence of DIBL-parameter" unit="");
        parameter real    CFLEXP   =  2.0                             `P(info="Exponent for length dependence of CF" unit="");
        parameter real    CFW      =  0.0                             `P(info="Width dependence of CF" unit="");
        parameter real    CFBO     =  0.0                             `P(info="Back-bias dependence of CF" unit="V^-1");

        // Mobility Parameters
        parameter real    UO       =  5e-2                            `P(info="Zero-field mobility at TR" unit="m^2/V/s");
        parameter real    FBET1    =  0.0                             `P(info="Relative mobility decrease due to first lateral profile" unit="");
        parameter real    FBET1W   =  0.0                             `P(info="Width dependence of relative mobility decrease due to first lateral profile" unit="");
        parameter real    LP1      =  1e-8     `from(  1e-10,inf    ) `P(info="Mobility-related characteristic length of first lateral profile" unit="m");
        parameter real    LP1W     =  0.0                             `P(info="Width dependence of mobility-related characteristic length of first lateral profile" unit="");
        parameter real    FBET2    =  0.0                             `P(info="Relative mobility decrease due to second lateral profile" unit="");
        parameter real    LP2      =  1e-8     `from(  1e-10,inf    ) `P(info="Mobility-related characteristic length of second lateral profile" unit="m");
        parameter real    BETW1    =  0.0                             `P(info="First higher-order width scaling coefficient of BETN" unit="");
        parameter real    BETW2    =  0.0                             `P(info="Second higher-order width scaling coefficient of BETN" unit="");
        parameter real    WBET     =  1e-9     `from(  1e-10,inf    ) `P(info="Characteristic width for width scaling of BETN" unit="m");
        parameter real    STBETO   =  1.0                             `P(info="Geometry independent temperature dependence of BETN" unit="");
        parameter real    STBETL   =  0.0                             `P(info="Length dependence of temperature dependence of BETN" unit="");
        parameter real    STBETW   =  0.0                             `P(info="Width dependence of temperature dependence of BETN" unit="");
        parameter real    STBETLW  =  0.0                             `P(info="Area dependence of temperature dependence of BETN" unit="");
        parameter real    MUEO     =  0.5                             `P(info="Geometry independent mobility reduction coefficient at TR" unit="m/V");
        parameter real    MUEW     =  0.0                             `P(info="Width dependence of mobility reduction coefficient at TR" unit="");
        parameter real    STMUEO   =  0.0                             `P(info="Temperature dependence of MUE" unit="");
        parameter real    THEMUO   =  1.5                             `P(info="Mobility reduction exponent at TR" unit="");
        parameter real    STTHEMUO =  1.5                             `P(info="Temperature dependence of THEMU" unit="");
        parameter real    CSO      =  0.0                             `P(info="Geometry independent coulomb scattering parameter at TR" unit="");
        parameter real    CSL      =  0.0                             `P(info="Length dependence of CS" unit="");
        parameter real    CSLEXP   =  0.0                             `P(info="Exponent for length dependence of CS" unit="");
        parameter real    CSW      =  0.0                             `P(info="Width dependence of CS" unit="");
        parameter real    CSLW     =  0.0                             `P(info="Area dependence of CS" unit="");
        parameter real    STCSO    =  0.0                             `P(info="Temperature dependence of CS" unit="");
        parameter real    XCORO    =  0.0                             `P(info="Geometry independent non-universality parameter" unit="V^-1");
        parameter real    XCORL    =  0.0                             `P(info="Length dependence of non-universality parameter" unit="");
        parameter real    XCORW    =  0.0                             `P(info="Width dependence of non-universality parameter" unit="");
        parameter real    XCORLW   =  0.0                             `P(info="Area dependence of non-universality parameter" unit="");
        parameter real    STXCORO  =  0.0                             `P(info="Temperature dependence of XCOR" unit="");
        parameter real    FETAO    =  1.0                             `P(info="Effective field parameter" unit="");

        // Series Resistance
        parameter real    RSW1     =  2.5e3                           `P(info="Source/drain series resistance for 1 um wide channel at TR" unit="Ohm");
        parameter real    RSW2     =  0.0                             `P(info="Higher-order width scaling of RS" unit="");
        parameter real    STRSO    =  1.0                             `P(info="Temperature dependence of RS" unit="");
        parameter real    RSBO     =  0.0                             `P(info="Back-bias dependence of series resistance" unit="V^-1");
        parameter real    RSGO     =  0.0                             `P(info="Gate-bias dependence of series resistance" unit="V^-1");

        // Velocity Saturation
        parameter real    THESATO  =  0.0                             `P(info="Geometry independent velocity saturation parameter at TR" unit="V^-1");
        parameter real    THESATL  =  0.05                            `P(info="Length dependence of THESAT" unit="V^-1");
        parameter real    THESATLEXP= 1.0                             `P(info="Exponent for length dependence of THESAT" unit="");
        parameter real    THESATW  =  0.0                             `P(info="Width dependence of velocity saturation parameter" unit="");
        parameter real    THESATLW =  0.0                             `P(info="Area dependence of velocity saturation parameter" unit="");
        parameter real    STTHESATO=  1.0                             `P(info="Geometry independent temperature dependence of THESAT" unit="");
        parameter real    STTHESATL=  0.0                             `P(info="Length dependence of temperature dependence of THESAT" unit="");
        parameter real    STTHESATW=  0.0                             `P(info="Width dependence of temperature dependence of THESAT" unit="");
        parameter real    STTHESATLW= 0.0                             `P(info="Area dependence of temperature dependence of THESAT" unit="");
        parameter real    THESATBO =  0.0                             `P(info="Back-bias dependence of velocity saturation" unit="V^-1");
        parameter real    THESATGO =  0.0                             `P(info="Gate-bias dependence of velocity saturation" unit="V^-1");

        // Saturation Voltage
        parameter real    AXO      =  18                              `P(info="Geometry independent linear/saturation transition factor" unit="");
        parameter real    AXL      =  0.4      `from(    0.0,inf    ) `P(info="Length dependence of AX" unit="");

        // Channel Length Modulation
        parameter real    ALPL     =  5e-4                            `P(info="Length dependence of ALP" unit="");
        parameter real    ALPLEXP  =  1.0                             `P(info="Exponent for length dependence of ALP" unit="");
        parameter real    ALPW     =  0.0                             `P(info="Width dependence of ALP" unit="");
        parameter real    ALP1L1   =  0.0                             `P(info="Length dependence of CLM enhancement factor above threshold" unit="V");
        parameter real    ALP1LEXP =  0.5                             `P(info="Exponent for length dependence of ALP1" unit="");
        parameter real    ALP1L2   =  0.0      `from(    0.0,inf    ) `P(info="Second_order length dependence of ALP1" unit="");
        parameter real    ALP1W    =  0.0                             `P(info="Width dependence of ALP1" unit="");
        parameter real    ALP2L1   =  0.0                             `P(info="Length dependence of CLM enhancement factor below threshold" unit="V^-1");
        parameter real    ALP2LEXP =  0.5                             `P(info="Exponent for length dependence of ALP2" unit="");
        parameter real    ALP2L2   =  0.0      `from(    0.0,inf    ) `P(info="Second_order length dependence of ALP2" unit="");
        parameter real    ALP2W    =  0.0                             `P(info="Width dependence of ALP2" unit="");
        parameter real    VPO      =  0.05                            `P(info="CLM logarithmic dependence parameter" unit="V");

        // Weak-avalanche parameters
        parameter real    A1O      =  1.0                             `P(info="Geometry independent impact-ionization pre-factor" unit="");
        parameter real    A1L      =  0.0                             `P(info="Length dependence of A1" unit="");
        parameter real    A1W      =  0.0                             `P(info="Width dependence of A1" unit="");
        parameter real    A2O      =  10                              `P(info="Impact-ionization exponent at TR" unit="V");
        parameter real    STA2O    =  0.0                             `P(info="Temperature dependence of A2" unit="V");
        parameter real    A3O      =  1.0                             `P(info="Geometry independent saturation-voltage dependence of II" unit="");
        parameter real    A3L      =  0.0                             `P(info="Length dependence of A3" unit="");
        parameter real    A3W      =  0.0                             `P(info="Width dependence of A3" unit="");
        parameter real    A4O      =  0.0                             `P(info="Geometry independent back-bias dependence of II" unit="V^-0.5");
        parameter real    A4L      =  0.0                             `P(info="Length dependence of A4" unit="");
        parameter real    A4W      =  0.0                             `P(info="Width dependence of A4" unit="");

        // Gate current parameters
        parameter real    GCOO     =  0.0                             `P(info="Gate tunnelling energy adjustment" unit="");
        parameter real    IGINVLW  =  0.0                             `P(info="Gate channel current pre-factor for 1 um^2 channel area" unit="A");
        parameter real    IGOVW    =  0.0                             `P(info="Gate overlap current pre-factor for 1 um wide channel" unit="A");
        parameter real    IGOVDW   =  0.0                             `P(info="Gate overlap current pre-factor for 1 um wide channel for drain side" unit="A");
        parameter real    STIGO    =  2.0                             `P(info="Temperature dependence of IGINV and IGOV" unit="");
        parameter real    GC2O     =  0.375                           `P(info="Gate current slope factor" unit="");
        parameter real    GC3O     =  0.063                           `P(info="Gate current curvature factor" unit="");
        parameter real    CHIBO    =  3.1                             `P(info="Tunnelling barrier height" unit="V");

        // Gate-induced drain leakage parameters
        parameter real    AGIDLW   =  0.0                             `P(info="Width dependence of GIDL pre-factor" unit="A/V^3");
        parameter real    AGIDLDW  =  0.0                             `P(info="Width dependence of GIDL pre-factor for drain side" unit="A/V^3");
        parameter real    BGIDLO   =  41                              `P(info="GIDL probability factor at TR" unit="V");
        parameter real    BGIDLDO  =  41                              `P(info="GIDL probability factor at TR for drain side" unit="V");
        parameter real    STBGIDLO =  0.0                             `P(info="Temperature dependence of BGIDL" unit="V/K");
        parameter real    STBGIDLDO=  0.0                             `P(info="Temperature dependence of BGIDL for drain side" unit="V/K");
        parameter real    CGIDLO   =  0.0                             `P(info="Back-bias dependence of GIDL" unit="");
        parameter real    CGIDLDO  =  0.0                             `P(info="Back-bias dependence of GIDL for drain side" unit="");

        // Charge Model Parameters
        parameter real    CGBOVL   =  0.0                             `P(info="Oxide capacitance for gate-bulk overlap for 1 um long channel" unit="F");
        parameter real    CFRW     =  0.0                             `P(info="Outer fringe capacitance for 1 um wide channel" unit="F");
        parameter real    CFRDW    =  0.0                             `P(info="Outer fringe capacitance for 1 um wide channel for drain side" unit="F");

        // Noise Model Parameters
        parameter real    FNTO     =  1.0                             `P(info="Thermal noise coefficient" unit="");
        parameter real    FNTEXCL  =  0.0      `from(    0.0,inf    ) `P(info="Length dependence coefficient of excess noise" unit="");
        parameter real    NFALW    =  8e22                            `P(info="First coefficient of flicker noise for 1 um^2 channel area" unit="V^-1/m^4");
        parameter real    NFBLW    =  3e7                             `P(info="Second coefficient of flicker noise for 1 um^2 channel area" unit="V^-1/m^2");
        parameter real    NFCLW    =  0.0                             `P(info="Third coefficient of flicker noise for 1 um^2 channel area" unit="V^-1");
        parameter real    EFO      =  1.0                             `P(info="Flicker noise frequency exponent" unit="");
        parameter real    LINTNOI  =  0.0                             `P(info="Length offset for flicker noise" unit="m");
        parameter real    ALPNOI   =  2.0                             `P(info="Exponent for length offset for flicker noise" unit="");

        // Other Parameters
        parameter real    DTA      =  0                               `P(info="Temperature offset w.r.t. ambient circuit temperature" unit="K");

        // Well proximity effect Parameters
        parameter real    KVTHOWEO =  0                               `P(info="Geometrical independent threshold shift parameter" unit="");
        parameter real    KVTHOWEL =  0                               `P(info="Length dependent threshold shift parameter" unit="");
        parameter real    KVTHOWEW =  0                               `P(info="Width dependent threshold shift parameter" unit="");
        parameter real    KVTHOWELW=  0                               `P(info="Area dependent threshold shift parameter" unit="");
        parameter real    KUOWEO   =  0                               `P(info="Geometrical independent mobility degradation factor" unit="");
        parameter real    KUOWEL   =  0                               `P(info="Length dependent mobility degradation factor" unit="");
        parameter real    KUOWEW   =  0                               `P(info="Width dependent mobility degradation factor" unit="");
        parameter real    KUOWELW  =  0                               `P(info="Area dependent mobility degradation factor" unit="");
`endif // Binning
`ifdef NQSmodel
        parameter real    SWNQS    =  0.0      `from(    0.0,9.0    ) `P(info="Flag for NQS, 0=off, 1, 2, 3, 5, or 9=number of collocation points" unit="");
        parameter real    MUNQSO   =  1.0                             `P(info="Relative mobility for NQS modelling");
`endif // NQSmodel

    // Parasitic resistance parameters
    parameter real    RGO      =  0.0                             `P(info="Gate resistance" unit="Ohm");
    parameter real    RBULKO   =  0.0                             `P(info="Bulk resistance between node BP and BI" unit="Ohm");
    parameter real    RWELLO   =  0.0                             `P(info="Well resistance between node BI and B" unit="Ohm");
    parameter real    RJUNSO   =  0.0                             `P(info="Source-side bulk resistance between node BI and BS" unit="Ohm");
    parameter real    RJUNDO   =  0.0                             `P(info="Drain-side bulk resistance between node BI and BD" unit="Ohm");
    parameter real    RINT     =  0.0      `from(    0.0,inf    ) `P(info="Contact resistance between silicide and ploy" unit="Ohm/Sqr");
    parameter real    RVPOLY   =  0.0      `from(    0.0,inf    ) `P(info="Vertical poly resistance" unit="Ohm/Sqr");
    parameter real    RSHG     =  0.0      `from(    0.0,inf    ) `P(info="Gate electrode diffusion sheet resistance" unit="Ohm/Sqr");
    parameter real    DLSIL    =  0.0                             `P(info="Silicide extension over the physical gate length" unit="m");

    // Stress Model Parameters
    parameter real    SAREF    =  1.0e-6   `from(   1e-9,inf    ) `P(info="Reference distance between OD-edge and poly from one side" unit="m");
    parameter real    SBREF    =  1.0e-6   `from(   1e-9,inf    ) `P(info="Reference distance between OD-edge and poly from other side" unit="m");
    parameter real    WLOD     =  0                               `P(info="Width parameter" unit="m");
    parameter real    KUO      =  0                               `P(info="Mobility degradation/enhancement coefficient" unit="m");
    parameter real    KVSAT    =  0        `from(   -1.0,1.0    ) `P(info="Saturation velocity degradation/enhancement coefficient" unit="m");
    parameter real    TKUO     =  0                               `P(info="Temperature dependence of KUO" unit="");
    parameter real    LKUO     =  0                               `P(info="Length dependence of KUO" unit="m^LLODKUO");
    parameter real    WKUO     =  0                               `P(info="Width dependence of KUO" unit="m^WLODKUO");
    parameter real    PKUO     =  0                               `P(info="Cross-term dependence of KUO" unit="m^(LLODKUO+WLODKUO)");
    parameter real    LLODKUO  =  0        `from(    0.0,inf    ) `P(info="Length parameter for UO stress effect" unit="");
    parameter real    WLODKUO  =  0        `from(    0.0,inf    ) `P(info="Width parameter for UO stress effect" unit="");
    parameter real    KVTHO    =  0                               `P(info="Threshold shift parameter" unit="Vm");
    parameter real    LKVTHO   =  0                               `P(info="Length dependence of KVTHO" unit="m^LLODVTH");
    parameter real    WKVTHO   =  0                               `P(info="Width dependence of KVTHO" unit="m^WLODVTH");
    parameter real    PKVTHO   =  0                               `P(info="Cross-term dependence of KVTHO" unit="m^(LLODVTH+WLODVTH)");
    parameter real    LLODVTH  =  0        `from(    0.0,inf    ) `P(info="Length parameter for VTH-stress effect" unit="");
    parameter real    WLODVTH  =  0        `from(    0.0,inf    ) `P(info="Width parameter for VTH-stress effect" unit="");
    parameter real    STETAO   =  0                               `P(info="eta0 shift factor related to VTHO change" unit="m");
    parameter real    LODETAO  =  1.0      `from(    0.0,inf    ) `P(info="eta0 shift modification factor for stress effect" unit="");

    // Well proximity effect Parameters
    parameter real    SCREF    =  1e-6     `from(    0.0,inf    ) `P(info="Distance between OD-edge and well edge of a reference device" unit="m");
    parameter real    WEB      =  0                               `P(info="Coefficient for SCB" unit="");
    parameter real    WEC      =  0                               `P(info="Coefficient for SCC" unit="");

    // JUNCAP Parameters
    `include "JUNCAP200_parlist.include"

    // Instance parameters
    parameter real    L        =  10e-6    `from(   1e-9,inf    ) `P(type="instance" info="Design length" unit="m");
    parameter real    W        =  10e-6    `from(   1e-9,inf    ) `P(type="instance" info="Design width" unit="m");
    parameter real    SA       =  0.0                             `P(type="instance" info="Distance between OD-edge and poly from one side" unit="m");
    parameter real    SB       =  0.0                             `P(type="instance" info="Distance between OD-edge and poly from other side" unit="m");
    parameter real    SD       =  0.0                             `P(type="instance" info="Distance between neighbouring fingers" unit="m");
    parameter real    SCA      =  0.0      `from(    0.0,inf    ) `P(type="instance" info="Integral of the first distribution function for scattered well dopants" unit="");
    parameter real    SCB      =  0.0      `from(    0.0,inf    ) `P(type="instance" info="Integral of the second distribution function for scattered well dopants" unit="");
    parameter real    SCC      =  0.0      `from(    0.0,inf    ) `P(type="instance" info="Integral of the third distribution function for scattered well dopants" unit="");
    parameter real    SC       =  0.0                             `P(type="instance" info="Distance between OD-edge and nearest well edge" unit="m");
    parameter real    DELVTO   =  0.0                             `P(type="instance" info="Threshold voltage shift parameter" unit="V");
    parameter real    FACTUO   =  1.0      `from(    0.0,inf    ) `P(type="instance" info="Zero-field mobility pre-factor" unit="");
    parameter real    ABSOURCE = 1E-12     `from(`AB_cliplow,inf) `P(type="instance" info="Bottom area of source junction" unit="m^2");
    parameter real    LSSOURCE = 1E-6      `from(`LS_cliplow,inf) `P(type="instance" info="STI-edge length of source junction" unit="m");
    parameter real    LGSOURCE = 1E-6      `from(`LG_cliplow,inf) `P(type="instance" info="Gate-edge length of source junction" unit="m");
    parameter real    ABDRAIN  = 1E-12     `from(`AB_cliplow,inf) `P(type="instance" info="Bottom area of drain junction" unit="m^2");
    parameter real    LSDRAIN  = 1E-6      `from(`LS_cliplow,inf) `P(type="instance" info="STI-edge length of drain junction" unit="m");
    parameter real    LGDRAIN  = 1E-6      `from(`LG_cliplow,inf) `P(type="instance" info="Gate-edge length of drain junction" unit="m");
    parameter real    AS       = 1E-12     `from(`AB_cliplow,inf) `P(type="instance" info="Bottom area of source junction" unit="m^2");
    parameter real    PS       = 1E-6      `from(`LS_cliplow,inf) `P(type="instance" info="Perimeter of source junction" unit="m");
    parameter real    AD       = 1E-12     `from(`AB_cliplow,inf) `P(type="instance" info="Bottom area of drain junction" unit="m^2");
    parameter real    PD       = 1E-6      `from(`LS_cliplow,inf) `P(type="instance" info="Perimeter of drain junction" unit="m");
    parameter real    MULT     =  1.0      `from(    0.0,inf    ) `P(type="instance" info="Number of devices in parallel" unit="");
    parameter real    NF       =  1.0      `from(    1.0,inf    ) `P(type="instance" info="Number of fingers" unit="");
    parameter real    NGCON    =  1.0      `from(    1.0,2.0    ) `P(type="instance" info="Number of gate contacts" unit="");
    parameter real    XGW      =  1.0E-7                          `P(type="instance" info="Distance from the gate contact to the channel edge" unit="m");

    //////////////////////////
    //
    //  Variables
    //
    //////////////////////////

    // Variables for geometrical scaling rules
    real    L_i, W_i, L_f, L_slif, W_f, SA_i, SB_i, SD_i, SC_i;
    real    LEN, WEN, iL, iW, delLPS, delWOD, LE, WE, iLE, iWE, Lcv, Wcv, LEcv, WEcv;

`ifdef Binning
        // Auxiliary variables for binning-rules
        real    iLEWE, iiLE, iiWE, iiLEWE, iiiLEWE;
        real    iLEcv, iiLEcv, iiWEcv, iiLEWEcv, iiiLEWEcv;
        real    iLcv, iiLcv, iiWcv, iiLWcv, iiiLWcv;
`else // Binning
        // Intermediate variables used for geometry-scaling
        real    NSUBO_i, WSEG_i, NPCK_i, WSEGP_i, LPCK_i, LOV_i, LOVD_i;
        real    LP1_i, LP2_i, WBET_i, AXL_i, ALP1L2_i, ALP2L2_i;
        real    NSUB, AA, BB, NSUB0e, NPCKe, LPCKe;
        real    FBET1e, LP1e, GPE, GWE, tmpx;
        real    Lnoi, Lred;
`endif // Binning

    // List of local parameters
    real    VFB, STVFB, TOX, EPSROX, NEFF, VNSUB, NSLP, DNSUB, DPHIB, NP, CT;
    real    TOXOV, TOXOVD, NOV, NOVD, CF, CFB;
    real    BETN, STBET, MUE, STMUE, THEMU, STTHEMU, CS, STCS, XCOR, STXCOR, FETA;
    real    RS, STRS, RSB, RSG;
    real    THESAT, STTHESAT, THESATB, THESATG;
    real    AX;
    real    ALP, ALP1, ALP2, VP;
    real    A1, A2, STA2, A3, A4;
    real    GCO, IGINV, IGOV, IGOVD, STIG, GC2, GC3, CHIB;
    real    AGIDL, BGIDL, STBGIDL, CGIDL;
    real    AGIDLD, BGIDLD, CGIDLD, STBGIDLD;
    real    COX, CGOV, CGBOV, CGOVD, CFR, CFRD;
    real    FNT, NFA, NFB, NFC, EF;
    real    FNTEXC;
`ifdef NQSmodel
        real    MUNQS;
`endif // NQSmodel
    real    RG, RBULK, RWELL, RJUNS, RJUND;

    real    NF_i, invNF;
    real    RINT_i, RVPOLY_i, RSHG_i, NGCON_i, XGW_i, XGWE;

    // Variables for stress-model
    real    SAREF_i, SBREF_i, KVSAT_i, LLODKUO_i, WLODKUO_i, LLODVTH_i, WLODVTH_i, LODETAO_i;
    real    Invsa, Invsb, Invsaref, Invsbref, Kstressu0, rhobeta, rhobetaref, Kstressvth0;
    real    temp0, temp00, templ, tempw, Lx, Wx;
    real    loop, tmpa, tmpb;

    // Variables for well proximity effect model
    real    SCA_i, SCB_i, SCC_i, SCREF_i, WEB_i, WEC_i, KVTHOWE, KUOWE;
`endif // LocalModel

real ABS_i, LSS_i, LGS_i, ABD_i, LSD_i, LGD_i;

// Variables used in electrical equations
real FACTUO_i, DELVTO_i, VFB_i, STVFB_i, TOX_i, EPSROX_i, NEFF_i, VNSUB_i, NSLP_i, DNSUB_i, NP_i, QMC_i, CT_i, TOXOV_i, TOXOVD_i, NOV_i, NOVD_i;
real CF_i, CFB_i, DPHIB_i;
real BETN_i, STBET_i, MUE_i, STMUE_i, THEMU_i, STTHEMU_i, CS_i, STCS_i, XCOR_i, STXCOR_i, FETA_i;
real RS_i, THER_i, STRS_i, RSB_i, RSG_i;
real THESAT_i, STTHESAT_i, THESATB_i, THESATG_i;
real AX_i, ALP_i, ALP1_i, ALP2_i, VP_i;
real A1_i, A2_i, STA2_i, A3_i, A4_i;
real GCO_i, IGINV_i, IGOV_i, IGOVD_i, STIG_i, GC2_i, GC3_i, CHIB_i;
real AGIDL_i, BGIDL_i, STBGIDL_i, CGIDL_i;
real AGIDLD_i, BGIDLD_i, STBGIDLD_i, CGIDLD_i;
real COX_i, CGOV_i, CGOVD_i, CGBOV_i, CFR_i, CFRD_i;
real FNT_i, NFA_i, NFB_i, NFC_i, EF_i;
real FNTEXC_i;
real TR_i, MULT_i;
real vth_i, vts_i;

real temp, temp1, temp2, tempM;
real help;
real jwcorr;

real TKR, TKD, TKD_sq, dT, rT, rTn;
real EPSOX, EPSSI;
real BET_i, phit, inv_phit, Eg, phibFac, CoxPrime, tox_sq;
real delVg, CoxovPrime, GOV, GOV2;
real np, kp, qq, qb0, dphibq, qlim2;
real E_eff0, eta_mu, BCH, BOV, inv_CHIB, GCQ, Dch, Dov;
real tf_bet, tf_mue, tf_cs, tf_xcor, tf_ther, tf_thesat, tf_ig;
real xi_ov, inv_xi_ov, x_mrg_ov, x1, inv_xg1, Vdsat_lim;
real nt, Cox_over_q;

real CoxovPrime_d, GOV_s, GOV_d, GOV2_s, GOV2_d;
real BOV_d, xi_ov_s, xi_ov_d, inv_xi_ov_s, inv_xi_ov_d, x_mrg_ov_s, x_mrg_ov_d, inv_xg1_s, inv_xg1_d;
real AGIDLs, AGIDLDs, BGIDLs, BGIDLDs, BGIDL_T, BGIDLD_T;

real phib, sqrt_phib, phix, aphi, bphi, phix1, phix2, G_0, phit1, inv_phit1, alpha_b;
real inv_VP, inv_AX, Sfl_prefac;

real Vgs, Vds, Vsb, Vsbstar;
real Vgb, Vgb1, Vgbstar, Vdb, Vdsx, Vsbx;

real VgsPrime, VgdPrime, VsbPrime, VdbPrime;

real Dnsub;
real Igidl, Igisl, Vtovd, Vtovs;
real x_s, sqm, alpha1, phi_inf, za, rhob;
real thesat1, wsat, ysat, zsat, r1, r2, dL, GdL, dL1, GR, QCLM;
real xgm, Voxm, qim1_1, xgs_ov, xgd_ov;
real Ux;
real mutau, nu, xn_s, delta_ns;
real Gf, Gf2, inv_Gf2, xi, inv_xi, margin;
real qeff;

real SP_xg1, SP_S_temp,SP_S_temp1,SP_S_temp2;
real SP_S_yg, SP_S_ysub, SP_S_y0, SP_S_a, SP_S_b, SP_S_c;
real SP_S_bx, SP_S_tau, SP_S_eta, SP_S_delta0, SP_S_delta1;
real SP_S_pC, SP_S_qC, SP_S_A_fac;
real SP_S_x1, SP_S_w, SP_S_xbar, SP_S_x0;
real SP_S_xi0, SP_S_xi1, SP_S_xi2;
real SP_OV_yg, SP_OV_z, SP_OV_eta, SP_OV_a, SP_OV_c;
real SP_OV_tau, SP_OV_D0, SP_OV_y0, SP_OV_xi, SP_OV_temp;
real SP_OV_p, SP_OV_q, SP_OV_w, SP_OV_Afac, SP_OV_xbar;
real SP_OV_x0, SP_OV_u;

real x_d, x_m, x_ds, Rxcor, delta_1s, xi0s, xi1s, xi2s, xi0d;
real Es, Em, Ed, Ds, Dm, Dd, Ps, xgs, qis, qbs, qbm, Eeffm, Vm;
real Phi_0, Phi_2, asat, Phi_0_2, Phi0_Phi2;
real Vdsat, xn_d, k_ds, Udse;
real Mutmp, Phi_sat, delta_nd;
real pC, qC, Pm;
real d0, D_bar, km, x_pm, xi_pd, p_pd, u_pd, q_pd;
real xs_ov, xd_ov, Vovs, Vovd, psi_t;
real zg, delVsat, TP, Dsi, Dgate, u0, u0_div_H, x, xsq, inv_x, ex, inv_ex, Ag, Bg, Sg;
real Fj, Fj2;
real N1, Nm1, Delta_N1, Sfl;
real t1, t2, sqt2, r, lc, lcinv2, g_ideal, CGeff, mid, mig, migid, c_igid, sqid, sqig;
real H0;

// excess noise variables
real fac_exc, temp2_exc, wsat_exc, temp_exc, thesat1_exc, zsat_exc, Gvsat_exc, gfac, Sidexc;

real shot_igs, shot_igd, shot_iavl;
real shot_igcsx, shot_igcdx, shot_igsov, shot_igdov;

real Igc0, igc, igcd_h;
real Igc, Igb;
real QI, QD, QB, QG, Qg, Qd, Qb, Qs, Qgs_ov, Qgd_ov;
real Qfgs, Qfgd, Qgb_ov;

real arg1, arg2max, arg2mina;

real RG_i, RBULK_i, RWELL_i, RJUNS_i, RJUND_i;
real ggate, gbulk, gwell, gjund, gjuns, nt0;
real rgatenoise, rbulknoise, rwellnoise, rjundnoise, rjunsnoise;

integer CHNL_TYPE;

`ifdef NQSmodel
    // Variables used in NQS-calculations
    real SWNQS_i, MUNQS_i;
    real Qp1_0, Qp2_0, Qp3_0, Qp4_0, Qp5_0, Qp6_0, Qp7_0, Qp8_0, Qp9_0;
    real fk1, fk2, fk3, fk4, fk5, fk6, fk7, fk8, fk9;

    real phi_p1, phi_p2, phi_p3;
    real phi_p4, phi_p5, phi_p6;
    real phi_p7, phi_p8, phi_p9;

    real Qp1, Qp2, Qp3;
    real Qp4, Qp5, Qp6;
    real Qp7, Qp8, Qp9;
    real Qp0, QpN;

    real QG_NQS, QS_NQS, QD_NQS;
    real pd, Gp, Gp2, a_factrp, marginp, x_sp, x_dp;

    real dfQi, fQi, dQis, dQis_1, d2Qis, dQbs, dQy, d2Qy, dpsy2;
    real ym, inorm, Tnorm, Qb_tmp, QbSIGN;
    real r_nqs, vnorm, vnorm_inv;
    real NQS_xg1, NQS_yg, NQS_z, NQS_eta, NQS_a, NQS_c, NQS_tau, NQS_D0, NQS_xi, NQS_p;
    real NQS_q, NQS_temp, NQS_A_fac, NQS_xbar, NQS_w, NQS_x0, NQS_u, NQS_y0;
    real xphi, fk0, thesat2, Fvsat;
    real temp3, temp4, temp5, temp6, temp7, temp8, temp9;
`endif // NQSmodel

// JUNCAP2 variables
`include "JUNCAP200_varlist1.include"
`include "JUNCAP200_varlist2.include"
real jnoisex_s, jnoisex_d;
real Vjun_s, Vjun_d;
`ifdef LocalModel
    real JW_i;
`endif // LocalModel

real ijunbot, ijunsti, ijungat, qjunbot, qjunsti, qjungat;
real alpha, eta_p, xitsb, FdL, Gmob, Gmob_dL, Gvsat, Gvsatinv, dps, qim, qim1;
real H, Iimpact, mavl, Ids, Vdse, Igsov, Igdov, Igcs, Igcd, xg, sigVds;
real COX_qm;
real ijun_s, ijunbot_s, ijunsti_s, ijungat_s;
real ijun_d, ijunbot_d, ijunsti_d, ijungat_d;
real qjun_s, qjunbot_s, qjunsti_s, qjungat_s;
real qjun_d, qjunbot_d, qjunsti_d, qjungat_d;
real jnoise_s, jnoise_d;
real tmpv, vjv;
real sf;

/////////////////////////////////////////////////////////////////////////////
//
// Variables for operating point info
//
/////////////////////////////////////////////////////////////////////////////

real id_op, is, ig, ib, P_D, facvsb, facvsb0, sig1k;

`OPdef(ctype     ,"Flag for channel type" ,"")
`OPdef(sdint     ,"Flag for source-drain interchange" ,"")

`OPdef(ise       ,"Total source current" ,"A")
`OPdef(ige       ,"Total gate current" ,"A")
`OPdef(ide       ,"Total drain current" ,"A")
`OPdef(ibe       ,"Total bulk current" ,"A")
`OPdef(ids       ,"Drain current, excl. avalanche, tunnel, GISL, GIDL, and junction currents" ,"A")
`OPdef(idb       ,"Drain to bulk current" ,"A")
`OPdef(isb       ,"Source to bulk current" ,"A")
`OPdef(igs       ,"Gate-source tunneling current" ,"A")
`OPdef(igd       ,"Gate-drain tunneling current" ,"A")
`OPdef(igb       ,"Gate-bulk tunneling current" ,"A")
`OPdef(igcs      ,"Gate-channel tunneling current (source component)" ,"A")
`OPdef(igcd      ,"Gate-channel tunneling current (drain component)" ,"A")
`OPdef(iavl      ,"Substrate current due to weak avelanche" ,"A")
`OPdef(igisl     ,"Gate-induced source leakage current" ,"A")
`OPdef(igidl     ,"Gate-induced drain leakage current" ,"A")

`OPdef(ijs       ,"Total source junction current" ,"A")
`OPdef(ijsbot    ,"Source junction current (bottom component)" ,"A")
`OPdef(ijsgat    ,"Source junction current (gate-edge component)" ,"A")
`OPdef(ijssti    ,"Source junction current (STI-edge component)" ,"A")
`OPdef(ijd       ,"Total drain junction current" ,"A")
`OPdef(ijdbot    ,"Drain junction current (bottom component)" ,"A")
`OPdef(ijdgat    ,"Drain junction current (gate-edge component)" ,"A")
`OPdef(ijdsti    ,"Drain junction current (STI-edge component)" ,"A")

`OPdef(vds       ,"Drain-source voltage" ,"V")
`OPdef(vgs       ,"Gate-source voltage" ,"V")
`OPdef(vsb       ,"Source-bulk voltage" ,"V")
`OPdef(vto       ,"Zero-bias threshold voltage" ,"V")
`OPdef(vts       ,"Threshold voltage including back bias effects" ,"V")
`OPdef(vth       ,"Threshold voltage including back bias and drain bias effects" ,"V")
`OPdef(vgt       ,"Effective gate drive voltage including back bias and drain bias effects" ,"V")
`OPdef(vdss      ,"Drain saturation voltage at actual bias" ,"V")
`OPdef(vsat      ,"Saturation limit" ,"V")

`OPdef(gm        ,"Transconductance" ,"1/Ohm")
`OPdef(gmb       ,"Substrate transconductance" ,"1/Ohm")
`OPdef(gds       ,"Output conductance" ,"1/Ohm")
`OPdef(gjs       ,"Source junction conductance" ,"1/Ohm")
`OPdef(gjd       ,"Drain junction conductance" ,"1/Ohm")

`OPdef(cdd       ,"Drain capacitance" ,"F")
`OPdef(cdg       ,"Drain-gate capacitance" ,"F")
`OPdef(cds       ,"Drain-source capacitance" ,"F")
`OPdef(cdb       ,"Drain-bulk capacitance" ,"F")
`OPdef(cgd       ,"Gate-drain capacitance" ,"F")
`OPdef(cgg       ,"Gate capacitance" ,"F")
`OPdef(cgs       ,"Gate-source capacitance" ,"F")
`OPdef(cgb       ,"Gate-bulk capacitance" ,"F")
`OPdef(csd       ,"Source-drain capacitance" ,"F")
`OPdef(csg       ,"Source-gate capacitance" ,"F")
`OPdef(css       ,"Source capacitance" ,"F")
`OPdef(csb       ,"Source-bulk capacitance" ,"F")
`OPdef(cbd       ,"Bulk-drain capacitance" ,"F")
`OPdef(cbg       ,"Bulk-gate capacitance" ,"F")
`OPdef(cbs       ,"Bulk-source capacitance" ,"F")
`OPdef(cbb       ,"Bulk capacitance" ,"F")
`OPdef(cgsol     ,"Total gate-source overlap capacitance" ,"F")
`OPdef(cgdol     ,"Total gate-drain overlap capacitance" ,"F")

`OPdef(cjs       ,"Total source junction capacitance" ,"F")
`OPdef(cjsbot    ,"Source junction capacitance (bottom component)" ,"F")
`OPdef(cjsgat    ,"Source junction capacitance (gate-edge component)" ,"F")
`OPdef(cjssti    ,"Source junction capacitance (STI-edge component)" ,"F")
`OPdef(cjd       ,"Total drain junction capacitance" ,"F")
`OPdef(cjdbot    ,"Drain junction capacitance (bottom component)" ,"F")
`OPdef(cjdgat    ,"Drain junction capacitance (gate-edge component)" ,"F")
`OPdef(cjdsti    ,"Drain junction capacitance (STI-edge component)" ,"F")

`OPdef(weff      ,"Effective channel width for geometrical models" ,"m")
`OPdef(leff      ,"Effective channel length for geometrical models" ,"m")
`OPdef(u         ,"Transistor gain" ,"")
`OPdef(rout      ,"Small-signal output resistance" ,"Ohm")
`OPdef(vearly    ,"Equivalent Early voltage" ,"V")
`OPdef(beff      ,"Gain factor" ,"A/V^2")
`OPdef(fug       ,"Unity gain frequency at actual bias" ,"Hz")
`OPdef(rg_op     ,"Gate resistance" ,"Ohm")

`OPdef(sfl       ,"Flicker noise current spectral density at 1 Hz" ,"A^2/Hz")
`OPdef(sqrtsff   ,"Input-referred RMS white noise voltage spectral density at 1 kHz" ,"V/sqrt(Hz)")
`OPdef(sqrtsfw   ,"Input-referred RMS white noise voltage spectral density" ,"V/sqrt(Hz)")
`OPdef(sid       ,"White noise current spectral density" ,"A^2/Hz")
`OPdef(sig       ,"Induced gate noise current spectral density at 1 Hz" ,"A^2/Hz")
`OPdef(cigid     ,"Imaginary part of correlation coefficient between Sig and Sid" ,"")
`OPdef(fknee     ,"Cross-over frequency above which white noise is dominant" ,"Hz")
`OPdef(sigs      ,"Gate-source current noise spectral density" ,"A^2/Hz")
`OPdef(sigd      ,"Gate-drain current noise spectral density" ,"A^2/Hz")
`OPdef(siavl     ,"Impact ionization current noise spectral density" ,"A^2/Hz")
`OPdef(ssi       ,"Total source junction current noise spectral density" ,"A^2/Hz")
`OPdef(sdi       ,"Total drain junction current noise spectral density" ,"A^2/Hz")

// local parameters after scaling, T-scaling, and clipping
`OPdef(lp_vfb      , "Local parameter VFB after T-scaling and clipping", "V")
`OPdef(lp_stvfb    , "Local parameter STVFB after clipping", "V/K")
`OPdef(lp_tox      , "Local parameter TOX after clipping", "m")
`OPdef(lp_epsrox   , "Local parameter EPSROX after clipping", "")
`OPdef(lp_neff     , "Local parameter NEFF after clipping", "m^-3")
`OPdef(lp_vnsub    , "Local parameter VNSUB after clipping", "V")
`OPdef(lp_nslp     , "Local parameter NSLP after clipping", "V")
`OPdef(lp_dnsub    , "Local parameter DNSUB after clipping", "V^-1")
`OPdef(lp_dphib    , "Local parameter DPHIB after clipping", "V")
`OPdef(lp_np       , "Local parameter NP after clipping", "m^-3")
`OPdef(lp_ct       , "Local parameter CT after clipping", "")
`OPdef(lp_toxov    , "Local parameter TOXOV after clipping", "m")
`OPdef(lp_toxovd   , "Local parameter TOXOVD after clipping", "m")
`OPdef(lp_nov      , "Local parameter NOV after clipping", "m^-3")
`OPdef(lp_novd     , "Local parameter NOVD after clipping", "m^-3")
`OPdef(lp_cf       , "Local parameter CF after clipping", "")
`OPdef(lp_cfb      , "Local parameter CFB after clipping", "V^-1")
`OPdef(lp_betn     , "Local parameter BETN after T-scaling and clipping", "m^2/(V s)")
`OPdef(lp_stbet    , "Local parameter STBET after clipping", "")
`OPdef(lp_mue      , "Local parameter MUE after T-scaling and clipping", "m/V")
`OPdef(lp_stmue    , "Local parameter STMUE after clipping", "")
`OPdef(lp_themu    , "Local parameter THEMU after T-scaling and clipping", "")
`OPdef(lp_stthemu  , "Local parameter STTHEMU after clipping", "")
`OPdef(lp_cs       , "Local parameter CS after T-scaling and clipping", "")
`OPdef(lp_stcs     , "Local parameter STCS after clipping", "")
`OPdef(lp_xcor     , "Local parameter XCOR after T-scaling and clipping", "V^-1")
`OPdef(lp_stxcor   , "Local parameter STXCOR after clipping", "")
`OPdef(lp_feta     , "Local parameter FETA after clipping", "")
`OPdef(lp_rs       , "Local parameter RS after T-scaling and clipping", "Ohm")
`OPdef(lp_strs     , "Local parameter STRS after clipping", "")
`OPdef(lp_rsb      , "Local parameter RSB after clipping", "V^-1")
`OPdef(lp_rsg      , "Local parameter RSG after clipping", "V^-1")
`OPdef(lp_thesat   , "Local parameter THESAT after T-scaling and clipping", "V^-1")
`OPdef(lp_stthesat , "Local parameter STTHESAT after clipping", "")
`OPdef(lp_thesatb  , "Local parameter THESATB after clipping", "V^-1")
`OPdef(lp_thesatg  , "Local parameter THESATG after clipping", "V^-1")
`OPdef(lp_ax       , "Local parameter AX after clipping", "")
`OPdef(lp_alp      , "Local parameter ALP after clipping", "")
`OPdef(lp_alp1     , "Local parameter ALP1 after clipping", "V")
`OPdef(lp_alp2     , "Local parameter ALP2 after clipping", "V^-1")
`OPdef(lp_vp       , "Local parameter VP after clipping", "V")
`OPdef(lp_a1       , "Local parameter A1 after clipping", "")
`OPdef(lp_a2       , "Local parameter A2 after T-scaling and clipping", "V")
`OPdef(lp_sta2     , "Local parameter STA2 after clipping", "")
`OPdef(lp_a3       , "Local parameter A3 after clipping", "")
`OPdef(lp_a4       , "Local parameter A4 after clipping", "1/sqrt(V)")
`OPdef(lp_gco      , "Local parameter GCO after clipping", "")
`OPdef(lp_iginv    , "Local parameter IGINV after T-scaling and clipping", "A")
`OPdef(lp_igov     , "Local parameter IGOV after T-scaling and clipping", "A")
`OPdef(lp_igovd    , "Local parameter IGOVD after T-scaling and clipping", "A")
`OPdef(lp_stig     , "Local parameter STIG after clipping", "")
`OPdef(lp_gc2      , "Local parameter GC2 after clipping", "")
`OPdef(lp_gc3      , "Local parameter GC3 after clipping", "")
`OPdef(lp_chib     , "Local parameter CHIB after clipping", "V")
`OPdef(lp_agidl    , "Local parameter AGIDL after clipping", "A/V^3")
`OPdef(lp_agidld   , "Local parameter AGIDLD after clipping", "A/V^3")
`OPdef(lp_bgidl    , "Local parameter BGIDL after T-scaling and clipping", "V")
`OPdef(lp_bgidld   , "Local parameter BGIDLD after T-scaling and clipping", "V")
`OPdef(lp_stbgidl  , "Local parameter STBGIDL after clipping", "V/K")
`OPdef(lp_stbgidld , "Local parameter STBGIDLD after clipping", "V/K")
`OPdef(lp_cgidl    , "Local parameter CGIDL after clipping", "")
`OPdef(lp_cgidld   , "Local parameter CGIDLD after clipping", "")
`OPdef(lp_cox      , "Local parameter COX after clipping", "F")
`OPdef(lp_cgov     , "Local parameter CGOV after clipping", "F")
`OPdef(lp_cgovd    , "Local parameter CGOVD after clipping", "F")
`OPdef(lp_cgbov    , "Local parameter CGBOV after clipping", "F")
`OPdef(lp_cfr      , "Local parameter CFR after clipping", "F")
`OPdef(lp_cfrd     , "Local parameter CFRD after clipping", "F")
`OPdef(lp_fnt      , "Local parameter FNT after clipping", "")
`OPdef(lp_fntexc   , "Local parameter FNTEXC after clipping", "")
`OPdef(lp_nfa      , "Local parameter NFA after clipping", "1/(V m^4)")
`OPdef(lp_nfb      , "Local parameter NFB after clipping", "1/(V m^2)")
`OPdef(lp_nfc      , "Local parameter NFC after clipping", "V^-1")
`OPdef(lp_ef       , "Local parameter EF after clipping", "")
`OPdef(lp_rg       , "Local parameter RG after clipping", "Ohm")
`OPdef(lp_rbulk    , "Local parameter RBULK after clipping", "Ohm")
`OPdef(lp_rwell    , "Local parameter RWELL after clipping", "Ohm")
`OPdef(lp_rjuns    , "Local parameter RJUNS after clipping", "Ohm")
`OPdef(lp_rjund    , "Local parameter RJUND after clipping", "Ohm")
`OPdef(tk          , "Device Temperature", "K")
`OPdef(cjosbot     , "Bottom component of total zero-bias source junction capacitance at device temperature", "F")
`OPdef(cjossti     , "STI-edge component of total zero-bias source junction capacitance at device temperature", "F")
`OPdef(cjosgat     , "Gate-edge component of total zero-bias source junction capacitance at device temperature", "F")
`OPdef(vbisbot     , "Built-in voltage of source-side bottom junction at device temperature", "V")
`OPdef(vbissti     , "Built-in voltage of source-side STI-edge junction at device temperature", "V")
`OPdef(vbisgat     , "Built-in voltage of source-side gate-edge junction at device temperature", "V")
`OPdef(idsatsbot   , "Total source-side bottom junction saturation current", "A")
`OPdef(idsatssti   , "Total source-side STI-edge junction saturation current", "A")
`OPdef(idsatsgat   , "Total source-side gate-edge junction saturation current", "A")
`OPdef(cjosbotd    , "Bottom component of total zero-bias drain junction capacitance at device temperature", "F")
`OPdef(cjosstid    , "STI-edge component of total zero-bias drain junction capacitance at device temperature", "F")
`OPdef(cjosgatd    , "Gate-edge component of total zero-bias drain junction capacitance at device temperature", "F")
`OPdef(vbisbotd    , "Built-in voltage of drain-side bottom junction at device temperature", "V")
`OPdef(vbisstid    , "Built-in voltage of drain-side STI-edge junction at device temperature", "V")
`OPdef(vbisgatd    , "Built-in voltage of drain-side gate-edge junction at device temperature", "V")
`OPdef(idsatsbotd  , "Total drain-side bottom junction saturation current", "A")
`OPdef(idsatsstid  , "Total drain-side STI-edge junction saturation current", "A")
`OPdef(idsatsgatd  , "Total drain-side gate-edge junction saturation current", "A")
`ifdef NQSmodel
    `OPdef(lp_munqs    , "Local parameter MUNQS after clipping", "")
`endif // NQSmodel

/////////////////////////////////////////////////////////////////////////////
//
//  Analog block with all calculations and contribs
//
/////////////////////////////////////////////////////////////////////////////

analog begin

    `INITIAL_MODEL
    begin : initial_model
        // Code independent of bias or instance parameters
        // This block needs to be evaluated only once

`ifdef LocalModel
            // Do nothing
`else // LocalModel
`ifdef Binning
                // There are no binning parameters that need clipping
`else // Binning
                // Clipping of global model parameters
                TOX_i      = `CLIP_LOW(TOXO, 1e-10);
                TOXOVD_i   = `CLIP_LOW(TOXOVDO, 1e-10);
                EPSROX_i   = `CLIP_LOW(EPSROXO, 1.0);
                NSUBO_i    = `CLIP_LOW(NSUBO, 1e20);
                WSEG_i     = `CLIP_LOW(WSEG, 1e-10);
                NPCK_i     = `CLIP_LOW(NPCK, 0.0);
                WSEGP_i    = `CLIP_LOW(WSEGP, 1e-10);
                LPCK_i     = `CLIP_LOW(LPCK, 1e-10);
                TOXOV_i    = `CLIP_LOW(TOXOVO, 1e-10);
                LOV_i      = `CLIP_LOW(LOV, 0.0);
                LOVD_i     = `CLIP_LOW(LOVD, 0.0);
                LP1_i      = `CLIP_LOW(LP1, 1e-10);
                LP2_i      = `CLIP_LOW(LP2, 1e-10);
                WBET_i     = `CLIP_LOW(WBET, 1e-10);
                AXL_i      = `CLIP_LOW(AXL, 0.0);
                ALP1L2_i   = `CLIP_LOW(ALP1L2, 0.0);
                ALP2L2_i   = `CLIP_LOW(ALP2L2, 0.0);
`endif // Binning

            SAREF_i    = `CLIP_LOW(SAREF, 1e-9);
            SBREF_i    = `CLIP_LOW(SBREF, 1e-9);
            KVSAT_i    = `CLIP_BOTH(KVSAT, -1.0, 1.0);
            LLODKUO_i  = `CLIP_LOW(LLODKUO, 0.0);
            WLODKUO_i  = `CLIP_LOW(WLODKUO, 0.0);
            LLODVTH_i  = `CLIP_LOW(LLODVTH, 0.0);
            WLODVTH_i  = `CLIP_LOW(WLODVTH, 0.0);
            LODETAO_i  = `CLIP_LOW(LODETAO, 0.0);
            SCREF_i    = `CLIP_LOW(SCREF, 0.0);
            WEB_i      = WEB;
            WEC_i      = WEC;
            RINT_i     = `CLIP_LOW(RINT, 0.0);
            RVPOLY_i   = `CLIP_LOW(RVPOLY, 0.0);
`endif // LocalModel

        // 4.1 Internal parameters (including temperature scaling)
        // (only internal parameters independent on instance parameters
        //  are calculated in this section)
        if (`PGIVEN(nmos)) begin
            CHNL_TYPE = `NMOS;
        end else if (`PGIVEN(pmos)) begin
            CHNL_TYPE = `PMOS;
        end else begin
            CHNL_TYPE = (`PGIVEN(TYPE)) ? TYPE : `NMOS;
        end
        //$strobe("CHNL_TYPE %d", CHNL_TYPE);

        // Transistor temperature
        TR_i       =  `CLIP_LOW(TR, -273);
        TKR        =  `KELVINCONVERSION + TR_i;
        TKD   = $temperature + DTA;
        TKD_sq     =  TKD * TKD;
        dT         =  TKD - TKR;
        rT         =  TKD / TKR;
        rTn        =  TKR / TKD;
        phit       =  TKD * `KBOL / `QELE;
        inv_phit   =  1.0 / phit;

        // Local process parameters
        Eg         =  1.179 - 9.025e-5 * TKD - 3.05e-7 * TKD_sq;
        phibFac    =  (1.045 + 4.5e-4 * TKD) * (0.523 + 1.4e-3 * TKD - 1.48e-6 * TKD_sq) * TKD_sq / 9.0E4;
        phibFac    =  `MAX(phibFac, 1.0E-3);
        EPSSI      =  `EPSO * `EPSRSI;

`ifdef NQSmodel
            // Round SWNQS to nearest allowed value
            if (SWNQS < 0.5) begin
                SWNQS_i = 0.0;
            end else begin
                if (SWNQS < 1.5) begin
                    SWNQS_i = 1.0;
                end else begin
                    if (SWNQS < 2.5) begin
                        SWNQS_i = 2.0;
                    end else begin
                        if (SWNQS < 4.0) begin
                            SWNQS_i = 3.0;
                        end else begin
                            if (SWNQS < 7.0) begin
                                SWNQS_i = 5.0;
                            end else begin
                                SWNQS_i = 9.0;
                            end
                        end
                    end
                end
            end
            inorm      =  1.0e-12;
            r_nqs      =  1.0e+3;
            vnorm      =  10.0;
            vnorm_inv  =  1.0 / vnorm;

`endif // NQSmodel

        nt0        =  4 * `KBOL * TKD; // parameter for white noise of parasitic resistances

        // JUNCAP2
        `include "JUNCAP200_InitModel.include"

    end // initial_model

    `INITIAL_INSTANCE
    begin : initial_instance
        // Code independent of bias, but dependent on instance parameters,
        //   (including code dependent on parameters which could IN PRINCIPLE be scaled)
        // This block needs to be evaluated only once for each instance


        ABSOURCE_i = ABSOURCE;
        LSSOURCE_i = LSSOURCE;
        LGSOURCE_i = LGSOURCE;
        ABDRAIN_i  = ABDRAIN;
        LSDRAIN_i  = LSDRAIN;
        LGDRAIN_i  = LGDRAIN;
        AS_i       = AS;
        PS_i       = PS;
        AD_i       = AD;
        PD_i       = PD;


`ifdef LocalModel
            JW_i       = JW;
`else // LocalModel

            L_i        = L;
            W_i        = W;
            SA_i       = SA;
            SB_i       = SB;
            SD_i       = SD;
            SC_i       = SC;
            XGW_i      = XGW;

            // Clipping of the instance parameters
            NF_i       = `CLIP_LOW(NF, 1.0);
            NF_i       = floor(NF_i + 0.5); // round to nearest integer
            invNF      = 1.0 / NF_i;
            L_i        = `CLIP_LOW(L_i, 1e-9);
            W_i        = `CLIP_LOW(W_i * invNF, 1e-9);
            SCA_i      = `CLIP_LOW(SCA, 0.0);
            SCB_i      = `CLIP_LOW(SCB, 0.0);
            SCC_i      = `CLIP_LOW(SCC, 0.0);
            RSHG_i     = `CLIP_LOW(RSHG, 0.0);
            NGCON_i    = (NGCON < 1.5) ? 1.0 : 2.0;

            ///////////////////////////////////////////
            //  GEOMETRICAL PARAMETERSCALING
            ///////////////////////////////////////////

            // 3.2 Transistor geometry
            LEN        = 1e-6;
            WEN        = 1e-6;
            iL         = LEN / L_i;
            iW         = WEN / W_i;
`ifdef Binning
                delLPS     = LVARO * (1.0 + LVARL * iL);
                delWOD     = WVARO * (1.0 + WVARW * iW);
`else // Binning
                delLPS     = LVARO * (1.0 + LVARL * iL) * (1.0 + LVARW * iW);
                delWOD     = WVARO * (1.0 + WVARL * iL) * (1.0 + WVARW * iW);
`endif // Binning
            LE         = `CLIP_LOW(L_i + delLPS - 2.0 * LAP, 1e-9);
            WE         = `CLIP_LOW(W_i + delWOD - 2.0 * WOT, 1e-9);
            LEcv       = `CLIP_LOW(L_i + delLPS - 2.0 * LAP + DLQ, 1e-9);
            WEcv       = `CLIP_LOW(W_i + delWOD - 2.0 * WOT + DWQ, 1e-9);
            Lcv        = `CLIP_LOW(L_i + delLPS + DLQ, 1e-9);
            Wcv        = `CLIP_LOW(W_i + delWOD + DWQ, 1e-9);
            iLE        = LEN / LE;
            iWE        = WEN / WE;
            jwcorr     = 0.0;
            ABS_i      = 0.0;
            LSS_i      = 0.0;
            LGS_i      = 0.0;
            ABD_i      = 0.0;
            LSD_i      = 0.0;
            LGD_i      = 0.0;
            if (SWJUNCAP == 3.0) begin
                jwcorr = WE;
            end
            if ((SWJUNCAP == 2.0) || (SWJUNCAP == 3.0)) begin
                ABS_i = AS_i * invNF;
                LSS_i = PS_i * invNF - jwcorr;
                LGS_i = WE;
                ABD_i = AD_i * invNF;
                LSD_i = PD_i * invNF - jwcorr;
                LGD_i = WE;
            end else begin
                ABS_i = ABSOURCE_i * invNF;
                LSS_i = LSSOURCE_i * invNF;
                LGS_i = LGSOURCE_i * invNF;
                ABD_i = ABDRAIN_i * invNF;
                LSD_i = LSDRAIN_i * invNF;
                LGD_i = LGDRAIN_i * invNF;
            end

            // Geometry for multi-finger devices
            L_f        = `CLIP_LOW(L_i + delLPS, 1e-9);
            L_slif     = `CLIP_LOW(L_f + DLSIL, 1e-9);
            W_f        = `CLIP_LOW(W_i + delWOD, 1e-9);
            XGWE       = `CLIP_LOW(XGW_i - 0.5 * delWOD, 1e-9);

`ifdef Binning
                // 3.3 Geometry scaling with binning scaling rules
            `include "PSP102_binning.include"

`else // Binning
                // 3.2 Geometry scaling with physical scaling rules

                // Process parameters
                VFB        = VFBO * (1.0 + VFBL * iLE) * (1.0 + VFBW * iWE) * (1.0 + VFBLW * iLE * iWE);
                STVFB      = STVFBO * (1.0 + STVFBL * iLE) * (1.0 + STVFBW * iWE) * (1.0 + STVFBLW * iLE * iWE);
                TOX        = TOXO;
                EPSROX     = EPSROXO;
                NSUB0e     = NSUBO_i * `MAX(( 1.0 + NSUBW * iWE * ln( 1.0 + WE / WSEG_i )), 1.0E-03);
                NPCKe      = NPCK_i * `MAX(( 1.0 + NPCKW * iWE * ln( 1.0 + WE / WSEGP_i )), 1.0E-03);
                LPCKe      = LPCK_i * `MAX(( 1.0 + LPCKW * iWE * ln( 1.0 + WE / WSEGP_i )), 1.0E-03);
                if (LE > (2 * LPCKe)) begin
                    AA         = 7.5e10;
                    BB         = sqrt(NSUB0e + 0.5 * NPCKe) - sqrt(NSUB0e);
                    NSUB       = sqrt(NSUB0e) + AA * ln(1 + 2 * LPCKe / LE * (exp(BB / AA) - 1));
                    NSUB       = NSUB * NSUB;
                end else begin
                    if (LE >= LPCKe) begin
                        NSUB       = NSUB0e + NPCKe * LPCKe / LE;
                    end else begin // LE < LPCK
                        NSUB       = NSUB0e + NPCKe * (2 - LE / LPCKe);
                    end
                end
                NEFF       = NSUB * (1 - FOL1 * iLE - FOL2 * iLE * iLE);
                VNSUB      = VNSUBO;
                NSLP       = NSLPO;
                DNSUB      = DNSUBO;
                DPHIB      = (DPHIBO + DPHIBL * pow(iLE, DPHIBLEXP)) * (1.0 + DPHIBW * iWE) * (1.0 + DPHIBLW * iLE * iWE);
                NP         = NPO * `MAX(1e-6, (1.0 + NPL * iLE));
                CT         = (CTO + CTL * pow(iLE, CTLEXP)) * (1.0 + CTW * iWE) * (1.0 + CTLW * iLE * iWE);
                TOXOV      = TOXOVO;
                TOXOVD     = TOXOVDO;
                NOV        = NOVO;
                NOVD       = NOVDO;

                // DIBL parameters
                CF         = CFL * pow(iLE, CFLEXP) * (1.0 + CFW * iWE);
                CFB        = CFBO;

                // Mobility parameters
                FBET1e     = FBET1 * (1.0 + FBET1W * iWE);
                LP1e       = LP1_i * `MAX(1.0 + LP1W * iWE, 1.0E-03);
                GPE        = 1.0 + FBET1e * LP1e / LE * (1.0 - exp(-LE / LP1e)) + FBET2 * LP2_i / LE * (1.0 - exp(-LE / LP2_i));
                GPE        = `MAX(GPE, 1e-15);
                GWE        = 1.0 + BETW1 * iWE + BETW2 * iWE * ln(1.0 + WE / WBET_i);
                BETN       = UO * WE / (GPE * LE) * GWE;
                STBET      = STBETO * (1.0 + STBETL * iLE) * (1.0 + STBETW * iWE) * (1.0 + STBETLW * iLE * iWE);
                MUE        = MUEO * (1.0 + MUEW * iWE);
                STMUE      = STMUEO;
                THEMU      = THEMUO;
                STTHEMU    = STTHEMUO;
                CS         = (CSO + CSL * pow(iLE, CSLEXP)) * (1.0 + CSW * iWE) * (1.0 + CSLW * iLE * iWE);
                STCS       = STCSO;
                XCOR       = XCORO * (1.0 + XCORL * iLE) * (1.0 + XCORW * iWE) * (1.0 + XCORLW * iLE * iWE);
                STXCOR     = STXCORO;
                FETA       = FETAO;

                // Series resistance
                RS         = RSW1 * iWE * (1.0 + RSW2 * iWE);
                STRS       = STRSO;
                RSB        = RSBO;
                RSG        = RSGO;

                // Velocity saturation
                THESAT     = (THESATO + THESATL* GWE / GPE * pow(iLE, THESATLEXP)) * (1.0 + THESATW * iWE) * (1.0 + THESATLW * iLE * iWE);
                STTHESAT   = STTHESATO * (1.0 + STTHESATL * iLE) * (1.0 + STTHESATW * iWE) * (1.0 + STTHESATLW * iLE * iWE);
                THESATB    = THESATBO;
                THESATG    = THESATGO;

                // Saturation voltage
                AX         = AXO / (1.0 + AXL_i * iLE);

                // Channel length modulation
                ALP        = ALPL * pow(iLE, ALPLEXP) * (1.0 + ALPW * iWE);
                tmpx       = pow(iLE, ALP1LEXP);
                ALP1       = ALP1L1 * tmpx * (1.0 + ALP1W * iWE) / (1.0 + ALP1L2_i * iLE * tmpx);
                tmpx       = pow(iLE, ALP2LEXP);
                ALP2       = ALP2L1 * tmpx * (1.0 + ALP2W * iWE) / (1.0 + ALP2L2_i * iLE * tmpx);
                VP         = VPO;

                // Impact ionization
                A1         = A1O * (1.0 + A1L * iLE) * (1.0 + A1W * iWE);
                A2         = A2O;
                STA2       = STA2O;
                A3         = A3O * (1.0 + A3L * iLE) * (1.0 + A3W * iWE);
                A4         = A4O * (1.0 + A4L * iLE) * (1.0 + A4W * iWE);

                // Gate current
                GCO        = GCOO;
                IGINV      = IGINVLW / (iWE * iLE);
                IGOV       = IGOVW * LOV_i / (LEN * iWE);
                IGOVD      = IGOVDW * LOVD_i / (LEN * iWE);
                STIG       = STIGO;
                GC2        = GC2O;
                GC3        = GC3O;
                CHIB       = CHIBO;

                // GIDL
                AGIDL      = AGIDLW * LOV_i / (LEN * iWE);
                AGIDLD     = AGIDLDW * LOVD_i / (LEN * iWE);
                BGIDL      = BGIDLO;
                BGIDLD     = BGIDLDO;
                STBGIDL    = STBGIDLO;
                STBGIDLD   = STBGIDLDO;
                CGIDL      = CGIDLO;
                CGIDLD     = CGIDLDO;

                // Charge model parameters
                COX        = `EPSO * EPSROX_i * WEcv * LEcv / TOX_i;
                CGOV       = `EPSO * EPSROX_i * WEcv * LOV_i / TOXOV_i;
                CGOVD      = `EPSO * EPSROX_i * WEcv * LOVD_i / TOXOVD_i;
                CGBOV      = CGBOVL * Lcv / LEN;
                CFR        = CFRW * Wcv / WEN;
                CFRD       = CFRDW * Wcv / WEN;

                // Noise model parameters
                temp0      = 1.0 - 2.0 * LINTNOI * iLE / LEN;
                Lnoi       = `MAX(temp0, 1.0e-3);
                Lred       = 1.0 / pow(Lnoi, ALPNOI);

                FNT        = FNTO;
                FNTEXC     = FNTEXCL * BETN * BETN * iWE * iWE;
                NFA        = Lred * iWE * iLE * NFALW;
                NFB        = Lred * iWE * iLE * NFBLW;
                NFC        = Lred * iWE * iLE * NFCLW;
                EF         = EFO;

                // Well proximity effect parameters
                KVTHOWE    = KVTHOWEO + KVTHOWEL * iLE + KVTHOWEW * iWE + KVTHOWELW * iLE * iWE;
                KUOWE      = KUOWEO + KUOWEL * iLE + KUOWEW * iWE + KUOWELW * iLE * iWE;
`endif // Binning

            // Parasitic resistance parameters
            RG         = RSHG_i * (`oneThird * W_f / NGCON_i + XGWE) / (NGCON_i * L_slif)
                         + (RINT_i + RVPOLY_i) / (W_f * L_f) + NF_i * RGO;
            RWELL      = NF_i * RWELLO;
            RBULK      = NF_i * RBULKO;
            RJUNS      = NF_i * RJUNSO;
            RJUND      = NF_i * RJUNDO;

`ifdef NQSmodel
                MUNQS      = MUNQSO;
`endif // NQSModel

            ///////////////////////////////////////////
            //  STRESSMODEL
            ///////////////////////////////////////////

            // 3.4 Stress equations
            tmpa     = 0.0;
            tmpb     = 0.0;
            loop     = 0.0;
            if ((SA_i > 0.0) && (SB_i > 0.0) && ((NF_i == 1.0) || ((NF_i > 1.0) && (SD_i > 0.0)))) begin
                // Auxiliary variables

                // Note: some verilog-A compilers will (unnecesarily) cause the while-loop
                // below to be executed at every bias step; this has a negative impact on
                // the simulation speed of PSP.
                while (loop < (NF_i - 0.5)) begin
                    tmpa    = tmpa + 1.0 / (SA_i + 0.5 * L_i + loop * (SD_i + L_i));
                    tmpb    = tmpb + 1.0 / (SB_i + 0.5 * L_i + loop * (SD_i + L_i));
                    loop    = loop + 1.0;
                end
                Invsa      = tmpa * invNF;
                Invsb      = tmpb * invNF;
                Invsaref   = 1.0 / (SAREF_i + 0.5 * L_i);
                Invsbref   = 1.0 / (SBREF_i + 0.5 * L_i);
                Lx         = `MAX(L_i + delLPS, 1e-9);
                Wx         = `MAX(W_i + delWOD + WLOD, 1e-9);
                templ      =  1.0 / pow(Lx, LLODKUO_i);
                tempw      =  1.0 / pow(Wx, WLODKUO_i);
                Kstressu0  = (1.0 + LKUO * templ + WKUO * tempw + PKUO * templ * tempw) * (1.0 + TKUO * (rT - 1.0));
                rhobeta    = KUO * (Invsa + Invsb) / Kstressu0;
                rhobetaref = KUO * (Invsaref + Invsbref) / Kstressu0;
                templ      = 1.0 / pow(Lx, LLODVTH_i);
                tempw      = 1.0 / pow(Wx, WLODVTH_i);
                Kstressvth0= 1.0 + LKVTHO * templ + WKVTHO * tempw + PKVTHO * templ * tempw;
                temp0      = Invsa + Invsb - Invsaref - Invsbref;

                //  Parameter adaptations
                BETN       = BETN * (1.0 + rhobeta) / (1.0 + rhobetaref);
                THESAT     = THESAT * (1.0 + rhobeta) * (1.0 + KVSAT_i * rhobetaref) / ((1.0 + rhobetaref) * (1.0  + KVSAT_i * rhobeta));
                VFB        = VFB + KVTHO * temp0 / Kstressvth0;
                CF         = CF + STETAO * temp0 / pow(Kstressvth0, LODETAO_i);
            end

            ///////////////////////////////////////////
            //  WELL PROXIMITY EFFECT MODEL
            ///////////////////////////////////////////

            // 3.5 Well proximity effect equations
            if ((SCA_i > 0.0) || (SCB_i > 0.0) || (SCC_i > 0.0) || (SC_i > 0.0)) begin
                if ((SCA_i == 0.0) && (SCB_i == 0.0) && (SCC_i == 0.0)) begin
                    temp0      = SC_i + W_i;
                    temp00     = 1.0 / SCREF_i;
                    SCA_i      = SCREF_i * SCREF_i / (SC_i * temp0);
                    SCB_i      = ((0.1 * SC_i + 0.01 * SCREF_i) * exp(-10.0 * SC_i * temp00)
                                 - (0.1 * temp0 + 0.01 * SCREF_i) * exp(-10.0 * temp0 * temp00)) / W_i;
                    SCC_i      = ((0.05 * SC_i + 0.0025 * SCREF_i) * exp(-20.0 * SC_i * temp00)
                                 - (0.05 * temp0 + 0.0025 * SCREF_i) * exp(-20.0 * temp0 * temp00)) / W_i;
                end

                // Parameter adaptations
                temp0     = SCA_i + WEB_i * SCB_i + WEC_i * SCC_i;
                VFB       = VFB + KVTHOWE * temp0;
                BETN      = BETN * (1.0 + KUOWE * temp0);
            end

            /////////////////////////////////////////////
            //  END OF SCALINGRULES AND STRESS/WPE MODELS
            /////////////////////////////////////////////

`endif // LocalModel
        // 4.1 Internal parameters (including temperature scaling)

        // Clipping of the local model parameters
        VFB_i      =  VFB;
        STVFB_i    =  STVFB;
        TOX_i      = `CLIP_LOW(TOX, 1e-10);
        EPSROX_i   = `CLIP_LOW(EPSROX, 1.0);
        NEFF_i     = `CLIP_BOTH(NEFF,  1e20,  1e26);
        VNSUB_i    =  VNSUB;
        NSLP_i     = `CLIP_LOW(NSLP, 1e-3);
        DNSUB_i    = `CLIP_BOTH(DNSUB, 0.0, 1.0);
        DPHIB_i    =  DPHIB;
        NP_i       = `CLIP_LOW(NP, 0.0);
        QMC_i      = `CLIP_LOW(QMC, 0.0);
        CT_i       = `CLIP_LOW(CT, 0.0);
        TOXOV_i    = `CLIP_LOW(TOXOV, 1e-10);
        TOXOVD_i   = `CLIP_LOW(TOXOVD, 1e-10);
        NOV_i      = `CLIP_BOTH(NOV, 1e20, 1e27);
        NOVD_i     = `CLIP_BOTH(NOVD, 1e20, 1e27);
        CF_i       = `CLIP_LOW(CF, 0.0);
        CFB_i      = `CLIP_BOTH(CFB, 0.0, 1.0);
        BETN_i      = `CLIP_LOW(BETN, 0.0);
        STBET_i    =  STBET;
        MUE_i      = `CLIP_LOW(MUE, 0.0);
        STMUE_i    =  STMUE;
        THEMU_i    = `CLIP_LOW(THEMU, 0.0);
        STTHEMU_i  =  STTHEMU;
        CS_i       = `CLIP_LOW(CS,  0.0);
        STCS_i     =  STCS;
        XCOR_i     = `CLIP_LOW(XCOR, 0.0);
        STXCOR_i   =  STXCOR;
        FETA_i     = `CLIP_LOW(FETA, 0.0);
        RS_i       = `CLIP_LOW(RS, 0.0);
        STRS_i     =  STRS;
        RSB_i      = `CLIP_BOTH(RSB, -0.5, 1.0);
        RSG_i      = `CLIP_LOW(RSG, -0.5);
        THESAT_i   = `CLIP_LOW(THESAT, 0.0);
        STTHESAT_i =  STTHESAT;
        THESATB_i  = `CLIP_BOTH(THESATB, -0.5, 1.0);
        THESATG_i  = `CLIP_LOW(THESATG, -0.5);
        AX_i       = `CLIP_LOW(AX, 2.0);
        ALP_i      = `CLIP_LOW(ALP, 0.0);
        ALP1_i     = `CLIP_LOW(ALP1, 0.0);
        ALP2_i     = `CLIP_LOW(ALP2, 0.0);
        VP_i       = `CLIP_LOW(VP, 1.0e-10);
        A1_i       = `CLIP_LOW(A1, 0.0);
        A2_i       = `CLIP_LOW(A2, 0.0);
        STA2_i     =  STA2;
        A3_i       = `CLIP_LOW(A3, 0.0);
        A4_i       = `CLIP_LOW(A4, 0.0);
        GCO_i      = `CLIP_BOTH(GCO, -10.0, 10.0);
        IGINV_i    = `CLIP_LOW(IGINV, 0.0);
        IGOV_i     = `CLIP_LOW(IGOV, 0.0);
        IGOVD_i    = `CLIP_LOW(IGOVD, 0.0);
        STIG_i     =  STIG;
        GC2_i      = `CLIP_BOTH(GC2, 0.0, 10.0);
        GC3_i      = `CLIP_BOTH(GC3, -10.0, 10.0);
        CHIB_i     = `CLIP_LOW(CHIB, 1.0);
        AGIDL_i    = `CLIP_LOW(AGIDL, 0.0);
        AGIDLD_i   = `CLIP_LOW(AGIDLD, 0.0);
        BGIDL_i    = `CLIP_LOW(BGIDL, 0.0);
        BGIDLD_i   = `CLIP_LOW(BGIDLD, 0.0);
        STBGIDL_i  =  STBGIDL;
        STBGIDLD_i =  STBGIDLD;
        CGIDL_i    =  CGIDL;
        CGIDLD_i   =  CGIDLD;
        COX_i      = `CLIP_LOW(COX, 0.0);
        CGOV_i     = `CLIP_LOW(CGOV, 0.0);
        CGOVD_i    = `CLIP_LOW(CGOVD, 0.0);
        CGBOV_i    = `CLIP_LOW(CGBOV, 0.0);
        CFR_i      = `CLIP_LOW(CFR, 0.0);
        CFRD_i     = `CLIP_LOW(CFRD, 0.0);
        FNT_i      = `CLIP_LOW(FNT, 0.0);
        FNTEXC_i   = `CLIP_LOW(FNTEXC, 0.0);
        NFA_i      = `CLIP_LOW(NFA, 0.0);
        NFB_i      = `CLIP_LOW(NFB, 0.0);
        NFC_i      = `CLIP_LOW(NFC, 0.0);
        EF_i       = `CLIP_LOW(EF,  0.0);
        FACTUO_i   = `CLIP_LOW(FACTUO, 0.0);
        DELVTO_i   =  DELVTO;
        // Parasitic resistances
        RG_i       = `CLIP_LOW(RG, 1.0e-3);
        RBULK_i    = `CLIP_LOW(RBULK, 1.0e-3);
        RJUNS_i    = `CLIP_LOW(RJUNS, 1.0e-3);
        RJUND_i    = `CLIP_LOW(RJUND, 1.0e-3);
        RWELL_i    = `CLIP_LOW(RWELL, 1.0e-3);
`ifdef LocalModel
            MULT_i     = `CLIP_LOW(MULT, 0.0);
`else // LocalModel
            MULT_i     = `CLIP_LOW(MULT * NF_i, 0.0);
`endif // LocalModel
`ifdef NQSmodel
            MUNQS_i    = `CLIP_LOW(MUNQS, 0.0);
`endif // NQSmodel

        // ignore drain-side values in case of symmetric junctions
        if (SWJUNASYM == 0.0) begin
            TOXOVD_i    =  TOXOV_i;
            NOVD_i      =  NOV_i;
            AGIDLD_i    =  AGIDL_i;
            BGIDLD_i    =  BGIDL_i;
            STBGIDLD_i  =  STBGIDL_i;
            CGIDLD_i    =  CGIDL_i;
            IGOVD_i     =  IGOV_i;
            CGOVD_i     =  CGOV_i;
            CFRD_i      =  CFR_i;
        end

`ifdef LocalModel
            // JUNCAP instance parameters (local model)
            jwcorr     = 0.0;
            ABS_i      = 0.0;
            LSS_i      = 0.0;
            LGS_i      = 0.0;
            ABD_i      = 0.0;
            LSD_i      = 0.0;
            LGD_i      = 0.0;
            if (SWJUNCAP == 3.0) begin
                jwcorr = JW_i;
            end
            if ((SWJUNCAP == 2.0) || (SWJUNCAP == 3.0)) begin
                ABS_i = AS_i;
                LSS_i = PS_i - jwcorr;
                LGS_i = JW_i;
                ABD_i = AD_i;
                LSD_i = PD_i - jwcorr;
                LGD_i = JW_i;
            end else begin
                ABS_i = ABSOURCE_i;
                LSS_i = LSSOURCE_i;
                LGS_i = LGSOURCE_i;
                ABD_i = ABDRAIN_i;
                LSD_i = LSDRAIN_i;
                LGD_i = LGDRAIN_i;
            end
`endif // LocalModel
        if ((SWJUNCAP == 1.0) || (SWJUNCAP == 2.0) || (SWJUNCAP == 3.0)) begin
            ABSOURCE_i = `CLIP_LOW(ABS_i, `AB_cliplow);
            LSSOURCE_i = `CLIP_LOW(LSS_i, `LS_cliplow);
            LGSOURCE_i = `CLIP_LOW(LGS_i, `LG_cliplow);
            ABDRAIN_i  = `CLIP_LOW(ABD_i, `AB_cliplow);
            LSDRAIN_i  = `CLIP_LOW(LSD_i, `LS_cliplow);
            LGDRAIN_i  = `CLIP_LOW(LGD_i, `LG_cliplow);
        end else begin
            ABSOURCE_i = 0.0;
            LSSOURCE_i = 0.0;
            LGSOURCE_i = 0.0;
            ABDRAIN_i  = 0.0;
            LSDRAIN_i  = 0.0;
            LGDRAIN_i  = 0.0;
        end

        // Local process parameters
        EPSOX      =  `EPSO * EPSROX_i;
        phit1      =  phit * (1.0 + CT_i * rTn);
        inv_phit1  =  1.0 / phit1;

        VFB_i      =  VFB_i + STVFB_i * dT + DELVTO_i;
        phib       =  Eg + DPHIB_i + 2.0 * phit * ln(NEFF_i * pow(phibFac, -0.75) * 4.0e-26);
        phib       =  `MAX(phib, 5.0E-2);
        CoxPrime   =  EPSOX / TOX_i;
        tox_sq     =  TOX_i * TOX_i;
        G_0        =  sqrt(2.0 * `QELE * NEFF_i * EPSSI * inv_phit) / CoxPrime;

        // Poly-silicon depletion
        kp         =  0.0;
        if (NP_i > 0.0) begin
            arg2max    =  8.0e7 / tox_sq;
            np         = `MAX(NP_i, arg2max);
            np         = `MAX(5.0e24, np);
            kp         =  2.0 * CoxPrime * CoxPrime * phit / (`QELE * np * EPSSI);
        end

        // QM corrections
        qlim2      =  100.0 * phit * phit;
        qq         =  0.0;
        if (QMC_i > 0.0) begin
            qq         =  0.4 * `QMN * QMC_i * pow(CoxPrime, `twoThirds);
            if (CHNL_TYPE==`PMOS) begin
                qq         =  `QMP / `QMN * qq;
            end
            qb0        =  sqrt(phit * G_0 * G_0 * phib);
            dphibq     =  0.75 * qq * pow(qb0, `twoThirds);
            phib       =  phib + dphibq;
            G_0        =  G_0 * (1.0 + 2.0 * `twoThirds * dphibq / qb0);
        end
        sqrt_phib  =  sqrt(phib);
        phix       =  0.95 * phib;
        aphi       =  0.0025 * phib * phib;
        bphi       =  aphi;
        phix2      =  0.5 * sqrt(bphi);
        phix1      =  `MINA(phix - phix2, 0, aphi);

        // Gate overlap
        CoxovPrime   =  EPSOX / TOXOV_i;
        CoxovPrime_d =  EPSOX / TOXOVD_i;
        GOV_s        =  sqrt(2.0 * `QELE * NOV_i * EPSSI * inv_phit) / CoxovPrime;
        GOV_d        =  sqrt(2.0 * `QELE * NOVD_i * EPSSI * inv_phit) / CoxovPrime_d;
        GOV2_s       =  GOV_s * GOV_s;
        GOV2_d       =  GOV_d * GOV_d;
        xi_ov_s      =  1.0 + GOV_s * `invSqrt2;
        xi_ov_d      =  1.0 + GOV_d * `invSqrt2;
        inv_xi_ov_s  =  1.0 / xi_ov_s;
        inv_xi_ov_d  =  1.0 / xi_ov_d;
        x_mrg_ov_s   =  1.0e-5 * xi_ov_s;
        x_mrg_ov_d   =  1.0e-5 * xi_ov_d;

        // Mobility parameters
        tf_bet     =  pow(rTn, STBET_i);
        BET_i      =  FACTUO_i * BETN_i * CoxPrime * tf_bet;
        BETN_i     =  BETN_i * tf_bet;
        THEMU_i    =  THEMU_i * pow(rTn, STTHEMU_i);
        tf_mue     =  pow(rTn, STMUE_i);
        MUE_i      =  MUE_i * tf_mue;
        tf_cs      =  pow(rTn, STCS_i);
        CS_i       =  CS_i * tf_cs;
        tf_xcor    =  pow(rTn, STXCOR_i);
        XCOR_i     =  XCOR_i * tf_xcor;
        E_eff0     =  1.0e-8 * CoxPrime / EPSSI;
        eta_mu     =  0.5 * FETA_i;
        if (CHNL_TYPE == `PMOS) begin
            eta_mu     =  `oneThird * FETA_i;
        end

        // Series resistance
        tf_ther    =  pow(rTn, STRS_i);
        RS_i       =  RS_i * tf_ther;
        THER_i     =  2 * BET_i * RS_i;

        // Velocity saturation
        tf_thesat  =  pow(rTn, STTHESAT_i);
        THESAT_i   =  THESAT_i * tf_thesat;
        Vdsat_lim  =  3.912023005 * phit1;

        inv_AX     =  1.0 / AX_i;
        inv_VP     =  1.0 / VP_i;

        // Impact ionization
        A2_i       =  A2_i * pow(rT, STA2_i);

        // Gate current
        tf_ig      =  pow(rT, STIG_i);
        IGINV_i    =  IGINV_i * tf_ig;
        IGOV_i     =  IGOV_i * tf_ig;
        IGOVD_i    =  IGOVD_i * tf_ig;
        inv_CHIB   =  1.0 / CHIB_i;
        tempM      =  4.0 * `oneThird * sqrt(2 * `QELE * `MELE * CHIB_i) / `HBAR;
        BCH        =  tempM * TOX_i;
        BOV        =  tempM * TOXOV_i;
        BOV_d      =  tempM * TOXOVD_i;
        GCQ        =  0;
        if (GC3_i < 0) begin
            GCQ        =  -0.495 * GC2_i / GC3_i;
        end
        alpha_b    =  0.5 * (phib + Eg);
        Dch        =  GCO_i * phit1;
        Dov        =  GCO_i * phit;

        // GIDL
        AGIDLs     =  AGIDL_i * 4e-18 / (TOXOV_i * TOXOV_i);
        AGIDLDs    =  AGIDLD_i * 4e-18 / (TOXOVD_i * TOXOVD_i);
        tempM      = `MAX(1.0 + STBGIDL_i * dT, 0);
        BGIDL_T    =  BGIDL_i * tempM;
        BGIDLs     =  BGIDL_T * TOXOV_i * 5e8;
        tempM      = `MAX(1.0 + STBGIDLD_i * dT, 0);
        BGIDLD_T   =  BGIDLD_i * tempM;
        BGIDLDs    =  BGIDLD_T * TOXOVD_i * 5e8;

        // Noise
        nt         =  FNT_i * 4 * `KBOL * TKD;
        Cox_over_q =  CoxPrime / `QELE;
        Sfl_prefac =  phit * phit * BET_i / Cox_over_q;
        fac_exc    =  `MELE * 1e9 * FNTEXC_i;

        // Additional internal parameters
        x1         =  1.25;
        inv_xg1_s  =  1.0 / (x1 + GOV_s * 7.324648775608221e-1); // =  1.0/(x1+GOV*sqrt(exp(-x1)+x1-1));
        inv_xg1_d  =  1.0 / (x1 + GOV_d * 7.324648775608221e-1);

        // Conductance of parasitic resistance
        if (RG_i > 0.0) begin
            ggate  =  1.0 / RG_i;
        end else begin
            ggate  =  0.0;
        end
        if (RBULK_i > 0.0) begin
            gbulk  =  1.0 / RBULK_i;
        end else begin
            gbulk  =  0.0;
        end
        if (RJUNS_i > 0.0) begin
            gjuns  =  1.0 / RJUNS_i;
        end else begin
            gjuns  =  0.0;
        end
        if (RJUND_i > 0.0) begin
            gjund  =  1.0 / RJUND_i;
        end else begin
            gjund  =  0.0;
        end
        if (RWELL_i > 0.0) begin
            gwell  =  1.0 / RWELL_i;
        end else begin
            gwell  =  0.0;
        end


        // JUNCAP2

        // Initialization of JUNCAP (global) variables; required for some verilog-A compilers
        vbimin_s   = 0.0;
        vbimin_d   = 0.0;
        vfmin_s    = 0.0;
        vfmin_d    = 0.0;
        vch_s      = 0.0;
        vch_d      = 0.0;
        vbbtlim_s  = 0.0;
        vbbtlim_d  = 0.0;
        VMAX_s     = 0.0;
        VMAX_d     = 0.0;
        exp_VMAX_over_phitd_s = 0.0;
        exp_VMAX_over_phitd_d = 0.0;
        ISATFOR1_s = 0.0;
        ISATFOR1_d = 0.0;
        MFOR1_s    = 1.0;
        MFOR1_d    = 1.0;
        ISATFOR2_s = 0.0;
        ISATFOR2_d = 0.0;
        MFOR2_s    = 1.0;
        MFOR2_d    = 1.0;
        ISATREV_s  = 0.0;
        ISATREV_d  = 0.0;
        MREV_s     = 1.0;
        MREV_d     = 1.0;
        m0flag_s   = 0.0;
        m0flag_d   = 0.0;
        xhighf1_s  = 0.0;
        xhighf1_d  = 0.0;
        expxhf1_s  = 0.0;
        expxhf1_d  = 0.0;
        xhighf2_s  = 0.0;
        xhighf2_d  = 0.0;
        expxhf2_s  = 0.0;
        expxhf2_d  = 0.0;
        xhighr_s   = 0.0;
        xhighr_d   = 0.0;
        expxhr_s   = 0.0;
        expxhr_d   = 0.0;
        zflagbot_s = 1.0;
        zflagbot_d = 1.0;
        zflagsti_s = 1.0;
        zflagsti_d = 1.0;
        zflaggat_s = 1.0;
        zflaggat_d = 1.0;

        m0_rev     = 0.0;
        mcor_rev   = 0.0;
        I1_cor     = 0.0;
        I2_cor     = 0.0;
        I3_cor     = 0.0;
        I4_cor     = 0.0;
        I5_cor     = 0.0;
        tt0        = 0.0;
        tt1        = 0.0;
        tt2        = 0.0;
        zfrac      = 0.0;
        alphaje    = 0.0;

        if (SWJUNCAP > 0) begin
            `JuncapInitInstance(ABSOURCE_i, LSSOURCE_i, LGSOURCE_i, idsatbot,   idsatsti,   idsatgat,   vbibot,   vbisti,   vbigat,   PBOT_i,  PSTI_i,  PGAT_i,  VBIRBOT_i,  VBIRSTI_i,  VBIRGAT_i,  VMAX_s, exp_VMAX_over_phitd_s, vbimin_s, vch_s, vfmin_s, vbbtlim_s)
            `JuncapInitInstance(ABDRAIN_i,  LSDRAIN_i,  LGDRAIN_i,  idsatbot_d, idsatsti_d, idsatgat_d, vbibot_d, vbisti_d, vbigat_d, PBOTD_i, PSTID_i, PGATD_i, VBIRBOTD_i, VBIRSTID_i, VBIRGATD_i, VMAX_d, exp_VMAX_over_phitd_d, vbimin_d, vch_d, vfmin_d, vbbtlim_d)

            if (SWJUNEXP_i == 1.0) begin : JUNCAPexpressInit
                // The variables in the macro below are (re-)declared LOCALLY, to keep
                // them separated from their globally declared counterparts. This trick
                // allows one to use the "juncapcommon" macro both in the JUNCAP-express
                // initialization and in the full-JUNCAP evaluation, while in the former
                // the verilog-A compiler can still consider the variables as
                // voltage-INdependent. This is essential to avoid recomputation of the
                // JUNCAP-express initialization at each bias-step.
//                `LocalGlobalVars
                // results computed here are not used elsewhere
//                real ijunbot, ijunsti, ijungat, qjunbot, qjunsti, qjungat;

                // Initialization of (local) variables; required for some verilog-A compilers
                ysq = 0.0;
                terfc = 0.0;
                erfcpos = 0.0;
                h1 = 0.0;
                h2 = 0.0;
                h2d = 0.0;
                h3 = 0.0;
                h4 = 0.0;
                h5 = 0.0;
                idmult = 0.0;
                vj = 0.0;
                z = 0.0;
                zinv = 0.0;
                two_psistar = 0.0;
                vjlim = 0.0;
                vjsrh = 0.0;
                vbbt = 0.0;
                vav = 0.0;
                tmp = 0.0;
                id = 0.0;
                isrh = 0.0;
                vbi_minus_vjsrh = 0.0;
                wsrhstep = 0.0;
                dwsrh = 0.0;
                wsrh = 0.0;
                wdep = 0.0;
                asrh = 0.0;
                itat = 0.0;
                btat = 0.0;
                twoatatoverthreebtat = 0.0;
                umaxbeforelimiting = 0.0;
                umax = 0.0;
                sqrtumax = 0.0;
                umaxpoweronepointfive = 0.0;
                wgamma = 0.0;
                wtat = 0.0;
                ktat = 0.0;
                ltat = 0.0;
                mtat = 0.0;
                xerfc = 0.0;
                erfctimesexpmtat = 0.0;
                gammamax = 0.0;
                ibbt = 0.0;
                Fmaxr = 0.0;
                fbreakdown = 0.0;
                qjunbot = 0.0;
                qjunsti = 0.0;
                qjungat = 0.0;

                // Computation of JUNCAP-express internal parameters
                `JuncapExpressInit1(ABSOURCE_i, LSSOURCE_i, LGSOURCE_i, VJUNREF_i, qprefbot, qpref2bot, vbiinvbot, one_minus_PBOT, idsatbot, CSRHBOT_i, CTATBOT_i, vbibot, wdepnulrbot, VBIRBOTinv, PBOT_i, ftdbot, btatpartbot, atatbot, one_over_one_minus_PBOT, CBBTBOT_i, VBIRBOT_i, wdepnulrinvbot, fbbtbot, VBRBOT_i, VBRinvbot, PBRBOT_i, fstopbot, slopebot, qprefsti, qpref2sti, vbiinvsti, one_minus_PSTI, idsatsti, CSRHSTI_i, CTATSTI_i, vbisti, wdepnulrsti, VBIRSTIinv, PSTI_i, ftdsti, btatpartsti, atatsti, one_over_one_minus_PSTI, CBBTSTI_i, VBIRSTI_i, wdepnulrinvsti, fbbtsti, VBRSTI_i, VBRinvsti, PBRSTI_i, fstopsti, slopesti, qprefgat, qpref2gat, vbiinvgat, one_minus_PGAT, idsatgat, CSRHGAT_i, CTATGAT_i, vbigat, wdepnulrgat, VBIRGATinv, PGAT_i, ftdgat, btatpartgat, atatgat, one_over_one_minus_PGAT, CBBTGAT_i, VBIRGAT_i, wdepnulrinvgat, fbbtgat, VBRGAT_i, VBRinvgat, PBRGAT_i, fstopgat, slopegat, VMAX_s, exp_VMAX_over_phitd_s, vbimin_s, vch_s, vfmin_s, vbbtlim_s)
                `JuncapExpressInit2(ABSOURCE_i, LSSOURCE_i, LGSOURCE_i, qprefbot, qpref2bot, vbiinvbot, one_minus_PBOT, idsatbot, CSRHBOT_i, CTATBOT_i, vbibot, wdepnulrbot, VBIRBOTinv, PBOT_i, ftdbot, btatpartbot, atatbot, one_over_one_minus_PBOT, CBBTBOT_i, VBIRBOT_i, wdepnulrinvbot, fbbtbot, VBRBOT_i, VBRinvbot, PBRBOT_i, fstopbot, slopebot, qprefsti, qpref2sti, vbiinvsti, one_minus_PSTI, idsatsti, CSRHSTI_i, CTATSTI_i, vbisti, wdepnulrsti, VBIRSTIinv, PSTI_i, ftdsti, btatpartsti, atatsti, one_over_one_minus_PSTI, CBBTSTI_i, VBIRSTI_i, wdepnulrinvsti, fbbtsti, VBRSTI_i, VBRinvsti, PBRSTI_i, fstopsti, slopesti, qprefgat, qpref2gat, vbiinvgat, one_minus_PGAT, idsatgat, CSRHGAT_i, CTATGAT_i, vbigat, wdepnulrgat, VBIRGATinv, PGAT_i, ftdgat, btatpartgat, atatgat, one_over_one_minus_PGAT, CBBTGAT_i, VBIRGAT_i, wdepnulrinvgat, fbbtgat, VBRGAT_i, VBRinvgat, PBRGAT_i, fstopgat, slopegat, VMAX_s, exp_VMAX_over_phitd_s, vbimin_s, vch_s, vfmin_s, vbbtlim_s)
                `JuncapExpressInit3(ABSOURCE_i, LSSOURCE_i, LGSOURCE_i, idsatbot, idsatsti, idsatgat, ISATFOR1_s, MFOR1_s, ISATFOR2_s, MFOR2_s, ISATREV_s, MREV_s, m0flag_s)
                `JuncapExpressInit4(ABSOURCE_i, LSSOURCE_i, LGSOURCE_i, FJUNQ_i, cjobot, cjosti, cjogat, zflagbot_s, zflagsti_s, zflaggat_s)
                `JuncapExpressInit5(ABSOURCE_i, LSSOURCE_i, LGSOURCE_i, ISATFOR1_s, ISATFOR2_s, ISATREV_s, xhighf1_s, expxhf1_s, xhighf2_s, expxhf2_s, xhighr_s, expxhr_s)

                `JuncapExpressInit1(ABDRAIN_i, LSDRAIN_i, LGDRAIN_i, VJUNREFD_i, qprefbot_d, qpref2bot_d, vbiinvbot_d, one_minus_PBOT_d, idsatbot_d, CSRHBOTD_i, CTATBOTD_i, vbibot_d, wdepnulrbot_d, VBIRBOTinv_d, PBOTD_i, ftdbot_d, btatpartbot_d, atatbot_d, one_over_one_minus_PBOT_d, CBBTBOTD_i, VBIRBOTD_i, wdepnulrinvbot_d, fbbtbot_d, VBRBOTD_i, VBRinvbot_d, PBRBOTD_i, fstopbot_d, slopebot_d, qprefsti_d, qpref2sti_d, vbiinvsti_d, one_minus_PSTI_d, idsatsti_d, CSRHSTID_i, CTATSTID_i, vbisti_d, wdepnulrsti_d, VBIRSTIinv_d, PSTID_i, ftdsti_d, btatpartsti_d, atatsti_d, one_over_one_minus_PSTI_d, CBBTSTID_i, VBIRSTID_i, wdepnulrinvsti_d, fbbtsti_d, VBRSTID_i, VBRinvsti_d, PBRSTID_i, fstopsti_d, slopesti_d, qprefgat_d, qpref2gat_d, vbiinvgat_d, one_minus_PGAT_d, idsatgat_d, CSRHGATD_i, CTATGATD_i, vbigat_d, wdepnulrgat_d, VBIRGATinv_d, PGATD_i, ftdgat_d, btatpartgat_d, atatgat_d, one_over_one_minus_PGAT_d, CBBTGATD_i, VBIRGATD_i, wdepnulrinvgat_d, fbbtgat_d, VBRGATD_i, VBRinvgat_d, PBRGATD_i, fstopgat_d, slopegat_d, VMAX_d, exp_VMAX_over_phitd_d, vbimin_d, vch_d, vfmin_d, vbbtlim_d)
                `JuncapExpressInit2(ABDRAIN_i, LSDRAIN_i, LGDRAIN_i, qprefbot_d, qpref2bot_d, vbiinvbot_d, one_minus_PBOT_d, idsatbot_d, CSRHBOTD_i, CTATBOTD_i, vbibot_d, wdepnulrbot_d, VBIRBOTinv_d, PBOTD_i, ftdbot_d, btatpartbot_d, atatbot_d, one_over_one_minus_PBOT_d, CBBTBOTD_i, VBIRBOTD_i, wdepnulrinvbot_d, fbbtbot_d, VBRBOTD_i, VBRinvbot_d, PBRBOTD_i, fstopbot_d, slopebot_d, qprefsti_d, qpref2sti_d, vbiinvsti_d, one_minus_PSTI_d, idsatsti_d, CSRHSTID_i, CTATSTID_i, vbisti_d, wdepnulrsti_d, VBIRSTIinv_d, PSTID_i, ftdsti_d, btatpartsti_d, atatsti_d, one_over_one_minus_PSTI_d, CBBTSTID_i, VBIRSTID_i, wdepnulrinvsti_d, fbbtsti_d, VBRSTID_i, VBRinvsti_d, PBRSTID_i, fstopsti_d, slopesti_d, qprefgat_d, qpref2gat_d, vbiinvgat_d, one_minus_PGAT_d, idsatgat_d, CSRHGATD_i, CTATGATD_i, vbigat_d, wdepnulrgat_d, VBIRGATinv_d, PGATD_i, ftdgat_d, btatpartgat_d, atatgat_d, one_over_one_minus_PGAT_d, CBBTGATD_i, VBIRGATD_i, wdepnulrinvgat_d, fbbtgat_d, VBRGATD_i, VBRinvgat_d, PBRGATD_i, fstopgat_d, slopegat_d, VMAX_d, exp_VMAX_over_phitd_d, vbimin_d, vch_d, vfmin_d, vbbtlim_d)
                `JuncapExpressInit3(ABDRAIN_i, LSDRAIN_i, LGDRAIN_i, idsatbot_d, idsatsti_d, idsatgat_d, ISATFOR1_d, MFOR1_d, ISATFOR2_d, MFOR2_d, ISATREV_d, MREV_d, m0flag_d)
                `JuncapExpressInit4(ABDRAIN_i, LSDRAIN_i, LGDRAIN_i, FJUNQD_i, cjobot_d, cjosti_d, cjogat_d, zflagbot_d, zflagsti_d, zflaggat_d)
                `JuncapExpressInit5(ABDRAIN_i, LSDRAIN_i, LGDRAIN_i, ISATFOR1_d, ISATFOR2_d, ISATREV_d, xhighf1_d, expxhf1_d, xhighf2_d, expxhf2_d, xhighr_d, expxhr_d)
            end // JUNCAPexpressInit

        end


    end // initial_instance

    /////////////////////////////////////////////////////////////////////////////
    //
    //      DC bias dependent quantities (calculations for current contribs)
    //
    /////////////////////////////////////////////////////////////////////////////

    begin : evaluateblock
//        real alpha, eta_p, xitsb, FdL, Gmob, Gmob_dL, Gvsat, Gvsatinv, dps, qim, qim1;
//        real H, Iimpact, mavl, Ids, Vdse, Igsov, Igdov, Igcs, Igcd, xg, sigVds;
//        real COX_qm;
//        real ijun_s, ijunbot_s, ijunsti_s, ijungat_s;
//        real ijun_d, ijunbot_d, ijunsti_d, ijungat_d;
//        real qjun_s, qjunbot_s, qjunsti_s, qjungat_s;
//        real qjun_d, qjunbot_d, qjunsti_d, qjungat_d;
//        real jnoise_s, jnoise_d;

        begin : evaluateStatic

            // Initialisation of some variables
            SP_S_x1    =  0.0;
            x_s        =  0.0;
            sqm        =  0.0;
            alpha      =  0.0;
            eta_p      =  1.0;
            xitsb      =  0.0;
            rhob       =  0.0;
            GdL        =  1.0;
            FdL        =  1.0;
            Gmob       =  1.0;
            Gmob_dL    =  1.0;
            Udse       =  0.0;
            QCLM       =  0.0;
            thesat1    =  0.0;
            Gvsat      =  1.0;
            Gvsatinv   =  1.0;
            xgm        =  0.0;
            dps        =  0.0;
            qim        =  0.0;
            qim1       =  0.0;
            H          =  1.0;
            xs_ov      =  0.0;
            xd_ov      =  0.0;
            Vovs       =  0.0;
            Vovd       =  0.0;
            Iimpact    =  0.0;
            mavl       =  0.0;
`ifdef NQSmodel
                // Initialization of variables for NQS model
                pd         =  1.0;
                ym         =  0.0;
`endif // NQSmodel

            if (CHNL_TYPE == `NMOS) begin
                Vgs        =  V(GP, S);
                Vds        =  V(D, S);
                Vsb        =  V(S, BP);
                Vjun_s     = -V(S, BS);
                Vjun_d     = -V(D, BD);
            end else begin
                Vgs        = -V(GP, S);
                Vds        = -V(D, S);
                Vsb        = -V(S, BP);
                Vjun_s     =  V(S, BS);
                Vjun_d     =  V(D, BD);
            end

            // Voltages NOT subject to S/D-interchange
            VgsPrime   = Vgs;
            VsbPrime   = Vsb;
            VdbPrime   = Vds + Vsb;
            VgdPrime   = Vgs - Vds;
            xgs_ov     = -VgsPrime * inv_phit;
            xgd_ov     = -VgdPrime * inv_phit;

            // Source-drain interchange
            sigVds     =  1.0;
            if (Vds < 0.0) begin
                sigVds     = -1.0;
                Vgs        =  Vgs - Vds;
                Vsb        =  Vsb + Vds;
                Vds        = -Vds;
            end

            Vdb        =  Vds + Vsb;

            // 4.2.1 Conditioning of terminal voltages
            temp       =  `MINA(Vdb, Vsb, bphi) + phix;
            Vsbstar    =  Vsb - `MINA(temp, 0, aphi) + phix1;
            Vgbstar    =  Vgs + Vsbstar;
            Vgb1       =  Vgbstar - VFB_i;
            Vdsx       =  sqrt(Vds * Vds + 0.01) - 0.1;
            Vsbx       =  Vsbstar + 0.5 * (Vds - Vdsx);
            delVg      =  CF_i * (Vdsx * (1 + CFB_i * Vsbx)); // DIBL
            Vgb1       =  Vgb1 + delVg;
            xg         =  Vgb1 * inv_phit1;


            // 4.2.2 Bias dependent body factor
            if (DNSUB_i > 0.0) begin
                Dnsub       =  DNSUB_i * `MAXA(0, Vgs + Vsb - VNSUB_i, NSLP_i);
                Gf          =  G_0 * sqrt(1.0 + Dnsub);
            end else begin
                Gf          =  G_0;
            end

            // 4.2.3 Surface potential at source side
            Gf2        =  Gf * Gf;
            inv_Gf2    =  1.0 / Gf2;
            xi         =  1.0 + Gf * `invSqrt2;
            inv_xi     =  1.0 /  xi;
            Ux         =  Vsbstar * inv_phit1;
            xn_s       =  phib * inv_phit1 + Ux;
            if (xn_s < `se)
                delta_ns   =  exp(-xn_s);
            else
                delta_ns   = `ke / `P3(xn_s - `se);
            margin     =  1e-5 * xi;

            `sp_s(x_s, xg, xn_s, delta_ns)
            x_d        =  x_s;
            x_m        =  x_s;
            x_ds       =  0.0;

            //
            // Core PSP current calculation
            //
            if (xg <= 0.0) begin
                qis        =  0.0;
                Ids        =  0.0;
                xgm        =  xg - x_s;
                Voxm       =  xgm * phit1;
                qeff       =  Voxm;
                Vdsat      =  Vdsat_lim;
                Vdse       =  Vds;
            end else begin // (xg > 0)
                delta_1s   =  0.0;
                temp       =  1.0 / (2.0 + x_s * x_s);
                xi0s       =  x_s * x_s * temp;
                xi1s       =  4.0 * (x_s * temp * temp);
                xi2s       =  (8.0 * temp - 12.0 * xi0s) * temp * temp;
                if (x_s < `se05) begin
                    delta_1s   =  exp(x_s);
                    Es         =  1.0 / delta_1s;
                    delta_1s   =  delta_ns * delta_1s;
                end else if (x_s > (xn_s - `se05)) begin
                    delta_1s   =  exp(x_s - xn_s);
                    Es         =  delta_ns / delta_1s;
                end else begin
                    delta_1s   = `ke05 / `P3(xn_s - x_s - `se05);
                    Es         = `ke05 / `P3(x_s - `se05);
                end
                Ds         =  delta_1s - delta_ns * (x_s + 1.0 + xi0s);
                if (x_s < 1.0e-5) begin
                    Ps         =  0.5 * (x_s * x_s * (1.0 - `oneThird * (x_s * (1.0 - 0.25 * x_s))));
                    Ds         =  `oneSixth * (delta_ns * x_s * x_s * x_s * (1.0 + 1.75 * x_s));
                    temp       =  sqrt(1.0 - `oneThird * (x_s * (1.0 - 0.25 * x_s)));
                    sqm        =  `invSqrt2 * (x_s * temp);
                    alpha      =  1.0 + Gf * `invSqrt2 * (1.0 - 0.5 * x_s + `oneSixth * (x_s * x_s)) / temp;
                end else begin
                    Ps         =  x_s - 1.0 + Es;
                    sqm        =  sqrt(Ps);
                    alpha      =  1.0 + 0.5 * (Gf * (1.0 - Es) / sqm);
                end
                Em     =  Es;
                Ed     =  Em;
                Dm     =  Ds;
                Dd     =  Dm;

                // 4.2.4 Drain saturation voltage
                Rxcor      =  (1.0 + 0.2 * XCOR_i * Vsbx) / (1.0 + XCOR_i * Vsbx);
                if ( Ds > `ke05) begin
                    xgs        =  Gf * sqrt(Ps + Ds);
                    qis        =  Gf2 * Ds * phit1 / (xgs + Gf * sqm);
                    qbs        =  sqm * Gf * phit1;
                    if (RSB_i < 0) begin
                        rhob       = 1.0 / (1.0 - RSB_i * Vsbx);
                    end else begin
                        rhob       = 1.0 + RSB_i * Vsbx;
                    end
                    if (RSG_i < 0) begin
                        temp       = 1.0 - RSG_i * qis;
                    end else begin
                        temp       = 1.0 / (1.0 + RSG_i * qis);
                    end
                    GR         =  THER_i * (rhob * temp * qis);
                    Eeffm      =  E_eff0 * (qbs + eta_mu * qis);
                    Mutmp      =  pow(Eeffm * MUE_i, THEMU_i) + CS_i * (Ps / (Ps + Ds + 1.0e-14));
                    Gmob       =  (1.0 + Mutmp + GR) * Rxcor;
                    if (THESATB_i < 0) begin
                        xitsb      = 1.0 / (1.0 - THESATB_i * Vsbx);
                    end else begin
                        xitsb      = 1.0 + THESATB_i * Vsbx;
                    end
                    temp2      =  qis * xitsb;
                    wsat       =  100.0 * (temp2 / (100.0 + temp2));
                    if (THESATG_i < 0) begin
                        temp       = 1 / (1 - THESATG_i * wsat);
                    end else begin
                        temp       = 1 + THESATG_i * wsat;
                    end
                    thesat1    =  THESAT_i * (temp / Gmob);
                    phi_inf    =  qis / alpha + phit1;
                    ysat       =  thesat1 * phi_inf * `invSqrt2;
                    if (CHNL_TYPE==`PMOS) begin
                        ysat       =  ysat / sqrt(1.0 + ysat);
                    end
                    za         =  2.0 / (1.0 + sqrt(1.0 + 4.0 * ysat));
                    temp1      =  za * ysat;
                    Phi_0      =  phi_inf * za * (1.0 + 0.86 * (temp1 * (1.0 - temp1 * za) / (1.0 + 4.0 * (temp1 * temp1 * za))));
                    asat       =  xgs + 0.5 * Gf2;
                    Phi_2      =  0.98 * (Gf2 * Ds * phit1 / (asat + sqrt(asat * asat - Gf2 * Ds * 0.98)));
                    Phi_0_2    =  Phi_0 + Phi_2;
                    Phi0_Phi2  =  2.0 * (Phi_0 * Phi_2);
                    Phi_sat    =  Phi0_Phi2 / (Phi_0_2 + sqrt(Phi_0_2 * Phi_0_2 - 1.98 * Phi0_Phi2));
                    Vdsat      =  Phi_sat - phit1 * ln(1.0 + Phi_sat * (Phi_sat - 2.0 * asat * phit1) * inv_Gf2 / (phit1 * phit1 * Ds));
                end else begin
                    Vdsat      =  Vdsat_lim;
                end
                temp      =  pow(Vds / Vdsat, AX_i);
                Vdse      =  Vds * pow(1.0 + temp, -inv_AX);

                // 4.2.5 Surface potential at drain side
                Udse       =  Vdse * inv_phit1;
                xn_d       =  xn_s + Udse;
                if (Udse < `se) begin
                    k_ds       =  exp(-Udse);
                end else begin
                    k_ds       = `ke / `P3(Udse - `se);
                end
                delta_nd   =  delta_ns * k_ds;

                `sp_s_d(x_d, xg, xn_d, delta_nd)
                x_ds       =  x_d - x_s;

                //
                //  Approximations for extremely small x_ds: capacitance calulation
                //
                if (x_ds < 1.0E-10) begin
                    pC          =  2.0 * (xg - x_s) + Gf2 * (1.0 - Es + delta_1s * k_ds - delta_nd * (1.0 + xi1s));
                    qC          =  Gf2 * (1.0 - k_ds) * Ds;
                    temp        =  2.0 - Gf2 * (Es + delta_1s * k_ds - delta_nd * xi2s);
                    temp        =  pC * pC - 2.0 * (temp * qC);
                    x_ds        =  2.0 * (qC / (pC + sqrt(temp)));
                    x_d         =  x_s + x_ds;
                end
                dps         =  x_ds * phit1; // deltaPsi

                xi0d        =  x_d * x_d / (2.0 + x_d * x_d);
                if (x_d < `se05) begin
                    Ed         =  exp(-x_d);
                    if (x_d < 1.0e-5) begin
                        Dd         =  `oneSixth * delta_nd * x_d * x_d * x_d * (1.0 + 1.75 * x_d);
                    end else begin
                        Dd         =  delta_nd * (1.0 / Ed - x_d - 1.0 - xi0d);
                    end
                end else begin
                    if (x_d > (xn_d - `se05)) begin
                        temp       =  exp(x_d - xn_d);
                        Ed         =  delta_nd / temp;
                        Dd         =  temp - delta_nd * (x_d + 1.0 + xi0d);
                    end else begin
                        Ed         = `ke05 / `P3(x_d - `se05);
                        temp       = `ke05 / `P3(xn_d - x_d - `se05);
                        Dd         =  temp - delta_nd * (x_d + 1.0 + xi0d);
                    end
                end

                // 4.2.6 Mid-point surface potential
                x_m        =  0.5 * (x_s + x_d);
                Em         =  0.0;
                temp = Ed * Es;
                if (temp > 0.0) begin
                    Em         =  sqrt(temp);
                end
                D_bar      =  0.5 * (Ds + Dd);
                Dm         =  D_bar + 0.125 * (x_ds * x_ds * (Em - 2.0 * inv_Gf2));

                if (x_m < 1.0e-5) begin
                    Pm         =  0.5 * (x_m * x_m * (1.0 - `oneThird * (x_m * (1.0 - 0.25 * x_m))));
                    xgm        =  Gf * sqrt(Dm + Pm);

                    // 4.2.7 Polysilicon depletion
                    if (kp > 0.0) begin
                        eta_p       =  1.0 / sqrt(1.0 + kp * xgm);
                    end // (kp > 0.0)
                    temp       =  sqrt(1.0 - `oneThird * (x_m * (1.0 - 0.25 * x_m)));
                    sqm        =  `invSqrt2 * (x_m * temp);
                    alpha      =  eta_p + `invSqrt2 * (Gf * (1.0 - 0.5 * x_m + `oneSixth * (x_m * x_m)) / temp);
                end else begin
                    Pm         =  x_m - 1.0 + Em;
                    xgm        =  Gf * sqrt(Dm + Pm);

                    // 4.2.7 Polysilicon depletion
                    if (kp > 0.0) begin
                        d0         =  1.0 - Em + 2.0 * (xgm * inv_Gf2);
                        eta_p      =  1.0 / sqrt(1.0 + kp * xgm);
                        temp       =  eta_p / (eta_p + 1.0);
                        x_pm       =  kp * (temp * temp * Gf2 * Dm);
                        p_pd       =  2.0 * (xgm - x_pm) + Gf2 * (1.0 - Em + Dm);
                        q_pd       =  x_pm * (x_pm - 2.0 * xgm);
                        xi_pd      =  1.0 - 0.5 * (Gf2 * (Em + Dm));
                        u_pd       =  q_pd * p_pd / (p_pd * p_pd - xi_pd * q_pd);
                        x_m        =  x_m + u_pd;
                        km         =  exp(u_pd);
                        Em         =  Em / km;
                        Dm         =  Dm * km;
                        Pm         =  x_m - 1.0 + Em;
                        xgm        =  Gf * sqrt(Dm + Pm);
                        help       =  1.0 - Em + 2.0 * (xgm * eta_p * inv_Gf2);
                        x_ds       =  x_ds * km * (d0 + D_bar) / (help + km * D_bar);
                        dps        =  x_ds * phit1;
                    end // (kp > 0.0)
                    sqm        =  sqrt(Pm);
                    alpha      =  eta_p + 0.5 * (Gf * (1.0 - Em) / sqm);
                end

                // 4.2.8 Potential midpoint inversion charge
                qim        =  phit1 * (Gf2 * Dm / (xgm + Gf * sqm));

                // 4.2.8 Potential midpoint inversion charge (continued)
                qim1       =  qim + phit1 * alpha;
                qim1_1     =  1.0 / qim1;
                qbm        =  sqm * Gf * phit1;
                // Series resistance
                if (RSG_i < 0) begin
                    temp       = 1.0 - RSG_i * qim;
                end else begin
                    temp       = 1.0 / (1.0 + RSG_i * qim);
                end
                GR         =  THER_i * (rhob * temp * qim);
                // Mobility reduction
                qeff       =  qbm + eta_mu * qim;
                Eeffm      =  E_eff0 * qeff;
                Mutmp      =  pow(Eeffm * MUE_i, THEMU_i) + CS_i * (Pm / (Pm + Dm + 1.0e-14));
                Gmob       =  (1.0 + Mutmp + GR) * Rxcor;

                // 4.2.9 Drain-source channel current
                // Channel length modulation
                r1         =  qim * qim1_1;
                r2         =  phit1 * (alpha * qim1_1);
                temp       =  ln((1.0 + (Vds - dps) * inv_VP) / (1.0 + (Vdse - dps) * inv_VP));
                temp1      =  ln(1.0 + Vdsx * inv_VP);
                dL         =  ALP_i * temp;
                GdL        =  1.0 / (1.0 + dL + dL * dL);
                dL1        =  dL + ALP1_i * (qim1_1 * r1 * temp) + ALP2_i * (qbm * r2 * r2 * temp1);
                FdL        =  (1.0 + dL1 + dL1 * dL1) * GdL;
                // Velocity saturation
                temp2      =  qim * xitsb;
                wsat       =  100.0 * (temp2 / (100.0 + temp2));
                Gmob_dL    =  Gmob * GdL;
                if (THESATG_i < 0) begin
                    temp       = 1 / (1 - THESATG_i * wsat);
                end else begin
                    temp       = 1 + THESATG_i * wsat;
                end
                thesat1    =  THESAT_i * (temp / Gmob_dL);
                zsat       =  thesat1 * thesat1 * dps * dps;
                if (CHNL_TYPE == `PMOS) begin
                    zsat       =  zsat / (1.0 + thesat1 * dps);
                end
                Gvsat      =  0.5 * (Gmob_dL * (1.0 + sqrt(1.0 + 2.0 * zsat)));
                Gvsatinv   =  1.0 / Gvsat;
                // Drain-source current
                Ids        =  BET_i * (FdL * qim1 * dps * Gvsatinv);

                // 4.2.10 Variables for calculation of intrinsic charges and gate current
                Voxm       =  xgm * phit1;
                temp       =  Gmob_dL * Gvsatinv;
                alpha1     =  alpha * (1.0 + 0.5 * (zsat * temp * temp));
                H          =  temp * qim1 / alpha1;

                // 4.2.11 Impact-Ionization
                if (SWIMPACT != 0) begin
                    delVsat       =  Vds - A3_i * dps;
                    if (delVsat > 0) begin
                        temp2        =  A2_i * ((1.0 + A4_i * (sqrt(phib + Vsbstar) - sqrt_phib)) / (delVsat + 1e-30));
                        `expl(-temp2, temp)
                        mavl         =  A1_i * (delVsat * temp);
                        Iimpact      =  Ids * mavl;
                    end
                end
            end // (xg > 0)

            // 4.2.12 Surface potential in gate overlap regions
            if (((SWIGATE != 0) && ((IGOV_i > 0) || (IGOVD_i > 0))) || ((SWGIDL != 0) && ((AGIDL_i > 0) || (AGIDLD_i > 0))) || (CGOV_i > 0) || (CGOVD_i > 0)) begin
                `ChangeToSource
                `sp_ov(xs_ov, xgs_ov)
                `ChangeToDrain
                `sp_ov(xd_ov, xgd_ov)
                Vovs        = -phit * (xgs_ov + xs_ov);
                Vovd        = -phit * (xgd_ov + xd_ov);
            end

            // 4.2.13 Gate current
            Igsov      =  0.0;
            Igdov      =  0.0;
            Igc        =  0.0;
            Igb        =  0.0;
            Igcs       =  0.0;
            Igcd       =  0.0;
            if (SWIGATE != 0) begin
                if (IGOV_i > 0) begin

                    // Gate-source overlap component of gate current
                    arg2mina   =  Vovs + Dov;
                    psi_t      = `MINA(0.0, arg2mina, 0.01);
                    zg         =  sqrt(Vovs * Vovs + 1.0e-6) * inv_CHIB;
                    if (GC3_i < 0) begin
                        zg         = `MINA(zg, GCQ, 1.0e-6);
                    end
                    arg1       =  (3.0 + xs_ov + psi_t * inv_phit);
                    `expl(arg1, Dsi)
                    arg1       =  -VgsPrime * inv_phit;
                    `expl(arg1, temp)
                    Dgate      =  Dsi * temp;
                    temp       =  BOV * (-1.5 + zg * (GC2_i + GC3_i * zg));
                    if (temp > 0) begin
                        TP         = `P3(temp);
                    end else begin
                        `expl_low(temp, TP)
                    end
                    Igsov      =  IGOV_i * (TP * ln((1.0 + Dsi) / (1.0 + Dgate)));
                end

                if (IGOVD_i > 0) begin

                    // Gate-drain overlap component of gate current
                    arg2mina   =  Vovd + Dov;
                    psi_t      = `MINA(0.0, arg2mina, 0.01);
                    zg         =  sqrt(Vovd * Vovd + 1.0e-6) * inv_CHIB;
                    if (GC3_i < 0) begin
                        zg         = `MINA(zg, GCQ, 1.0e-6);
                    end
                    arg1       =  (3.0 + xd_ov + psi_t * inv_phit);
                    `expl(arg1, Dsi)
                    arg1       =  -VgdPrime * inv_phit;
                    `expl(arg1, temp)
                    Dgate      =  Dsi * temp;
                    temp       =  BOV_d * (-1.5 + zg * (GC2_i + GC3_i * zg));
                    if (temp > 0) begin
                        TP         = `P3(temp);
                    end else begin
                        `expl_low(temp, TP)
                    end
                    Igdov      =  IGOVD_i * (TP * ln((1.0 + Dsi) / (1.0 + Dgate)));
                end

                // Gate-channel component of gate current
                if (IGINV_i > 0) begin
                    if (xg <= 0.0) begin
                        temp       =  pow(Vds / Vdsat_lim, AX_i);
                        Udse       =  Vds * pow(1.0 + temp, -inv_AX) * inv_phit1;
                    end
                    `expl_low(x_ds-Udse, temp)
                    Vm         =  Vsbstar + phit1 * (0.5 * x_ds - ln(0.5 * (1.0 + temp)));

                    arg2mina   =  Voxm + Dch;
                    psi_t      = `MINA(0.0, arg2mina, 0.01);
                    zg         =  sqrt(Voxm * Voxm + 1.0e-6) * inv_CHIB;
                    if (GC3_i < 0) begin
                        zg     = `MINA(zg, GCQ, 1.0e-06);
                    end
                    arg1       =  (x_m + (psi_t - alpha_b - Vm) * inv_phit1);
                    `expl(arg1,Dsi)
                    arg1       = -(Vgs + Vsbstar - Vm) * inv_phit1;
                    `expl(arg1,temp)
                    Dgate      =  Dsi * temp;
                    temp       = BCH * (-1.5 + zg * (GC2_i + GC3_i * zg));
                    if (temp > 0) begin
                        TP         = `P3(temp);
                    end else begin
                        `expl_low(temp, TP)
                    end
                    Igc0       =  IGINV_i * (TP * ln((1.0 + Dsi) / (1.0 + Dgate)));

                    // Source/drain partitioning of gate-channel current
                    if ((xg <= 0) || ((GC2_i == 0) && (GC3_i == 0))) begin
                        igc        =  1.0;
                        igcd_h     =  0.5;
                    end else begin
                        temp       =  GC2_i + 2.0 * GC3_i * zg;
                        u0         =  CHIB_i / (temp * BCH);
                        x          =  0.5 * (dps / u0);
                        u0_div_H   =  u0 / H;
                        Bg         =  u0_div_H * (1.0 - u0_div_H) * 0.5;
                        Ag         =  0.5 - 3.0 * Bg;
                        if (x < 1.0e-3) begin
                            xsq        =  x * x;
                            igc        =  1.0 + xsq * (`oneSixth + u0_div_H * `oneThird + `oneSixth * (xsq * (0.05 + 0.2 * u0_div_H)));
                            igcd_h     =  0.5 * igc - `oneSixth * (x * (1.0 + xsq * (0.4 * (Bg + 0.25) + 0.0285714285714 * (xsq * (0.125 + Bg)))));
                        end else begin
                            inv_x      =  1.0 / x;
                            `expl(x, ex)
                            inv_ex     =  1.0 / ex;
                            temp       =  ex - inv_ex;
                            temp2      =  ex + inv_ex;
                            igc        =  0.5 * ((1.0 - u0_div_H) * temp * inv_x + u0_div_H * temp2);
                            igcd_h     =  0.5 * (igc - temp * (Bg - Ag * inv_x * inv_x) - Ag * temp2 * inv_x);
                        end
                    end
                    Sg         =  0.5 * (1.0 + xg / sqrt(xg * xg + 1.0e-6));
                    Igc        =  Igc0 * igc * Sg;
                    Igcd       =  Igc0 * igcd_h * Sg;
                    Igcs       =  Igc - Igcd;
                    Igb        =  Igc0 * igc * (1.0 - Sg);
                end // (IGINV >0)
            end // (SWIGATE != 0)

            // 4.2.14 GIDL/GISL current
            Igidl        = 0.0;
            Igisl        = 0.0;
            if (SWGIDL != 0) begin

                // GIDL current computation
                if ((AGIDLD_i > 0) && (Vovd < 0)) begin
                    Vtovd        = sqrt(Vovd * Vovd + CGIDLD_i * CGIDLD_i * (VdbPrime * VdbPrime) + 1.0e-6);
                    temp = -BGIDLDs / Vtovd;
                    `expl_low(temp, temp2)
                    Igidl        = -AGIDLDs * (VdbPrime * Vovd * Vtovd * temp2);
                end

                // GISL current computation
                if ((AGIDL_i > 0) && (Vovs < 0)) begin
                    Vtovs        = sqrt(Vovs * Vovs + CGIDL_i * CGIDL_i * (VsbPrime * VsbPrime) + 1.0e-6);
                    temp = -BGIDLs / Vtovs;
                    `expl_low(temp, temp2)
                    Igisl        = -AGIDLs * (VsbPrime * Vovs * Vtovs * temp2);
                end
            end // (SWGIDL != 0)

            P_D        = 1 + 0.25 * (Gf * kp);
            facvsb0    = phib + 2 * phit1;
            facvsb     = Vsbstar + facvsb0;
            vts_i      = VFB_i + P_D * facvsb - Vsbstar + Gf * sqrt(phit1 * facvsb );
            vth_i      = vts_i - delVg;
        end // evaluateStatic


        /////////////////////////////////////////////////////////////////////////////
        //
        //      AC bias dependent quantities (calculations for charge contribs)
        //
        /////////////////////////////////////////////////////////////////////////////



        begin : evaluateDynamic

            // 4.2.16 Quantum mechanical corrections
            Vgb        =  Vgs + Vsb;
            COX_qm     =  COX_i;
            if (qq > 0.0) begin
                COX_qm     =  COX_i / (1.0 + qq * pow(qeff * qeff + qlim2, -1.0 * `oneSixth));
            end

            // 4.2.17 Intrinsic charge model
            if (xg <= 0.0) begin
                QG         =  Voxm;
                QI         =  0.0;
                QD         =  0.0;
                QB         =  QG;
            end else begin
                Fj         =  0.5 * (dps / H);
                Fj2        =  Fj * Fj;
                QCLM       =  (1.0 - GdL) * (qim - 0.5 * (alpha * dps));
                QG         =  Voxm + 0.5 * (eta_p * dps * (Fj * GdL * `oneThird - 1.0 + GdL));
                temp       =  alpha * dps * `oneSixth;
                QI         =  GdL * (qim + temp * Fj) + QCLM;
                QD         =  0.5 * (GdL * GdL * (qim - temp * (1.0 - Fj - 0.2 * Fj2)) + QCLM * (1.0 + GdL));
                QB         =  QG - QI;
            end
            Qg         =  QG * COX_qm;
            Qd         = -QD * COX_qm;
            Qb         = -QB * COX_qm;

            // 4.2.18 Extrinsic charge model
            Qgs_ov     =  CGOV_i * Vovs;
            Qgd_ov     =  CGOVD_i * Vovd;
            Qgb_ov     =  CGBOV_i * Vgb;

            // Outer fringe charge
            Qfgs       =  CFR_i * VgsPrime;
            Qfgd       =  CFRD_i * VgdPrime;

`ifdef NQSmodel
                // Variables for NQS model
                Gp         =  0.0;
                Gp2        =  0.0;
                a_factrp   =  0.0;
                marginp    =  0.0;
                if (SWNQS_i != 0) begin
                    if (xg <= 0.0) begin
                        ym         =  0.5;
                        pd         =  1.0;
                        Gp         =  Gf;
                    end else begin
                        ym         =  0.5 * ( 1.0 + 0.25 * (dps / H));
                        pd         =  xgm / (xg - x_m);
                        Gp         =  Gf / pd;
                    end
                    Gp2        =  Gp * Gp;
                    a_factrp   =  1.0 + Gp * `invSqrt2;
                    marginp    =  1e-5 * a_factrp;
                end
`endif // NQSmodel

        end // evaluateDynamic


        /////////////////////////////////////////////////////////////////////////////
        //
        //      JUNCAP2 contribs
        //
        /////////////////////////////////////////////////////////////////////////////


        ijun_s     = 0.0;
        ijunbot_s  = 0.0;
        ijunsti_s  = 0.0;
        ijungat_s  = 0.0;

        ijun_d     = 0.0;
        ijunbot_d  = 0.0;
        ijunsti_d  = 0.0;
        ijungat_d  = 0.0;

        qjun_s     = 0.0;
        qjunbot_s  = 0.0;
        qjunsti_s  = 0.0;
        qjungat_s  = 0.0;

        qjun_d     = 0.0;
        qjunbot_d  = 0.0;
        qjunsti_d  = 0.0;
        qjungat_d  = 0.0;

        begin : evaluateStaticDynamic

            if (SWJUNCAP > 0.0) begin
                if (SWJUNEXP_i == 1.0) begin
                    `JuncapExpressCurrent(Vjun_s, MFOR1_s, ISATFOR1_s, MFOR2_s, ISATFOR2_s, MREV_s, ISATREV_s, m0flag_s, xhighf1_s, expxhf1_s, xhighf2_s, expxhf2_s, xhighr_s, expxhr_s, ijun_s)
                    `JuncapExpressCurrent(Vjun_d, MFOR1_d, ISATFOR1_d, MFOR2_d, ISATFOR2_d, MREV_d, ISATREV_d, m0flag_d, xhighf1_d, expxhf1_d, xhighf2_d, expxhf2_d, xhighr_d, expxhr_d, ijun_d)
                    begin : evaluateDynamic
//                        real tmpv, vjv;
                        `JuncapExpressCharge(Vjun_s, ABSOURCE_i, LSSOURCE_i, LGSOURCE_i, qprefbot,   qprefsti,   qprefgat,   qpref2bot,   qpref2sti,   qpref2gat,   vbiinvbot,   vbiinvsti,   vbiinvgat,   one_minus_PBOT,   one_minus_PSTI,   one_minus_PGAT,   vfmin_s, vch_s, zflagbot_s, zflagsti_s, zflaggat_s, qjunbot_s, qjunsti_s, qjungat_s)
                        `JuncapExpressCharge(Vjun_d, ABDRAIN_i,  LSDRAIN_i,  LGDRAIN_i,  qprefbot_d, qprefsti_d, qprefgat_d, qpref2bot_d, qpref2sti_d, qpref2gat_d, vbiinvbot_d, vbiinvsti_d, vbiinvgat_d, one_minus_PBOT_d, one_minus_PSTI_d, one_minus_PGAT_d, vfmin_d, vch_d, zflagbot_d, zflagsti_d, zflaggat_d, qjunbot_d, qjunsti_d, qjungat_d)
                    end
                end else begin
                    `juncapcommon(Vjun_s, ABSOURCE_i, LSSOURCE_i, LGSOURCE_i, qprefbot, qpref2bot, vbiinvbot, one_minus_PBOT, idsatbot, CSRHBOT_i, CTATBOT_i, vbibot, wdepnulrbot, VBIRBOTinv, PBOT_i, ftdbot, btatpartbot, atatbot, one_over_one_minus_PBOT, CBBTBOT_i, VBIRBOT_i, wdepnulrinvbot, fbbtbot, VBRBOT_i, VBRinvbot, PBRBOT_i, fstopbot, slopebot, qprefsti, qpref2sti, vbiinvsti, one_minus_PSTI, idsatsti, CSRHSTI_i, CTATSTI_i, vbisti, wdepnulrsti, VBIRSTIinv, PSTI_i, ftdsti, btatpartsti, atatsti, one_over_one_minus_PSTI, CBBTSTI_i, VBIRSTI_i, wdepnulrinvsti, fbbtsti, VBRSTI_i, VBRinvsti, PBRSTI_i, fstopsti, slopesti, qprefgat, qpref2gat, vbiinvgat, one_minus_PGAT, idsatgat, CSRHGAT_i, CTATGAT_i, vbigat, wdepnulrgat, VBIRGATinv, PGAT_i, ftdgat, btatpartgat, atatgat, one_over_one_minus_PGAT, CBBTGAT_i, VBIRGAT_i, wdepnulrinvgat, fbbtgat, VBRGAT_i, VBRinvgat, PBRGAT_i, fstopgat, slopegat, VMAX_s, exp_VMAX_over_phitd_s, vbimin_s, vch_s, vfmin_s, vbbtlim_s, ijunbot_s, qjunbot_s, ijunsti_s, qjunsti_s, ijungat_s, qjungat_s)
                    ijun_s = ABSOURCE_i * ijunbot_s + LSSOURCE_i * ijunsti_s + LGSOURCE_i * ijungat_s;
                    `juncapcommon(Vjun_d, ABDRAIN_i,  LSDRAIN_i,  LGDRAIN_i,  qprefbot_d, qpref2bot_d, vbiinvbot_d, one_minus_PBOT_d, idsatbot_d, CSRHBOTD_i, CTATBOTD_i, vbibot_d, wdepnulrbot_d, VBIRBOTinv_d, PBOTD_i, ftdbot_d, btatpartbot_d, atatbot_d, one_over_one_minus_PBOT_d, CBBTBOTD_i, VBIRBOTD_i, wdepnulrinvbot_d, fbbtbot_d, VBRBOTD_i, VBRinvbot_d, PBRBOTD_i, fstopbot_d, slopebot_d, qprefsti_d, qpref2sti_d, vbiinvsti_d, one_minus_PSTI_d, idsatsti_d, CSRHSTID_i, CTATSTID_i, vbisti_d, wdepnulrsti_d, VBIRSTIinv_d, PSTID_i, ftdsti_d, btatpartsti_d, atatsti_d, one_over_one_minus_PSTI_d, CBBTSTID_i, VBIRSTID_i, wdepnulrinvsti_d, fbbtsti_d, VBRSTID_i, VBRinvsti_d, PBRSTID_i, fstopsti_d, slopesti_d, qprefgat_d, qpref2gat_d, vbiinvgat_d, one_minus_PGAT_d, idsatgat_d, CSRHGATD_i, CTATGATD_i, vbigat_d, wdepnulrgat_d, VBIRGATinv_d, PGATD_i, ftdgat_d, btatpartgat_d, atatgat_d, one_over_one_minus_PGAT_d, CBBTGATD_i, VBIRGATD_i, wdepnulrinvgat_d, fbbtgat_d, VBRGATD_i, VBRinvgat_d, PBRGATD_i, fstopgat_d, slopegat_d, VMAX_d, exp_VMAX_over_phitd_d, vbimin_d, vch_d, vfmin_d, vbbtlim_d, ijunbot_d, qjunbot_d, ijunsti_d, qjunsti_d, ijungat_d, qjungat_d)
                    ijun_d = ABDRAIN_i * ijunbot_d + LSDRAIN_i * ijunsti_d + LGDRAIN_i * ijungat_d;
                end
            end

`ifdef NQSmodel
                // Set initial conditions for NQS model
            `include "PSP102_InitNQS.include"

`endif // NQSmodel
            // Parasitic resistances (including noise)
            rgatenoise = nt0 * ggate;
            rbulknoise = nt0 * gbulk;
            rjunsnoise = nt0 * gjuns;
            rjundnoise = nt0 * gjund;
            rwellnoise = nt0 * gwell;
        end // evaluateStaticDynamic


        /////////////////////////////////////////////////////////////////////////////
        //
        //      Current contribs
        //
        /////////////////////////////////////////////////////////////////////////////

        begin : loadStatic

            // 4.2.15 Total terminal currents

            // Convert back for NMOS-PMOS and Source-Drain interchange
            if (sigVds > 0) begin
                I(D, BP)     <+  CHNL_TYPE * MULT_i * Iimpact;
                I(D, S)      <+  CHNL_TYPE * MULT_i * Ids;
                I(GP, S)     <+  CHNL_TYPE * MULT_i * Igcs;
                I(GP, D)     <+  CHNL_TYPE * MULT_i * Igcd;
            end else begin
                I(S, BP)     <+  CHNL_TYPE * MULT_i * Iimpact;
                I(S, D)      <+  CHNL_TYPE * MULT_i * Ids;
                I(GP, D)     <+  CHNL_TYPE * MULT_i * Igcs;
                I(GP, S)     <+  CHNL_TYPE * MULT_i * Igcd;
            end
            I(GP, BP)    <+  CHNL_TYPE * MULT_i * Igb;
            I(GP,  S)    <+  CHNL_TYPE * MULT_i * Igsov;
            I(GP,  D)    <+  CHNL_TYPE * MULT_i * Igdov;
            I(S,  BP)    <+  CHNL_TYPE * MULT_i * Igisl;
            I(D,  BP)    <+  CHNL_TYPE * MULT_i * Igidl;
            I(BS,  S)    <+  CHNL_TYPE * MULT_i * ijun_s;
            I(BD,  D)    <+  CHNL_TYPE * MULT_i * ijun_d;

            `NonCollapsableR(ggate, RG_i,    rgatenoise, G,  GP)
            `NonCollapsableR(gbulk, RBULK_i, rbulknoise, BP, BI)
            `NonCollapsableR(gjuns, RJUNS_i, rjunsnoise, BS, BI)
            `NonCollapsableR(gjund, RJUND_i, rjundnoise, BD, BI)
            `NonCollapsableR(gwell, RWELL_i, rwellnoise, B,  BI)

            I(D, S)  <+  `GMIN * V(D, S);

        end // loadStatic

        /////////////////////////////////////////////////////////////////////////////
        //
        //      ddt() contribs from charges (Note: MULT is handled explicitly)
        //
        /////////////////////////////////////////////////////////////////////////////

`ifdef NQSmodel
            begin : loadStaticDynamic
                // Calculate NQS charge contributions
            `include "PSP102_ChargesNQS.include"
            end
`endif // NQSmodel

        begin : loadDynamic

            // 4.2.19 Total terminal charges

            // Intrinsic MOSFET charges
            Qs         =  -(Qg + Qb + Qd);

            // Total outerFringe + overlap for
            //    gate-source and gate-drain.
            Qfgs       =  Qfgs + Qgs_ov;
            Qfgd       =  Qfgd + Qgd_ov;

            // JUNCAP2
            qjun_s = ABSOURCE_i * qjunbot_s + LSSOURCE_i * qjunsti_s + LGSOURCE_i * qjungat_s;
            qjun_d = ABDRAIN_i * qjunbot_d + LSDRAIN_i * qjunsti_d + LGDRAIN_i * qjungat_d;

            // Convert back (undo S-D interchange)
            if (sigVds < 0) begin
                temp       = Qd;    // Qd <--> Qs
                Qd         = Qs;
                Qs         = temp;
            end

            I(GP, S)    <+  ddt(CHNL_TYPE * MULT_i * Qg);
            I(BP, S)    <+  ddt(CHNL_TYPE * MULT_i * Qb);
            I(D,  S)    <+  ddt(CHNL_TYPE * MULT_i * Qd);
            I(GP, S)    <+  ddt(CHNL_TYPE * MULT_i * Qfgs);
            I(GP, D)    <+  ddt(CHNL_TYPE * MULT_i * Qfgd);
            I(GP, BP)   <+  ddt(CHNL_TYPE * MULT_i * Qgb_ov);
            I(BS, S)    <+  ddt(CHNL_TYPE * MULT_i * qjun_s);
            I(BD, D)    <+  ddt(CHNL_TYPE * MULT_i * qjun_d);

        end // loadDynamic


        /////////////////////////////////////////////////////////////////////////////
        //
        //  Noise
        //
        /////////////////////////////////////////////////////////////////////////////

        begin : noise
//            real sf;

            // 4.2.20 Noise variable calculation
            Sfl        =  0.0;
            mid        =  0.0;
            mig        =  0.0;
            migid      =  0.0;
            c_igid     =  0.0;
            CGeff      =  COX_qm * eta_p;
            sqid       =  0.0;
            sqig       =  0.0;
            if ((xg > 0.0) && (MULT_i > 0.0) && (BET_i > 0.0)) begin
                N1         =  Cox_over_q * alpha * phit;
                Nm1        =  Cox_over_q * qim1;
                Delta_N1   =  Cox_over_q * (alpha * dps);
                Sfl        =  (NFA_i - NFB_i * N1 + NFC_i * (N1 * N1)) * ln((Nm1 + 0.5 * Delta_N1) / (Nm1 - 0.5 * Delta_N1));
                Sfl        =  Sfl + (NFB_i + NFC_i * (Nm1 - 2.0 * N1)) * Delta_N1;
                Sfl        =  Sfl_prefac * Ids * Gvsatinv * Sfl / N1;
                Sfl        =  `CLIP_LOW(Sfl, 0.0);

                H0         =  qim1 / alpha;
                t1         =  qim / qim1;
                sqt2       =  0.5 * `oneSixth * (dps / H0);
                t2         =  sqt2 * sqt2;
                r          =  H0 / H - 1.0;
                lc         =  `CLIP_LOW(1.0 - 12 * (r * t2), 1e-20);
                lcinv2     =  1 / (lc * lc);
                g_ideal    =  BET_i * (FdL * qim1 * Gvsatinv);
                CGeff      =  Gvsat * Gvsat * COX_qm * eta_p / (Gmob_dL * Gmob_dL);
                mid        =  t1 + 12 * t2 - 24 * ((1.0 + t1) * t2 * r);
                mid        =  `CLIP_LOW(mid, 1e-40);
                mid        =  g_ideal * lcinv2 * mid;
                mig        =  t1 / 12.0 - t2 * (t1 + 0.2 - 12 * t2) - 1.6 * (t2 * (t1 + 1.0 - 12 * t2) * r);
                mig        =  `CLIP_LOW(mig, 1e-40);
                mig        =  lcinv2 / g_ideal * mig;
                migid      =  lcinv2 * sqt2 * (1.0 - 12 * t2 - (t1 + 19.2 * t2 - 12 * (t1 * t2)) * r);

                // excess noise calculation
                if (FNTEXC_i > 0) begin
                    // recalculate Gvsat, excluding Gmob-effect
                    temp2_exc  =  qim * xitsb;
                    wsat_exc   =  100.0 * (temp2_exc / (100.0 + temp2_exc));
                    if (THESATG_i < 0) begin
                        temp_exc   = 1 / (1 - THESATG_i * wsat_exc);
                    end else begin
                        temp_exc   = 1 + THESATG_i * wsat_exc;
                    end
                    thesat1_exc=  THESAT_i * (temp_exc / Gmob);
                    zsat_exc   =  thesat1_exc * thesat1_exc * dps * dps;
                    if (CHNL_TYPE == `PMOS) begin
                        zsat_exc   =  zsat_exc / (1.0 + thesat1_exc * dps);
                    end
                    Gvsat_exc  =  0.5 * (Gmob * (1.0 + sqrt(1.0 + 2.0 * zsat_exc)));

                    gfac       =  Gmob / (Gvsat_exc * lc);
                    Sidexc     =  fac_exc * Ids * Vdse * gfac * gfac;

                    // add excess noise to conventional thermal noise
                    mid        =  mid + Sidexc / nt0;
                    mig        =  mig + Sidexc * (1 + 12 * t2) / (12 * g_ideal * g_ideal * nt0);
                    migid      =  migid - Sidexc * sqt2 * (1 + r) / (g_ideal * nt0);
                end
                sqid       =  sqrt(nt * mid);
                sqig       =  sqrt(nt / mig);
                if (sqid == 0) begin
                    c_igid     =  0.0;
                end else begin
                    c_igid     =  migid * sqig / sqid; // = migid / sqrt(mig * mid);
                end
                c_igid     =  `CLIP_BOTH(c_igid, 0.0, 1.0);
            end
            shot_igcsx  = 2.0 * `QELE * abs(Igcs);
            shot_igcdx  = 2.0 * `QELE * abs(Igcd);
            shot_igsov  = 2.0 * `QELE * abs(Igsov);
            shot_igdov  = 2.0 * `QELE * abs(Igdov);
            shot_iavl   = 2.0 * `QELE * ((mavl + 1) * abs(Iimpact));
            // JUNCAP2
            jnoisex_s  = 2.0 * `QELE * abs(ijun_s);
            jnoisex_d  = 2.0 * `QELE * abs(ijun_d);
            if (sigVds > 0) begin
                shot_igs   =  shot_igcsx + shot_igsov;
                shot_igd   =  shot_igcdx + shot_igdov;
                jnoise_s   =  jnoisex_s;
                jnoise_d   =  jnoisex_d + shot_iavl;
            end else begin
                shot_igs   =  shot_igcdx + shot_igsov;
                shot_igd   =  shot_igcsx + shot_igdov;
                jnoise_s   =  jnoisex_s + shot_iavl;
                jnoise_d   =  jnoisex_d;
            end

            sf          = 1e6 * CGeff * mig; // scale-factor to keep noise at NOI2 at reasonable level

            // Important note:
            // In Verilog-A, correlated noise sources can only be implemented by using two additional
            // internal nodes (NOI and NOI2). When implementing PSP in a circuit simlutor, it is
            // generally not necessary to retain these internal nodes and therefore (for execution
            // speed reasons) should be avoided.

            // Noise contribs
            I(NOI2)   <+  V(NOI2);
            I(NOI2)   <+  white_noise(sqig * sqig * sf * sf, "igig");
            I(NOII)   <+  -V(NOI2) / sf;
            I(NOIR)   <+  V(NOIR) / mig;
            I(NOIC)   <+  ddt(CGeff * V(NOIC));
            I(D,S)    <+  flicker_noise(MULT_i * Sfl, EF_i, "flicker");
            I(D,S)    <+  white_noise(MULT_i * sqid * sqid * (1.0 - c_igid * c_igid), "idid");
//            I(D,S)    <+  sigVds * sqrt(MULT_i) * migid * V(NOI2) / sf;
//            I(GP,S)   <+  ddt(sqrt(MULT_i) * 0.5 * CGeff * V(NOIC));
//            I(GP,D)   <+  ddt(sqrt(MULT_i) * 0.5 * CGeff * V(NOIC));
            I(GP,S)   <+  white_noise(MULT_i * shot_igs, "igs");
            I(GP,D)   <+  white_noise(MULT_i * shot_igd, "igd");
            // JUNCAP2
            I(BS,S)   <+  white_noise(MULT_i * jnoise_s, "ibs");
            I(BD,D)   <+  white_noise(MULT_i * jnoise_d, "ibd");
        end // noise


        /////////////////////////////////////////////////////////////////////////////
        //
        // Operating point info
        //
        /////////////////////////////////////////////////////////////////////////////

        begin : OPinfo
//            real temp;

            // Auxiliary variables
            id_op     = Ids + Iimpact - Igcd;
            is        = -Ids - Igcs;
            ig        = Igcs + Igcd + Igsov + Igdov + Igb;
            ib        = -Iimpact - Igb - Igidl - Igisl;

            sig1k      = 2 * `PI * 1000 * CGeff;
            sig1k      = sig1k * sig1k * mig;


            ////////////////////////////////////////////////////////////////////////////////////
            //
            // Actual operation point output variables
            //
            ////////////////////////////////////////////////////////////////////////////////////

            // Note: In this section (and ONLY in this section) `drain' always refers to
            //       the highest-potential end of the channel. Therefore, care has to be
            //       taken for derivatives w.r.t. terminal voltages when sigVds == -1.

            sdint      = sigVds;

            ctype      = CHNL_TYPE;

            if (sigVds < 0) begin
                // All variables in the actual model refering to junctions are
                // not subject to SD-interchange. In the OP-output variables,
                // SD-interchange is also done for the junctions, so that's
                // what is happening here. Similar precautions have to be taken
                // for those variables that are derivatives w.r.t. voltage branches
                ise        = MULT_i * (is - Igdov + Igidl - ijun_d);
                ige        = MULT_i * ig;
                ide        = MULT_i * (id_op - Igsov + Igisl - ijun_s);
                ibe        = MULT_i * (ib + ijun_s + ijun_d);
                ids        = MULT_i * Ids;
                idb        = MULT_i * (Iimpact + Igisl - ijun_s);
                isb        = MULT_i * (Igidl - ijun_d);
                igs        = MULT_i * (Igcs + Igdov);
                igd        = MULT_i * (Igcd + Igsov);
                igb        = MULT_i * Igb;
                igcs       = MULT_i * Igcs;
                igcd       = MULT_i * Igcd;
                iavl       = MULT_i * Iimpact;
                igisl      = MULT_i * Igidl;
                igidl      = MULT_i * Igisl;
                if (SWJUNEXP_i == 1) begin
                    ijsbot     = 0.0;
                    ijsgat     = 0.0;
                    ijssti     = 0.0;
                    ijdbot     = 0.0;
                    ijdgat     = 0.0;
                    ijdsti     = 0.0;
                    idsatsbot  = 0.0;
                    idsatssti  = 0.0;
                    idsatsgat  = 0.0;
                    idsatsbotd = 0.0;
                    idsatsstid = 0.0;
                    idsatsgatd = 0.0;
                end else begin
                    ijsbot     = MULT_i * ABDRAIN_i * ijunbot_d;
                    ijsgat     = MULT_i * LGDRAIN_i * ijungat_d;
                    ijssti     = MULT_i * LSDRAIN_i * ijunsti_d;
                    ijdbot     = MULT_i * ABSOURCE_i * ijunbot_s;
                    ijdgat     = MULT_i * LGSOURCE_i * ijungat_s;
                    ijdsti     = MULT_i * LSSOURCE_i * ijunsti_s;
                    idsatsbot  = MULT_i * ABSOURCE_i * idsatbot;
                    idsatssti  = MULT_i * LSSOURCE_i * idsatsti;
                    idsatsgat  = MULT_i * LGSOURCE_i * idsatgat;
                    idsatsbotd = MULT_i * ABDRAIN_i * idsatbot_d;
                    idsatsstid = MULT_i * LSDRAIN_i * idsatsti_d;
                    idsatsgatd = MULT_i * LGDRAIN_i * idsatgat_d;
                end
                ijs        = MULT_i * ijun_d;
                ijd        = MULT_i * ijun_s;

                vds        = Vds;
                vgs        = Vgs;
                vsb        = Vsb;
                vto        = VFB_i + P_D * facvsb0 + Gf * sqrt(phit1 * facvsb0);
                vts        = vts_i;
                vth        = vth_i;
                vgt        = vgs - vth;
                vdss       = Vdsat;
                vsat       = Vds - vdss;

                temp       = Ids + Iimpact + Igisl - Igcd - Igsov - ijun_s; // Total drain-current
`ifdef OPderiv
                    gm         =  CHNL_TYPE * MULT_i * ddx(temp, V(GP));
                    gmb        =  CHNL_TYPE * MULT_i * ddx(temp, V(BP));
                    gds        =  CHNL_TYPE * MULT_i * ddx(temp, V(S));

                    gjs        =  MULT_i * ddx(ijun_d, V(BD));
                    gjd        =  MULT_i * ddx(ijun_s, V(BS));

                    css        =  CHNL_TYPE * MULT_i * ddx(Qd, V(D));
                    csg        = -CHNL_TYPE * MULT_i * ddx(Qd, V(GP));
                    csb        = -CHNL_TYPE * MULT_i * ddx(Qd, V(BP));
                    csd        =  css - csg - csb;
                    cgs        = -CHNL_TYPE * MULT_i * ddx(Qg, V(D));
                    cgg        =  CHNL_TYPE * MULT_i * ddx(Qg, V(GP));
                    cgb        = -CHNL_TYPE * MULT_i * ddx(Qg, V(BP));
                    cgd        =  cgg - cgs - cgb;
                    cds        = -CHNL_TYPE * MULT_i * ddx(Qs, V(D));
                    cdg        = -CHNL_TYPE * MULT_i * ddx(Qs, V(GP));
                    cdb        = -CHNL_TYPE * MULT_i * ddx(Qs, V(BP));
                    cdd        =  cdg + cds + cdb;
                    cbs        = -CHNL_TYPE * MULT_i * ddx(Qb, V(D));
                    cbg        = -CHNL_TYPE * MULT_i * ddx(Qb, V(GP));
                    cbb        =  CHNL_TYPE * MULT_i * ddx(Qb, V(BP));
                    cbd        =  cbb - cbs - cbg;
                    cgsol      =  CHNL_TYPE * MULT_i * ddx(Qfgd, V(GP));
                    cgdol      =  CHNL_TYPE * MULT_i * ddx(Qfgs, V(GP));

                    cjsbot     = -MULT_i * CHNL_TYPE * ABDRAIN_i * ddx(qjunbot_d, V(D));
                    cjsgat     = -MULT_i * CHNL_TYPE * LGDRAIN_i * ddx(qjungat_d, V(D));
                    cjssti     = -MULT_i * CHNL_TYPE * LSDRAIN_i * ddx(qjunsti_d, V(D));
                    cjs        =  cjsbot + cjsgat + cjssti;
                    cjdbot     = -MULT_i * CHNL_TYPE * ABSOURCE_i * ddx(qjunbot_s, V(S));
                    cjdgat     = -MULT_i * CHNL_TYPE * LGSOURCE_i * ddx(qjungat_s, V(S));
                    cjdsti     = -MULT_i * CHNL_TYPE * LSSOURCE_i * ddx(qjunsti_s, V(S));
                    cjd        =  cjdbot + cjdgat + cjdsti;
`endif // OPderiv
            end else begin
                ise        = MULT_i * (is - Igsov + Igisl - ijun_s);
                ige        = MULT_i * ig;
                ide        = MULT_i * (id_op - Igdov + Igidl - ijun_d);
                ibe        = MULT_i * (ib + ijun_s + ijun_d);
                ids        = MULT_i * Ids;
                idb        = MULT_i * (Iimpact + Igidl - ijun_d);
                isb        = MULT_i * (Igisl - ijun_s);
                igs        = MULT_i * (Igcs + Igsov);
                igd        = MULT_i * (Igcd + Igdov);
                igb        = MULT_i * Igb;
                igcs       = MULT_i * Igcs;
                igcd       = MULT_i * Igcd;
                iavl       = MULT_i * Iimpact;
                igisl      = MULT_i * Igisl;
                igidl      = MULT_i * Igidl;
                if (SWJUNEXP_i == 1) begin
                    ijsbot     = 0.0;
                    ijsgat     = 0.0;
                    ijssti     = 0.0;
                    ijdbot     = 0.0;
                    ijdgat     = 0.0;
                    ijdsti     = 0.0;
                    idsatsbot  = 0.0;
                    idsatssti  = 0.0;
                    idsatsgat  = 0.0;
                    idsatsbotd = 0.0;
                    idsatsstid = 0.0;
                    idsatsgatd = 0.0;
                end else begin
                    ijsbot     = MULT_i * ABSOURCE_i * ijunbot_s;
                    ijsgat     = MULT_i * LGSOURCE_i * ijungat_s;
                    ijssti     = MULT_i * LSSOURCE_i * ijunsti_s;
                    ijdbot     = MULT_i * ABDRAIN_i * ijunbot_d;
                    ijdgat     = MULT_i * LGDRAIN_i * ijungat_d;
                    ijdsti     = MULT_i * LSDRAIN_i * ijunsti_d;
                    idsatsbot  = MULT_i * ABSOURCE_i * idsatbot;
                    idsatssti  = MULT_i * LSSOURCE_i * idsatsti;
                    idsatsgat  = MULT_i * LGSOURCE_i * idsatgat;
                    idsatsbotd = MULT_i * ABDRAIN_i * idsatbot_d;
                    idsatsstid = MULT_i * LSDRAIN_i * idsatsti_d;
                    idsatsgatd = MULT_i * LGDRAIN_i * idsatgat_d;
                end
                ijs        = MULT_i * ijun_s;
                ijd        = MULT_i * ijun_d;

                vds        = Vds;
                vgs        = Vgs;
                vsb        = Vsb;
                vto        = VFB_i + P_D * facvsb0 + Gf * sqrt(phit1 * facvsb0);
                vts        = vts_i;
                vth        = vth_i;
                vgt        = vgs - vth;
                vdss       = Vdsat;
                vsat       = Vds - vdss;

                temp       = Ids + Iimpact + Igidl - Igcd - Igdov - ijun_d;
`ifdef OPderiv
                    gm         =  CHNL_TYPE * MULT_i * ddx(temp, V(GP));
                    gmb        =  CHNL_TYPE * MULT_i * ddx(temp, V(BP));
                    gds        =  CHNL_TYPE * MULT_i * ddx(temp, V(D));

                    gjs        = -MULT_i * ddx(ijun_s, V(S));
                    gjd        = -MULT_i * ddx(ijun_d, V(D));

                    cdd        =  CHNL_TYPE * MULT_i * ddx(Qd, V(D));
                    cdg        = -CHNL_TYPE * MULT_i * ddx(Qd, V(GP));
                    cdb        = -CHNL_TYPE * MULT_i * ddx(Qd, V(BP));
                    cds        =  cdd - cdg - cdb;
                    cgd        = -CHNL_TYPE * MULT_i * ddx(Qg, V(D));
                    cgg        =  CHNL_TYPE * MULT_i * ddx(Qg, V(GP));
                    cgb        = -CHNL_TYPE * MULT_i * ddx(Qg, V(BP));
                    cgs        =  cgg - cgd - cgb;
                    csd        = -CHNL_TYPE * MULT_i * ddx(Qs, V(D));
                    csg        = -CHNL_TYPE * MULT_i * ddx(Qs, V(GP));
                    csb        = -CHNL_TYPE * MULT_i * ddx(Qs, V(BP));
                    css        =  csg + csd + csb;
                    cbd        = -CHNL_TYPE * MULT_i * ddx(Qb, V(D));
                    cbg        = -CHNL_TYPE * MULT_i * ddx(Qb, V(GP));
                    cbb        =  CHNL_TYPE * MULT_i * ddx(Qb, V(BP));
                    cbs        =  cbb - cbd - cbg;
                    cgsol      =  CHNL_TYPE * MULT_i * ddx(Qfgs, V(GP));
                    cgdol      =  CHNL_TYPE * MULT_i * ddx(Qfgd, V(GP));

                    cjsbot     = -MULT_i * CHNL_TYPE * ABSOURCE_i * ddx(qjunbot_s, V(S));
                    cjsgat     = -MULT_i * CHNL_TYPE * LGSOURCE_i * ddx(qjungat_s, V(S));
                    cjssti     = -MULT_i * CHNL_TYPE * LSSOURCE_i * ddx(qjunsti_s, V(S));
                    cjs        =  cjsbot + cjsgat + cjssti;
                    cjdbot     = -MULT_i * CHNL_TYPE * ABDRAIN_i * ddx(qjunbot_d, V(D));
                    cjdgat     = -MULT_i * CHNL_TYPE * LGDRAIN_i * ddx(qjungat_d, V(D));
                    cjdsti     = -MULT_i * CHNL_TYPE * LSDRAIN_i * ddx(qjunsti_d, V(D));
                    cjd        =  cjdbot + cjdgat + cjdsti;
`endif // OPderiv
            end
`ifdef LocalModel
                weff       = 0;
                leff       = 0;
`else // LocalModel
                weff       = WE;
                leff       = LE;
`endif // LocalModel

`ifdef OPderiv
                if (abs(gds) < 1e-18) begin
                    u          = 0;
                    rout       = 0;
                    vearly     = 0;
                end else begin
                    u          = gm / gds;
                    rout       = 1.0 / gds;
                    vearly     = ide / gds;
                end
                if (abs(vgt) < 1e-12) begin
                    beff       = 0;
                end else begin
                    beff       = 2 * abs(ide) / (vgt * vgt);
                end
                if (abs(cgg + cgsol + cgdol) < 1e-30) begin
                    fug        = 0.0;
                end else begin
                    fug        = gm / (2 * `PI * (cgg + cgsol + cgdol));
                end
                rg_op         = RG_i / MULT_i;

                sfl        = MULT_i * Sfl;
                if (abs(gm) < 1e-18) begin
                    sqrtsff    = 0;
                    sqrtsfw    = 0;
                end else begin
                    sqrtsff    = sqrt(MULT_i * Sfl / 1000.0) / gm;
                    sqrtsfw    = sqrt(MULT_i) * sqid / gm;
                end
                sid        = MULT_i * sqid * sqid;
                sig        = MULT_i * nt * sig1k / (1 + sig1k * mig);
                cigid      = c_igid;
                if (sid == 0) begin
                    fknee      = 0;
                end else begin
                    fknee      = sfl / sid;
                end
                siavl      = MULT_i * shot_iavl;
                if (sigVds < 0) begin
                    sigs       = MULT_i * (shot_igcsx + shot_igdov);
                    sigd       = MULT_i * (shot_igcdx + shot_igsov);
                    ssi        = MULT_i * jnoisex_d;
                    sdi        = MULT_i * jnoisex_s;
                end else begin
                    sigs       = MULT_i * (shot_igcsx + shot_igsov);
                    sigd       = MULT_i * (shot_igcdx + shot_igdov);
                    ssi        = MULT_i * jnoisex_s;
                    sdi        = MULT_i * jnoisex_d;
                end
`endif // OPderiv

            lp_vfb       = VFB_i;
            lp_stvfb     = STVFB_i;
            lp_tox       = TOX_i;
            lp_epsrox    = EPSROX_i;
            lp_neff      = NEFF_i;
            lp_vnsub     = VNSUB_i;
            lp_nslp      = NSLP_i;
            lp_dnsub     = DNSUB_i;
            lp_dphib     = DPHIB_i;
            lp_np        = NP_i;
            lp_ct        = CT_i;
            lp_toxov     = TOXOV_i;
            lp_toxovd    = TOXOVD_i;
            lp_nov       = NOV_i;
            lp_novd      = NOVD_i;
            lp_cf        = CF_i;
            lp_cfb       = CFB_i;
            lp_betn      = BETN_i;
            lp_stbet     = STBET_i;
            lp_mue       = MUE_i;
            lp_stmue     = STMUE_i;
            lp_themu     = THEMU_i;
            lp_stthemu   = STTHEMU_i;
            lp_cs        = CS_i;
            lp_stcs      = STCS_i;
            lp_xcor      = XCOR_i;
            lp_stxcor    = STXCOR_i;
            lp_feta      = FETA_i;
            lp_rs        = RS_i;
            lp_strs      = STRS_i;
            lp_rsb       = RSB_i;
            lp_rsg       = RSG_i;
            lp_thesat    = THESAT_i;
            lp_stthesat  = STTHESAT_i;
            lp_thesatb   = THESATB_i;
            lp_thesatg   = THESATG_i;
            lp_ax        = AX_i;
            lp_alp       = ALP_i;
            lp_alp1      = ALP1_i;
            lp_alp2      = ALP2_i;
            lp_vp        = VP_i;
            lp_a1        = A1_i;
            lp_a2        = A2_i;
            lp_sta2      = STA2_i;
            lp_a3        = A3_i;
            lp_a4        = A4_i;
            lp_gco       = GCO_i;
            lp_iginv     = IGINV_i;
            lp_igov      = IGOV_i;
            lp_igovd     = IGOVD_i;
            lp_stig      = STIG_i;
            lp_gc2       = GC2_i;
            lp_gc3       = GC3_i;
            lp_chib      = CHIB_i;
            lp_agidl     = AGIDL_i;
            lp_agidld    = AGIDLD_i;
            lp_bgidl     = BGIDL_T;
            lp_bgidld    = BGIDLD_T;
            lp_stbgidl   = STBGIDL_i;
            lp_stbgidld  = STBGIDLD_i;
            lp_cgidl     = CGIDL_i;
            lp_cgidld    = CGIDLD_i;
            lp_cox       = COX_i;
            lp_cgov      = CGOV_i;
            lp_cgovd     = CGOVD_i;
            lp_cgbov     = CGBOV_i;
            lp_cfr       = CFR_i;
            lp_cfrd      = CFRD_i;
            lp_fnt       = FNT_i;
            lp_fntexc    = FNTEXC_i;
            lp_nfa       = NFA_i;
            lp_nfb       = NFB_i;
            lp_nfc       = NFC_i;
            lp_ef        = EF_i;
            lp_rg        = RG_i;
            lp_rbulk     = RBULK_i;
            lp_rwell     = RWELL_i;
            lp_rjuns     = RJUNS_i;
            lp_rjund     = RJUND_i;
            tk           = TKD;
            cjosbot      = MULT_i * ABSOURCE_i * cjobot;
            cjossti      = MULT_i * LSSOURCE_i * cjosti;
            cjosgat      = MULT_i * LGSOURCE_i * cjogat;
            vbisbot      = vbibot;
            vbissti      = vbisti;
            vbisgat      = vbigat;
            cjosbotd     = MULT_i * ABDRAIN_i * cjobot_d;
            cjosstid     = MULT_i * LSDRAIN_i * cjosti_d;
            cjosgatd     = MULT_i * LGDRAIN_i * cjogat_d;
            vbisbotd     = vbibot_d;
            vbisstid     = vbisti_d;
            vbisgatd     = vbigat_d;
`ifdef NQSmodel
                lp_munqs     = MUNQS_i;
`endif // NQSmodel
        end // OPinfo

    end // evaluateblock

end // analogBlock
