`undef P
`define P(txt) (*txt*) 
//======================================================================================
//======================================================================================
// Filename: PSP102_module.include
//======================================================================================
//======================================================================================
//
//  (c) Copyright 2007, All Rights Reserved, NXP Semiconductors
//
//
//  Version: 102.1, April 2007 (Simkit 2.5)
//
//======================================================================================
//======================================================================================
//
// Further information can be found in the file readme.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;

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

    // 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 real    TYPE     =  1.0      `from(   -1.0,1.0    ) `P(info="Channel type parameter, +1=NMOS -1=PMOS" unit="" );
    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    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    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    NOV      =  5.0e+25  `from(   1e20,1e27   ) `P(info="Effective doping of overlap region" unit="m^-3" );

    //  DIBL parameters
    parameter real    CF       =  0.0      `from(    0.0,inf    ) `P(info="DIBL-parameter" unit="V^-1" );
    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    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    BGIDL    = 41.0      `from(   0.0,inf     ) `P(info="GIDL probability factor at TR" unit="V" );
    parameter real    STBGIDL  =  0.0                             `P(info="Temperature dependence of BGIDL" unit="V/K" );
    parameter real    CGIDL    =  0.0                             `P(info="Back-bias dependence of GIDL" 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    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" );

    //  Noise parameters
    parameter real    FNT      =  1.0      `from(    0.0,inf    ) `P(info="Thermal 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" );
`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" );
    parameter real    RG       =  1.0e-3   `from( 1.0e-6,inf    ) `P(info="Gate resistance" unit="Ohm" );
    parameter real    RBULK    =  1.0e-3   `from( 1.0e-6,inf    ) `P(info="Bulk resistance between node BP and BI" unit="Ohm" );
    parameter real    RWELL    =  1.0e-3   `from( 1.0e-6,inf    ) `P(info="Well resistance between node BI and B" unit="Ohm" );
    parameter real    RJUNS    =  1.0e-3   `from( 1.0e-6,inf    ) `P(info="Source-side bulk resistance between node BI and BS" unit="Ohm" );
    parameter real    RJUND    =  1.0e-3   `from( 1.0e-6,inf    ) `P(info="Drain-side bulk resistance between node BI and BD" unit="Ohm" );
`endif // NQSmodel

    // JUNCAP Parameters
    parameter real    TRJ      = 21        `from(`TRJ_cliplow,inf) `P(info="reference temperature" unit="C" );
    `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    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 real    TYPE     =  1.0      `from(     -1,1      ) `P(info="Channel type parameter, +1=NMOS -1=PMOS" unit="" );

    // 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    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    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    LOV      =  0        `from(    0.0,inf    ) `P(info="Overlap length for gate/drain and gate/source overlap capacitance" unit="m" );
    parameter real    NOVO     =  5e25                            `P(info="Effective doping of overlap region" unit="m^-3" );

    // DIBL Parameters
    parameter real    CFL      =  0.0                             `P(info="Length dependence of DIBL-parameter" unit="V^-1" );
    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    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    BGIDLO   =  41                              `P(info="GIDL probability factor at TR" unit="V" );
    parameter real    STBGIDLO =  0.0                             `P(info="Temperature dependence of BGIDL" unit="V/K" );
    parameter real    CGIDLO   =  0.0                             `P(info="Back-bias dependence of GIDL" unit="" );

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

    // Noise Model Parameters
    parameter real    FNTO     =  1.0                             `P(info="Thermal noise coefficient" 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" );

    // Other Parameters
    parameter real    DTA      =  0                               `P(info="Temperature offset w.r.t. ambient circuit temperature" unit="K" );
`endif // Binning
`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    MUNQSO   =  1.0                             `P(info="Relative mobility for NQS modelling" );    
    parameter real    RGO      =  1.0e-3                          `P(info="Gate resistance" );                        
    parameter real    RBULKO   =  1.0e-3                          `P(info="Bulk resistance between node BP and BI" unit="Ohm" ); 
    parameter real    RWELLO   =  1.0e-3                          `P(info="Well resistance between node BI and B" unit="Ohm" );  
    parameter real    RJUNSO   =  1.0e-3                          `P(info="Source-side bulk resistance between node BI and BS" unit="Ohm" ); 
    parameter real    RJUNDO   =  1.0e-3                          `P(info="Drain-side bulk resistance between node BI and BD" unit="Ohm" ); 
`endif // NQSmodel

    // Stress Model Parameters
    parameter real    SAREF    =  1.0e-6   `from(   1e-9,inf    ) `P(info="Reference distance beteen OD-edge to poly from one side" unit="m" );
    parameter real    SBREF    =  1.0e-6   `from(   1e-9,inf    ) `P(info="Reference distance beteen OD-edge to 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 modifaction factor for stress effect" unit="" );

    // JUNCAP Parameters
    parameter real    TRJ      = 21        `from(`TRJ_cliplow,inf) `P(info="reference temperature" unit="C");
    `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 beteen OD-edge to poly from one side" unit="m" );
    parameter real    SB       =  0.0                             `P(type="instance" info="Distance beteen OD-edge to poly from other side" unit="m" );
    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="" );

    //////////////////////////
    //
    //  Variables
    //
    //////////////////////////
    
    // Variables for geometrical scaling rules
    real    L_i, W_i, SA_i, SB_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;
    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;
`endif // Binning    

    // List of local parameters
    real    VFB, STVFB, TOX, NEFF, VNSUB, NSLP, DNSUB, DPHIB, NP, CT;
    real    TOXOV, NOV, 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, STIG, GC2, GC3, CHIB;
    real    AGIDL, BGIDL, STBGIDL, CGIDL;
    real    COX, CGOV, CGBOV, CFR;
    real    FNT, NFA, NFB, NFC;
`ifdef NQSmodel
    real    MUNQS, RG, RBULK, RWELL, RJUNS, RJUND;
`endif // NQSmodel

    // 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, templ, tempw, Lx, Wx;
`endif // LocalModel

    // Variables used in electrical equations
    real VFB_i, STVFB_i, TOX_i, NEFF_i, VNSUB_i, NSLP_i, DNSUB_i, NP_i, QMC_i, CT_i, TOXOV_i, NOV_i;
    real CF_i, CFB_i, DPHIB_i;
    real BET_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, STIG_i, GC2_i, GC3_i, CHIB_i;
    real AGIDL_i, BGIDL_i, STBGIDL_i, CGIDL_i;
    real COX_i, CGOV_i, CGBOV_i, CFR_i;
    real FNT_i, NFA_i, NFB_i, NFC_i;
    real TR_i, MULT_i;

    real temp, temp1, temp2, tempM;
    real help;
    
    real TKR, TKD, TKD_sq, dT, rT, rTn;
    real 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 phib, sqrt_phib, phix, aphi, bphi, phix1, phix2, G_0, phit1, inv_phit1, alpha_b;
    real inv_VP, inv_AX, Sfl_prefac;

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

    real Dnsub;
    real Igidl, Igisl, Vtovd, Vtovs;
    real x_s, sqm, alpha, alpha1, eta_p, phi_inf, za, xitsb, rhob;
    real thesat1, wsat, ysat, zsat, r1, r2, dL, GdL, dL1, FdL, GR, Gmob, Gmob_dL, Gvsat, Gvsatinv, QCLM;
    real xgm, Voxm, dps, qim, qim1, qim1_1, xgs_ov, xgd_ov, sigVds;
    real Ux, xg;
    real mu, nu, xn_s, delta_ns;
    real Gf, Gf2, inv_Gf2, xi, inv_xi, margin;
    real qeff, COX_qm;

    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 Vdse, 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 H, Fj, Fj2;
    real N1, Nm1, Delta_N1, Sfl;
    real H0, t1, t2, sqt2, r, lc, lcinv2, g_ideal, CGeff, mid, mig, migid, c_igid, sqid, sqig;
    real shot_igs, shot_igsx, shot_igd, shot_igdx, shot_iavl;

    real Ids, Iimpact, mavl, Igdov, Igsov, Igc0, igc, igcd_h;
    real Igc, Igcd, Igcs, Igb, Igs, Igd;
    real Idse, Igbe, Igse, Igde, Igidle, Igisle, Iimpacte;
    real QI, QD, QB, QG, Qg, Qd, Qb, Qs, Qgs_ov, Qgd_ov;
    real Qfgs, Qfgd, Qgb_ov;

    real arg1, arg2max, arg2mina;

    integer CHNL_TYPE;
    
`ifdef NQSmodel
    // Variables used in NQS-calculations
    real SWNQS_i, MUNQS_i, RG_i, RBULK_i, RWELL_i, RJUNS_i, RJUND_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 Vrg, Vrbulk, Vrwell, Vrjund, Vrjuns;
    real ggate, gbulk, gwell, gjund, gjuns, nt0;
    real rgatenoise, rbulknoise, rwellnoise, rjundnoise, rjunsnoise;
    real temp3, temp4, temp5, temp6, temp7, temp8, temp9;
`endif // NQSmodel

    // JUNCAP2 variables
    `include "JUNCAP200_varlist.include"
    real isjunbot, qsjunbot, isjunsti, qsjunsti, isjungat, qsjungat, isjun, qsjun, sjnoise, sjnoisex;
    real idjunbot, qdjunbot, idjunsti, qdjunsti, idjungat, qdjungat, idjun, qdjun, djnoise, djnoisex;
    real Vjuns, Vjund, VMAXS, VMAXD;
    real vbimins, vchs, vfmins, vbbtlims, vbimind, vchd, vfmind, vbbtlimd;
    real ABSOURCE_i, LSSOURCE_i, LGSOURCE_i;
    real ABDRAIN_i, LSDRAIN_i, LGDRAIN_i, juncapwidth;



`ifdef insideADMS // OPinfo
    /////////////////////////////////////////////////////////////////////////////
    //
    // Variables for operating point info
    //
    /////////////////////////////////////////////////////////////////////////////

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

    real ctype     `P(ask="yes" info="Flag for channel type" unit="");
    real sdint     `P(ask="yes" info="Flag for source-drain interchange" unit="");
    
    real ise       `P(ask="yes" info="Total source current" unit="A");
    real ige       `P(ask="yes" info="Total gate current" unit="A");
    real ide       `P(ask="yes" info="Total drain current" unit="A");
    real ibe       `P(ask="yes" info="Total bulk current" unit="A");
    real ids       `P(ask="yes" info="Drain current, excl. avalanche, tunnel, GISL, GIDL, and junction currents" unit="A");
    real idb       `P(ask="yes" info="Drain to bulk current" unit="A");
    real isb       `P(ask="yes" info="Source to bulk current" unit="A");
    real igs       `P(ask="yes" info="Gate-source tunneling current" unit="A");
    real igd       `P(ask="yes" info="Gate-drain tunneling current" unit="A");
    real igb       `P(ask="yes" info="Gate-bulk tunneling current" unit="A");
    real igcs      `P(ask="yes" info="Gate-channel tunneling current (source component)" unit="A");
    real igcd      `P(ask="yes" info="Gate-channel tunneling current (drain component)" unit="A");
    real iavl      `P(ask="yes" info="Substrate current due to weak avelanche" unit="A");
    real igisl     `P(ask="yes" info="Gate-induced source leakage current" unit="A");
    real igidl     `P(ask="yes" info="Gate-induced drain leakage current" unit="A");

    real ijs       `P(ask="yes" info="Total source junction current" unit="A");
    real ijsbot    `P(ask="yes" info="Source junction current (bottom component)" unit="A");
    real ijsgat    `P(ask="yes" info="Source junction current (gate-edge component)" unit="A");
    real ijssti    `P(ask="yes" info="Source junction current (STI-edge component)" unit="A");
    real ijd       `P(ask="yes" info="Total drain junction current" unit="A");
    real ijdbot    `P(ask="yes" info="Drain junction current (bottom component)" unit="A");
    real ijdgat    `P(ask="yes" info="Drain junction current (gate-edge component)" unit="A");
    real ijdsti    `P(ask="yes" info="Drain junction current (STI-edge component)" unit="A");

    real vds       `P(ask="yes" info="Drain-source voltage" unit="V");
    real vgs       `P(ask="yes" info="Gate-source voltage" unit="V");
    real vsb       `P(ask="yes" info="Source-bulk voltage" unit="V");
    real vto       `P(ask="yes" info="Zero-bias threshold voltage" unit="V");
    real vts       `P(ask="yes" info="Threshold voltage including back bias effects" unit="V");
    real vth       `P(ask="yes" info="Threshold voltage including back bias and drain bias effects" unit="V");
    real vgt       `P(ask="yes" info="Effective gate drive voltage including back bias and drain bias effects" unit="V");
    real vdss      `P(ask="yes" info="Drain saturation voltage at actual bias" unit="V");
    real vsat      `P(ask="yes" info="Saturation limit" unit="V");
    
    real gm        `P(ask="yes" info="Transconductance" unit="1/Ohm");
    real gmb       `P(ask="yes" info="Substrate transconductance" unit="1/Ohm");
    real gds       `P(ask="yes" info="Output conductance" unit="1/Ohm");
    real gjs       `P(ask="yes" info="Source junction conductance" unit="1/Ohm");
    real gjd       `P(ask="yes" info="Drain junction conductance" unit="1/Ohm");
    
    real cdd       `P(ask="yes" info="Drain capacitance" unit="F");
    real cdg       `P(ask="yes" info="Drain-gate capacitance" unit="F");
    real cds       `P(ask="yes" info="Drain-source capacitance" unit="F");
    real cdb       `P(ask="yes" info="Drain-bulk capacitance" unit="F");
    real cgd       `P(ask="yes" info="Gate-drain capacitance" unit="F");
    real cgg       `P(ask="yes" info="Gate capacitance" unit="F");
    real cgs       `P(ask="yes" info="Gate-source capacitance" unit="F");
    real cgb       `P(ask="yes" info="Gate-bulk capacitance" unit="F");
    real csd       `P(ask="yes" info="Source-drain capacitance" unit="F");
    real csg       `P(ask="yes" info="Source-gate capacitance" unit="F");
    real css       `P(ask="yes" info="Source capacitance" unit="F");
    real csb       `P(ask="yes" info="Source-bulk capacitance" unit="F");
    real cbd       `P(ask="yes" info="Bulk-drain capacitance" unit="F");
    real cbg       `P(ask="yes" info="Bulk-gate capacitance" unit="F");
    real cbs       `P(ask="yes" info="Bulk-source capacitance" unit="F");
    real cbb       `P(ask="yes" info="Bulk capacitance" unit="F");
    real cgsol     `P(ask="yes" info="Total gate-source overlap capacitance" unit="F");
    real cgdol     `P(ask="yes" info="Total gate-drain overlap capacitance" unit="F");
    
    real cjs       `P(ask="yes" info="Total source junction capacitance" unit="F");
    real cjsbot    `P(ask="yes" info="Source junction capacitance (bottom component)" unit="F");
    real cjsgat    `P(ask="yes" info="Source junction capacitance (gate-edge component)" unit="F");
    real cjssti    `P(ask="yes" info="Source junction capacitance (STI-edge component)" unit="F");
    real cjd       `P(ask="yes" info="Total drain junction capacitance" unit="F");
    real cjdbot    `P(ask="yes" info="Drain junction capacitance (bottom component)" unit="F");
    real cjdgat    `P(ask="yes" info="Drain junction capacitance (gate-edge component)" unit="F");
    real cjdsti    `P(ask="yes" info="Drain junction capacitance (STI-edge component)" unit="F");
    
    real weff      `P(ask="yes" info="Effective channel width for geometrical models" unit="m");
    real leff      `P(ask="yes" info="Effective channel length for geometrical models" unit="m");
    real u         `P(ask="yes" info="Transistor gain" unit="");
    real rout      `P(ask="yes" info="Small-signal output resistance" unit="Ohm");
    real vearly    `P(ask="yes" info="Equivalent Early voltage" unit="V");
    real beff      `P(ask="yes" info="Gain factor" unit="A/V^2");
    real fug       `P(ask="yes" info="Unity gain frequency at actual bias" unit="Hz");

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

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

    analog begin

        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);
            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);
            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

            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);
`endif // LocalModel

            // 4.1 Internal parameters (including temperature scaling)
            // (only internal parameters independent on instance parameters
            //  are calculated in this section)
            if (TYPE >= 0) begin
                CHNL_TYPE  = `NMOS;
            end else begin
                CHNL_TYPE  = `PMOS;
            end

            // 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);

`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;

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

            // JUNCAP2
            `include "JUNCAP200_InitModel.include"

        end // initial_model

        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

`ifdef LocalModel
            juncapwidth= JW;

`else // LocalModel
            // Clipping of the instance parameters
            SAREF_i    = `CLIP_LOW(SAREF, 1e-9);
            SBREF_i    = `CLIP_LOW(SBREF, 1e-9);
            L_i        = `CLIP_LOW(L, 1e-9);
            W_i        = `CLIP_LOW(W, 1e-9);
            SA_i       = SA;
            SB_i       = SB;

            ///////////////////////////////////////////
            //  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;
            juncapwidth= WE;

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

`else // Binning
            // 3.3 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;
            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;
            NOV        = NOVO;

            // 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);
            STIG       = STIGO;
            GC2        = GC2O;
            GC3        = GC3O;
            CHIB       = CHIBO;
            
            // GIDL
            AGIDL      = AGIDLW * LOV_i / (LEN * iWE);
            BGIDL      = BGIDLO;
            STBGIDL    = STBGIDLO;
            CGIDL      = CGIDLO;
            
            // Charge model parameters
            COX        = `EPSOX * WEcv * LEcv / TOX_i;
            CGOV       = `EPSOX * WEcv * LOV_i / TOXOV_i;
            CGBOV      = CGBOVL * Lcv / LEN;
            CFR        = CFRW * Wcv / WEN;
            FNT        = FNTO;
            
            // Noise model parameters
            NFA        = iWE * iLE * NFALW;
            NFB        = iWE * iLE * NFBLW;
            NFC        = iWE * iLE * NFCLW;
`endif // Binning

`ifdef NQSmodel
            MUNQS      = MUNQSO;
            RG         = RGO;
            RWELL      = RWELLO;
            RBULK      = RBULKO;
            RJUNS      = RJUNSO;
            RJUND      = RJUNDO;
`endif // NQSModel

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

            // 3.5 Stress equations
            if ((SA_i > 0) && (SB_i > 0)) begin
                // Auxiliary variables
                Invsa      = 1.0 / (SA_i + 0.5 * L_i);
                Invsb      = 1.0 / (SB_i + 0.5 * L_i);
                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

            ///////////////////////////////////////////
            //  END OF SCALINGRULES AND STRESSMODEL
            ///////////////////////////////////////////

`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);
            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);
            NOV_i      = `CLIP_BOTH(NOV, 1e20, 1e27);
            CF_i       = `CLIP_LOW(CF, 0.0);
            CFB_i      = `CLIP_BOTH(CFB, 0.0, 1.0);
            BET_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);
            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);
            BGIDL_i    = `CLIP_LOW(BGIDL, 0.0);
            STBGIDL_i  =  STBGIDL;
            CGIDL_i    =  CGIDL;
            COX_i      = `CLIP_LOW(COX, 0.0);
            CGOV_i     = `CLIP_LOW(CGOV, 0.0);
            CGBOV_i    = `CLIP_LOW(CGBOV, 0.0);
            CFR_i      = `CLIP_LOW(CFR, 0.0);
            FNT_i      = `CLIP_LOW(FNT, 0.0);
            NFA_i      = `CLIP_LOW(NFA, 0.0);
            NFB_i      = `CLIP_LOW(NFB, 0.0);
            NFC_i      = `CLIP_LOW(NFC, 0.0);
            MULT_i     = `CLIP_LOW(MULT, 0.0);


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

            VFB_i      =  VFB_i + STVFB_i * dT;
            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(3.0e25, 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;
            GOV        =  sqrt(2.0 * `QELE * NOV_i * `EPSSI * inv_phit) / CoxovPrime;
            GOV2       =  GOV * GOV;
            xi_ov      =  1.0 + GOV * `invSqrt2;
            inv_xi_ov  =  1.0 / xi_ov;
            x_mrg_ov   =  1.0e-5 * xi_ov;

            // Mobility parameters
            tf_bet     =  pow(rTn, STBET_i);
            BET_i      =  BET_i * CoxPrime * 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;
            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;
            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
            AGIDL_i    =  AGIDL_i * 4e-18 / (TOXOV_i * TOXOV_i);
            tempM      = `MAX(1.0 + STBGIDL_i * dT, 0);
            BGIDL_i    =  BGIDL_i * tempM * TOXOV_i * 5e8;

            // Noise
            nt         =  FNT_i * 4 * `KBOL * TKD;
            Cox_over_q =  CoxPrime / `QELE;
            Sfl_prefac =  phit * phit * BET_i / Cox_over_q;

            // Additional internal parameters
            x1         =  1.25;
            inv_xg1    =  1.0 / (x1 + GOV * 7.324648775608221e-1); // =  1.0/(x1+GOV*sqrt(exp(-x1)+x1-1));
`ifdef NQSmodel

            // NQS parameters and variables
            MUNQS_i    = `CLIP_LOW(MUNQS, 0.0);
            RG_i       = `CLIP_LOW(RG, 1.0e-6);
            RBULK_i    = `CLIP_LOW(RBULK, 1.0e-6);
            RJUNS_i    = `CLIP_LOW(RJUNS, 1.0e-6);
            RJUND_i    = `CLIP_LOW(RJUND, 1.0e-6);
            RWELL_i    = `CLIP_LOW(RWELL, 1.0e-6);
            
            // Conductance of parasitic resistance
            ggate      =  MULT_i / RG_i;
            gbulk      =  MULT_i / RBULK_i;
            gjuns      =  MULT_i / RJUNS_i;
            gjund      =  MULT_i / RJUND_i;
            gwell      =  MULT_i / RWELL_i;
`endif // NQSmodel

            // JUNCAP2
            vbimins    = 0.0;
            vfmins     = 0.0;
            vchs       = 0.0;
            vbbtlims   = 0.0;
            vbimind    = 0.0;
            vfmind     = 0.0;
            vchd       = 0.0;
            vbbtlimd   = 0.0;
            vj         = 0.0;
            idmult     = 0.0;
            vjsrh      = 0.0;
            zinv       = 0.0;
            wdep       = 0.0;
            wsrh       = 0.0;
            asrh       = 0.0;
            vav        = 0.0;
            vbi_minus_vjsrh = 0.0;
            
            if (SWJUNCAP == 0.0) 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;
                VMAXS      = 0.0;
                VMAXD      = 0.0;
            end else begin
                if (SWJUNCAP == 2.0) begin
                    ABSOURCE_i = `CLIP_LOW(AS, `AB_cliplow); 
                    LSSOURCE_i = `CLIP_LOW(PS, `LS_cliplow); 
                    LGSOURCE_i = juncapwidth; 
                    ABDRAIN_i  = `CLIP_LOW(AD, `AB_cliplow); 
                    LSDRAIN_i  = `CLIP_LOW(PD, `LS_cliplow); 
                    LGDRAIN_i  = juncapwidth;
                end else begin
                    if (SWJUNCAP == 3.0) begin
                        ABSOURCE_i = `CLIP_LOW(AS, `AB_cliplow); 
                        ABDRAIN_i  = `CLIP_LOW(AD, `AB_cliplow); 
                        LSSOURCE_i = `CLIP_LOW(PS - juncapwidth, `LS_cliplow); 
                        LGSOURCE_i = juncapwidth; 
                        LSDRAIN_i  = `CLIP_LOW(PD - juncapwidth, `LS_cliplow); 
                        LGDRAIN_i  = juncapwidth;
                    end else begin
                        ABSOURCE_i = `CLIP_LOW(ABSOURCE, `AB_cliplow); 
                        LSSOURCE_i = `CLIP_LOW(LSSOURCE, `LS_cliplow); 
                        LGSOURCE_i = `CLIP_LOW(LGSOURCE, `LG_cliplow); 
                        ABDRAIN_i  = `CLIP_LOW(ABDRAIN, `AB_cliplow); 
                        LSDRAIN_i  = `CLIP_LOW(LSDRAIN, `LS_cliplow); 
                        LGDRAIN_i  = `CLIP_LOW(LGDRAIN, `LG_cliplow);
                    end
                end
                `JuncapInitInstance(ABSOURCE_i, LSSOURCE_i, LGSOURCE_i, VMAXS, vbimins, vchs, vfmins, vbbtlims)
                `JuncapInitInstance(ABDRAIN_i,  LSDRAIN_i,  LGDRAIN_i,  VMAXD, vbimind, vchd, vfmind, vbbtlimd)
            end
            

        end // initial_instance

        /////////////////////////////////////////////////////////////////////////////
        //
        //      DC bias dependent quantities (calculations for current contribs)
        //
        /////////////////////////////////////////////////////////////////////////////
        
        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;

            Vrg        =  V(G ,GP);
            Vrjuns     =  V(BS,BI);
            Vrjund     =  V(BD,BI);
            Vrbulk     =  V(BP,BI);
            Vrwell     =  V(B ,BI);
`endif // NQSmodel
            if (CHNL_TYPE == `NMOS) begin
                Vgs        =  V(`Gint, S);
                Vds        =  V(D, S);
                Vsb        =  V(S, `Bint);
                Vjuns      = -V(S, `Bjs);
                Vjund      = -V(D, `Bjd);
            end else begin
                Vgs        = -V(`Gint, S);
                Vds        = -V(D, S);
                Vsb        = -V(S, `Bint);
                Vjuns      =  V(S, `Bjs);
                Vjund      =  V(D, `Bjd);
            end
            
            // Source-drain interchange
            sigVds     =  1.0;
            if (Vds < 0.0) begin
                sigVds     = -1.0;
                Vgs        =  Vgs - Vds;
                Vsb        =  Vsb + Vds;
                Vds        = -Vds;
            end

            Vgd        =  Vgs - Vds;
            Vdb        =  Vds + Vsb;
            Vgb        =  Vgs + Vsb;

            xgs_ov     = -Vgs * inv_phit;
            xgd_ov     = -Vgd * inv_phit;
            
            // 4.2.1 Conditioning of terminal voltages
            temp       =  `MINA(Vdb, Vsb, bphi) + phix;
            Vsbstar    =  Vsb - `MINA(temp, 0, aphi) + phix1;
            Vdbstar    =  Vds + Vsbstar;
            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;
            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);
                        `expl_low(-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)) || ((SWGIDL != 0) && (AGIDL_i > 0)) || (CGOV_i > 0)) begin
                `sp_ov(xs_ov, xgs_ov)
                `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;
            Igs        =  0.0;
            Igd        =  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       =  -Vgs * 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)));

                    // 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       =  -Vgd * 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
                    Igdov      =  IGOV_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)
                Igs        =  Igsov + Igcs;
                Igd        =  Igdov + Igcd;
            end // (SWIGATE != 0)

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

                // GIDL current computation
                if (Vovd < 0) begin                    
                    Vtovd        = sqrt(Vovd * Vovd + CGIDL_i * CGIDL_i * (Vdb * Vdb) + 1.0e-6);
                    temp = -BGIDL_i / Vtovd;
                    `expl_low(temp, temp2)
                    Igidl        = -AGIDL_i * (Vdb * Vovd * Vtovd * temp2);
                end

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


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

        begin : evaluateDynamic

            // 4.2.16 Quantum mechanical corrections
            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     =  CGOV_i * Vovd;
            Qgb_ov     =  CGBOV_i * Vgb;

            // Outer fringe charge
            Qfgs       =  CFR_i * Vgs;
            Qfgd       =  CFR_i * Vgd;
`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
        //
        /////////////////////////////////////////////////////////////////////////////

        begin : evaluateStaticDynamic
            
            // Source side
            VAK     = Vjuns;
            VMAX    = VMAXS;
            vbimin  = vbimins;
            vfmin   = vfmins;
            vch     = vchs;
            vbbtlim = vbbtlims;
            `juncapcommon(ABSOURCE_i,LSSOURCE_i,LGSOURCE_i,isjunbot,qsjunbot,isjunsti,qsjunsti,isjungat,qsjungat)
            
            // Drain side
            VAK     = Vjund;
            VMAX    = VMAXD;
            vbimin  = vbimind;
            vfmin   = vfmind;
            vch     = vchd;
            vbbtlim = vbbtlimd;
            `juncapcommon(ABDRAIN_i,LSDRAIN_i,LGDRAIN_i,idjunbot,qdjunbot,idjunsti,qdjunsti,idjungat,qdjungat)

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

`endif // NQSmodel
        end // evaluateStaticDynamic


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

        begin : loadStatic

            // 4.2.15 Total terminal currents

            // Intrinsic MOSFET current
            Idse       = MULT_i * Ids;

            // Gate (tunneling) current components
            Igbe       = MULT_i * Igb;
            Igse       = MULT_i * Igs;
            Igde       = MULT_i * Igd;
            
            // GIDL/GISL current
            Igidle     = MULT_i * Igidl;
            Igisle     = MULT_i * Igisl;

            // Impact ionization current
            Iimpacte   = MULT_i * Iimpact;

            // JUNCAP2
            isjun = MULT_i * (ABSOURCE_i * isjunbot + LSSOURCE_i * isjunsti + LGSOURCE_i * isjungat);
            idjun = MULT_i * (ABDRAIN_i  * idjunbot + LSDRAIN_i  * idjunsti + LGDRAIN_i  * idjungat);
            
            // Convert back for NMOS-PMOS and Source-Drain interchange
            if (sigVds > 0) begin
                I(D, `Bint)     <+  CHNL_TYPE * (Iimpacte + Igidle);
                I(D, S)         <+  CHNL_TYPE * Idse;
                I(`Gint, S)     <+  CHNL_TYPE * Igse;
                I(`Gint, D)     <+  CHNL_TYPE * Igde;
                I(S, `Bint)     <+  CHNL_TYPE * Igisle;
            end else begin
                I(S, `Bint)     <+  CHNL_TYPE * (Iimpacte + Igidle);
                I(S, D)         <+  CHNL_TYPE * Idse;
                I(`Gint, D)     <+  CHNL_TYPE * Igse;
                I(`Gint, S)     <+  CHNL_TYPE * Igde;
                I(D, `Bint)     <+  CHNL_TYPE * Igisle;
            end
            I(`Gint, `Bint) <+  CHNL_TYPE * Igbe;
            I(`Bjs, S)      <+  CHNL_TYPE * isjun;
            I(`Bjd, D)      <+  CHNL_TYPE * idjun;
`ifdef NQSmodel
            I(G, GP)  <+  Vrg * ggate;
            I(BP, BI) <+  Vrbulk * gbulk;
            I(BS, BI) <+  Vrjuns * gjuns;
            I(BD, BI) <+  Vrjund * gjund;
            I(B, BI)  <+  Vrwell * gwell;
`endif // NQSmodel

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

        end // loadStatic

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

        begin : loadDynamic
`ifdef NQSmodel

            // Calculate NQS charge contributions
            `include "PSP102_ChargesNQS.include"
`endif // NQSmodel

            // 4.2.19 Total terminal charges

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

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

            // Gate-bulk overlap charge
            Qgb_ov     =  MULT_i * Qgb_ov;

            // JUNCAP2
            qsjun      =  MULT_i * (ABSOURCE_i * qsjunbot + LSSOURCE_i * qsjunsti + LGSOURCE_i * qsjungat);
            qdjun      =  MULT_i * (ABDRAIN_i  * qdjunbot + LSDRAIN_i  * qdjunsti + LGDRAIN_i  * qdjungat);

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

            I(`Gint, S)     <+  ddt(CHNL_TYPE * Qg);
            I(`Bint, S)     <+  ddt(CHNL_TYPE * Qb);
            I(D, S)         <+  ddt(CHNL_TYPE * Qd);
            I(`Gint, S)     <+  ddt(CHNL_TYPE * Qfgs);
            I(`Gint, D)     <+  ddt(CHNL_TYPE * Qfgd);
            I(`Gint, `Bint) <+  ddt(CHNL_TYPE * Qgb_ov);
            I(`Bjs, S)      <+  ddt(CHNL_TYPE * qsjun);
            I(`Bjd, D)      <+  ddt(CHNL_TYPE * qdjun);

        end // loadDynamic


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

        begin : noise

            // 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) && (BET_i > 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;

                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 - 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);
                sqid       =  sqrt(MULT_i * nt * mid);
                sqig       =  sqrt(MULT_i * nt / mig);
                c_igid     =  (sqid == 0) ? 0.0 : (migid * sqig / sqid); // = migid / sqrt(mig * mid);
                c_igid     =  `CLIP_BOTH(c_igid, 0.0, 1.0);
            end
            shot_igsx  = 2.0 * `QELE * abs(Igse);
            shot_igdx  = 2.0 * `QELE * abs(Igde);
            shot_iavl  = 2.0 * `QELE * ((mavl + 1) * abs(Iimpacte));
            // JUNCAP2
            sjnoisex   = 2.0 * `QELE * abs(isjun);
            djnoisex   = 2.0 * `QELE * abs(idjun);
            if (sigVds > 0) begin
                shot_igs   =  shot_igsx;
                shot_igd   =  shot_igdx;
                sjnoise    =  sjnoisex;
                djnoise    =  djnoisex + shot_iavl;
            end else begin
                shot_igs   =  shot_igdx;
                shot_igd   =  shot_igsx;
                sjnoise    =  sjnoisex + shot_iavl;
                djnoise    =  djnoisex;
            end
`ifdef NQSmodel
            rgatenoise = nt0 * ggate;
            rbulknoise = nt0 * gbulk;
            rjunsnoise = nt0 * gjuns;
            rjundnoise = nt0 * gjund;
            rwellnoise = nt0 * gwell;
`endif // NQSmodel

            // 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(c_igid);
            I(NOII)   <+  white_noise(sqig * sqig * (1.0 - c_igid));
            I(NOII)   <+  -sqig * V(NOI2);
            I(NOIR)   <+  V(NOIR);
            I(NOIC)   <+  ddt(mig * CGeff * V(NOIC));
            I(D,S)    <+  flicker_noise(MULT_i * Sfl, 1.0);
            I(D,S)    <+  white_noise(sqid * sqid * (1.0 - c_igid));
            I(D,S)    <+  sqid * V(NOI2);
            I(`Gint,S)<+  ddt(0.5 * ((1.0 + sigVds) * mig * CGeff * V(NOIC)));
            I(`Gint,D)<+  ddt(0.5 * ((1.0 - sigVds) * mig * CGeff * V(NOIC)));
            I(`Gint,S)<+  white_noise(shot_igs);
            I(`Gint,D)<+  white_noise(shot_igd);
            // JUNCAP2
            I(`Bjs,S) <+  white_noise(sjnoise, "shot");
            I(`Bjd,D) <+  white_noise(djnoise, "shot");
`ifdef NQSmodel
            // Parasitic resistances
            I(GP,G)   <+  white_noise(rgatenoise);
            I(BP,BI)  <+  white_noise(rbulknoise);
            I(BS,BI)  <+  white_noise(rjunsnoise);
            I(BD,BI)  <+  white_noise(rjundnoise);
            I(B ,BI)  <+  white_noise(rwellnoise);
`endif // NQSmodel
        end // noise


`ifdef insideADMS // OPinfo
        /////////////////////////////////////////////////////////////////////////////
        //
        // Operating point info
        //
        /////////////////////////////////////////////////////////////////////////////

        begin : OPinfo
        
            // The output variables defined below are currently not available in
            // Verilog-A, but only in the SiMKit-C-code which was generated from
            // this source. Similar functionality will be available in Verilog-A
            // from Verilog-A version 2.2 onwards. However, a different syntax is
            // to be used (see Verilog AMS language reference manual, version 2.2,
            // november 2004, Accellera).
            
            // Auxiliary variables
            id_op     = Idse + Iimpacte + Igidle - Igde;
            is        = -Idse + Igisle - Igse;
            ig        = Igse + Igde + Igbe;
            ib        = -Iimpacte - Igbe - Igidle - Igisle;

            P_D         = 1 + 0.25 * (Gf * kp);
            facvsb0    = phib + 2 * phit1;
            facvsb     = Vsbstar + facvsb0;
            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        = is - idjun;
                ige        = ig;
                ide        = id_op - isjun;
                ibe        = ib + isjun + idjun;
                ids        = Idse;
                idb        = Iimpacte + Igidle - isjun;
                isb        = Igisle - idjun;
                igs        = Igse;
                igd        = Igde;
                igb        = Igbe;
                igcs       = MULT_i * Igcs;
                igcd       = MULT_i * Igcd;
                iavl       = Iimpacte;
                igisl      = Igisle;
                igidl      = Igidle;

                ijsbot     = MULT_i * ABDRAIN_i * idjunbot;
                ijsgat     = MULT_i * LGDRAIN_i * idjungat;
                ijssti     = MULT_i * LSDRAIN_i * idjunsti;
                ijs        = ijsbot + ijsgat + ijssti;
                ijdbot     = MULT_i * ABSOURCE_i * isjunbot;
                ijdgat     = MULT_i * LGSOURCE_i * isjungat;
                ijdsti     = MULT_i * LSSOURCE_i * isjunsti;
                ijd        = ijdbot + ijdgat + ijdsti;

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

                temp       = Idse + Iimpacte + Igidle - Igde - isjun; // Total drain-current
                gm         =  CHNL_TYPE * ddx(temp, V(`Gint, S));
                gmb        = -CHNL_TYPE * ddx(temp, V(S, `Bint));
                gds        = -CHNL_TYPE * ddx(temp, V(D, S)) - (gm + gmb);

                gjs        = -ddx(idjun, V(D, `Bjd));
                gjd        = -ddx(isjun, V(S, `Bjs));

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

                cjsbot     = -MULT_i * CHNL_TYPE * ABDRAIN_i * ddx(qdjunbot, V(D, `Bjd));
                cjsgat     = -MULT_i * CHNL_TYPE * LGDRAIN_i * ddx(qdjungat, V(D, `Bjd));
                cjssti     = -MULT_i * CHNL_TYPE * LSDRAIN_i * ddx(qdjunsti, V(D, `Bjd));
                cjs        =  cjsbot + cjsgat + cjssti;
                cjdbot     = -MULT_i * CHNL_TYPE * ABSOURCE_i * ddx(qsjunbot, V(S, `Bjs));
                cjdgat     = -MULT_i * CHNL_TYPE * LGSOURCE_i * ddx(qsjungat, V(S, `Bjs));
                cjdsti     = -MULT_i * CHNL_TYPE * LSSOURCE_i * ddx(qsjunsti, V(S, `Bjs));
                cjd        =  cjdbot + cjdgat + cjdsti;                                     
            end else begin
                ise        = is - isjun;
                ige        = ig;
                ide        = id_op - idjun;
                ibe        = ib + isjun + idjun;
                ids        = Idse;
                idb        = Iimpacte + Igidle - idjun;
                isb        = Igisle - isjun;
                igs        = Igse;
                igd        = Igde;
                igb        = Igbe;
                igcs       = MULT_i * Igcs;
                igcd       = MULT_i * Igcd;
                iavl       = Iimpacte;
                igisl      = Igisle;
                igidl      = Igidle;

                ijsbot     = MULT_i * ABSOURCE_i * isjunbot;
                ijsgat     = MULT_i * LGSOURCE_i * isjungat;
                ijssti     = MULT_i * LSSOURCE_i * isjunsti;
                ijs        = ijsbot + ijsgat + ijssti;
                ijdbot     = MULT_i * ABDRAIN_i * idjunbot;
                ijdgat     = MULT_i * LGDRAIN_i * idjungat;
                ijdsti     = MULT_i * LSDRAIN_i * idjunsti;
                ijd        = ijdbot + ijdgat + ijdsti;

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

                temp       = Idse + Iimpacte + Igidle - Igde - idjun;
                gm         =  CHNL_TYPE * ddx(temp, V(`Gint, S));
                gmb        = -CHNL_TYPE * ddx(temp, V(S, `Bint));
                gds        =  CHNL_TYPE * ddx(temp, V(D, S));

                gjs        = -ddx(isjun, V(S, `Bjs));
                gjd        = -ddx(idjun, V(D, `Bjd));

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

                cjsbot     = -MULT_i * CHNL_TYPE * ABSOURCE_i * ddx(qsjunbot, V(S, `Bjs));
                cjsgat     = -MULT_i * CHNL_TYPE * LGSOURCE_i * ddx(qsjungat, V(S, `Bjs));
                cjssti     = -MULT_i * CHNL_TYPE * LSSOURCE_i * ddx(qsjunsti, V(S, `Bjs));
                cjs        =  cjsbot + cjsgat + cjssti;                                     
                cjdbot     = -MULT_i * CHNL_TYPE * ABDRAIN_i * ddx(qdjunbot, V(D, `Bjd));
                cjdgat     = -MULT_i * CHNL_TYPE * LGDRAIN_i * ddx(qdjungat, V(D, `Bjd));
                cjdsti     = -MULT_i * CHNL_TYPE * LSDRAIN_i * ddx(qdjunsti, V(D, `Bjd));
                cjd        =  cjdbot + cjdgat + cjdsti;
            end
`ifdef LocalModel
            weff       = 0;
            leff       = 0;
`else
            weff       = WE;
            leff       = LE;
`endif
            u          = (abs(gds) < 1e-18) ? 0 : (gm / gds);
            rout       = (abs(gds) < 1e-18) ? 0 : (1.0 / gds);
            vearly     = (abs(gds) < 1e-18) ? 0 : (ide / gds);
            beff       = (abs(vgt) < 1e-12) ? 0 : (2 * abs(ide) / (vgt * vgt));
            fug        = (abs(cgg + cgsol + cgdol) < 1e-30) ? 0.0 : gm / (2 * `PI * (cgg + cgsol + cgdol));

            sfl        = Sfl;
            sqrtsff    = (abs(gm) < 1e-18) ? 0 : (sqrt(MULT_i * Sfl / 1000) / gm);
            sqrtsfw    = (abs(gm) < 1e-18) ? 0 : (sqid / gm);
            sid        = sqid * sqid;
            sig        = MULT_i * nt * sig1k / (1 + sig1k * mig);
            cigid      = c_igid;
            fknee      = (sid == 0) ? 0 : Sfl / sid;
            sigs       = shot_igsx;
            sigd       = shot_igdx;
            siavl      = shot_iavl;
            if (sigVds < 0) begin
                ssi        = djnoisex;
                sdi        = sjnoisex;
            end else begin
                ssi        = sjnoisex;
                sdi        = djnoisex;
            end
        end // OPinfo
`endif // OPinfo

    end // analogBlock
