update vdmos with self-heating network and tj and tcase terminal

src/frontend/inpcom.c

@@ -4009,7 +4009,7 @@ static int get_number_terminals(char *c)
                 char *inst = gettok_instance(&c);
                 strncpy(nam_buf, inst, sizeof(nam_buf) - 1);
                 txfree(inst);
-                if (strstr(nam_buf, "off") || strchr(nam_buf, '='))
+                if (strstr(nam_buf, "off") || strchr(nam_buf, '=') || strstr(nam_buf, "tnodeout"))
                     break;
                 i++;
             }

src/spicelib/devices/devsup.c

@@ -620,7 +620,7 @@ http://ltwiki.org/index.php5?title=Undocumented_LTspice#VDMOS:_Breakdown_and_Sub
 void
 DevCapVDMOS(double vgd, double cgdmin,
             double cgdmax, double a, double cgs,
-            double *capgs, double *capgd, double *capgb)
+            double *capgs, double *capgd)
 {
     double s = (cgdmax - cgdmin) / (1 + M_PI / 2);
     double y = cgdmax - s;
@@ -629,7 +629,6 @@ DevCapVDMOS(double vgd, double cgdmin,
     else
         *capgd = 0.5 * (s * atan(a * vgd) + y);
     *capgs = 0.5 * cgs;
-    *capgb = 0;
 }
 
 /* Compute the MOS overlap capacitances as functions of the device

src/spicelib/devices/vdmos/Makefile.am

@@ -22,6 +22,7 @@ libvdmos_la_SOURCES =	\
 	vdmospar.c	\
 	vdmospzld.c	\
 	vdmosset.c	\
+	vdmossoachk.c	\
 	vdmostemp.c	\
 	vdmostrun.c
 

src/spicelib/devices/vdmos/vdmos.c

@@ -2,7 +2,7 @@
 Copyright 1990 Regents of the University of California.  All rights reserved.
 Author: 1987 Thomas L. Quarles
 Modified: 2000 AlansFixes
-VDMOS Model: 2018 Holger Vogt
+VDMOS: 2018 Holger Vogt, 2020 Dietmar Warning
 **********/
 
 #include "ngspice/ngspice.h"
@@ -12,35 +12,36 @@ VDMOS Model: 2018 Holger Vogt
 #include "ngspice/suffix.h"
 
 IFparm VDMOSpTable[] = { /* parameters */ 
- IOPU("mu",           VDMOS_M,          IF_REAL,    "Multiplier"),
- IOPU("l",            VDMOS_L,          IF_REAL,    "Length"),
- IOPU("w",            VDMOS_W,          IF_REAL,    "Width"),
- IP("off",            VDMOS_OFF,        IF_FLAG,    "Device initially off"),
- IOPU("icvds",        VDMOS_IC_VDS,     IF_REAL,    "Initial D-S voltage"),
- IOPU("icvgs",        VDMOS_IC_VGS,     IF_REAL,    "Initial G-S voltage"),
- IOPU("temp",         VDMOS_TEMP,       IF_REAL,    "Instance temperature"),
- IOPU("dtemp",        VDMOS_DTEMP,      IF_REAL,    "Instance temperature difference"),
- IP( "ic",            VDMOS_IC,         IF_REALVEC, "Vector of D-S, G-S voltages"),
-
- OP( "id",           VDMOS_CD,         IF_REAL,    "Drain current"),
- OP( "is",           VDMOS_CS,         IF_REAL,    "Source current"),
- OP( "ig",           VDMOS_CG,         IF_REAL,    "Gate current "),
- OP( "vgs",          VDMOS_VGS,        IF_REAL,    "Gate-Source voltage"),
- OP( "vds",          VDMOS_VDS,        IF_REAL,    "Drain-Source voltage"),
- OP( "cgs",          VDMOS_CAPGS,      IF_REAL,    "Gate-Source capacitance"),
- OP( "cgd",          VDMOS_CAPGD,      IF_REAL,    "Gate-Drain capacitance"),
- OP( "cds",          VDMOS_CAPDS,      IF_REAL,    "Drain-Source capacitance"),
+ IOPU("m",          VDMOS_M,          IF_REAL,    "Multiplier"),
+ IOPU("l",          VDMOS_L,          IF_REAL,    "Length"),
+ IOPU("w",          VDMOS_W,          IF_REAL,    "Width"),
+ IP("off",          VDMOS_OFF,        IF_FLAG,    "Device initially off"),
+ IOPU("icvds",      VDMOS_IC_VDS,     IF_REAL,    "Initial D-S voltage"),
+ IOPU("icvgs",      VDMOS_IC_VGS,     IF_REAL,    "Initial G-S voltage"),
+ IOPU("temp",       VDMOS_TEMP,       IF_REAL,    "Instance temperature"),
+ IOPU("dtemp",      VDMOS_DTEMP,      IF_REAL,    "Instance temperature difference"),
+ IP( "ic",          VDMOS_IC,         IF_REALVEC, "Vector of D-S, G-S voltages"),
+
+ IOP("tnodeout",    VDMOS_TNODEOUT,   IF_FLAG,    "Thermal model switch on/off"),
+
+ OP( "id",          VDMOS_CD,         IF_REAL,    "Drain current"),
+ OP( "is",          VDMOS_CS,         IF_REAL,    "Source current"),
+ OP( "ig",          VDMOS_CG,         IF_REAL,    "Gate current"),
+ OP( "vgs",         VDMOS_VGS,        IF_REAL,    "Gate-Source voltage"),
+ OP( "vds",         VDMOS_VDS,        IF_REAL,    "Drain-Source voltage"),
+ OP( "cgs",         VDMOS_CAPGS,      IF_REAL,    "Gate-Source capacitance"),
+ OP( "cgd",         VDMOS_CAPGD,      IF_REAL,    "Gate-Drain capacitance"),
+ OP( "cds",         VDMOS_CAPDS,      IF_REAL,    "Drain-Source capacitance"),
 
  OPU( "dnode",      VDMOS_DNODE,      IF_INTEGER, "Number of the drain node "),
  OPU( "gnode",      VDMOS_GNODE,      IF_INTEGER, "Number of the gate node "),
  OPU( "snode",      VDMOS_SNODE,      IF_INTEGER, "Number of the source node "),
+ OPU( "tempnode",   VDMOS_TNODE,      IF_INTEGER, "Number of temperature node"),
+ OPU( "tcasenode",  VDMOS_TCASE,      IF_INTEGER, "Number of 2nd temperature node"),
  OPU( "dnodeprime", VDMOS_DNODEPRIME, IF_INTEGER, "Number of int. drain node"),
  OPU( "snodeprime", VDMOS_SNODEPRIME, IF_INTEGER, "Number of int. source node "),
 
- OP( "von",               VDMOS_VON,           IF_REAL, " "),
- OP( "vdsat",             VDMOS_VDSAT,         IF_REAL, "Saturation drain voltage"),
- OPU( "sourcevcrit",      VDMOS_SOURCEVCRIT,   IF_REAL, "Critical source voltage"),
- OPU( "drainvcrit",       VDMOS_DRAINVCRIT,    IF_REAL, "Critical drain voltage"),
+ OP( "von",               VDMOS_VON,           IF_REAL, "Device on state voltage"),
  OP( "rs",                VDMOS_SOURCERESIST,  IF_REAL, "Source resistance"),
  OPU("sourceconductance", VDMOS_SOURCECONDUCT, IF_REAL, "Conductance of source"),
  OP( "rd",                VDMOS_DRAINRESIST,   IF_REAL, "Drain conductance"),
@@ -54,13 +55,14 @@ IFparm VDMOSpTable[] = { /* parameters */
 
  OPU( "qgs",      VDMOS_QGS,        IF_REAL,    "Gate-Source charge storage"),
  OPU( "qgd",      VDMOS_QGD,        IF_REAL,    "Gate-Drain charge storage"),
- OPU( "p",        VDMOS_POWER,      IF_REAL,    "Instaneous power"),
+ OPU( "p",        VDMOS_POWER,      IF_REAL,    "Instantaneous power"),
 };
 
 IFparm VDMOSmPTable[] = { /* model parameters */
  /* basic device */
  OP("type",   VDMOS_MOD_TYPE,  IF_STRING, "N-channel or P-channel MOS"),
- IOP("vto",   VDMOS_MOD_VTO,   IF_REAL,   "Threshold voltage"),
+ IOP("vto",   VDMOS_MOD_VTH,   IF_REAL,   "Threshold voltage"),
+ IOPR("vth0", VDMOS_MOD_VTH,   IF_REAL,   "Threshold voltage"),
  IOP("kp",    VDMOS_MOD_KP,    IF_REAL,   "Transconductance parameter"),
  IOP("phi",   VDMOS_MOD_PHI,   IF_REAL,   "Surface potential"),
  IOP("lambda",VDMOS_MOD_LAMBDA,IF_REAL,   "Channel length modulation"),
@@ -80,10 +82,17 @@ IFparm VDMOSmPTable[] = { /* model parameters */
  IOP("vq",      VDMOS_MOD_VQ,      IF_REAL,   "Quasi saturation voltage fitting parameter"),
  IOP("mtriode", VDMOS_MOD_MTRIODE, IF_REAL,   "Conductance multiplier in triode region"),
 
+ /* temperature dependency */
+ IOP( "tcvth", VDMOS_MOD_TCVTH, IF_REAL,  "Linear Vth0 temperature coefficient"),
+ IOPR("vtotc", VDMOS_MOD_TCVTH, IF_REAL,  "Linear Vth0 temperature coefficient"),
+ IOP( "mu",    VDMOS_MOD_MU,    IF_REAL,  "Exponent of gain temperature dependency"),
+ IOPR("bex",   VDMOS_MOD_MU,    IF_REAL,  "Exponent of gain temperature dependency"),
+ IOP( "texp0", VDMOS_MOD_TEXP0, IF_REAL,  "Drain resistance rd0 temperature exponent"),
+ IOP( "texp1", VDMOS_MOD_TEXP1, IF_REAL,  "Drain resistance rd1 temperature exponent"),
+
  /* weak inversion */
- IOP("subslope",  VDMOS_MOD_SUBSLOPE,   IF_REAL, "Slope of weak inversion log current versus vgs - vth"),
  IOP("subshift",  VDMOS_MOD_SUBSHIFT,   IF_REAL, "Shift of weak inversion plot on the vgs axis"),
- IOP("ksubthres", VDMOS_MOD_KSUBTHRES,  IF_REAL, "Shift of weak inversion plot on the vgs axis"),
+ IOP("ksubthres", VDMOS_MOD_KSUBTHRES,  IF_REAL, "Slope of weak inversion log current versus vgs"),
 
 /* body diode */
  IOP("bv",    VDMOS_MOD_BV,    IF_REAL,   "Vds breakdown voltage"),
@@ -94,7 +103,7 @@ IFparm VDMOSmPTable[] = { /* model parameters */
  IOP("n",     VDMOS_MOD_N,     IF_REAL,   "Body diode emission coefficient"),
  IOP("tt",    VDMOS_MOD_TT,    IF_REAL,   "Body diode transit time"),
  IOP("eg",    VDMOS_MOD_EG,    IF_REAL,   "Body diode activation energy for temperature effect on Is"),
- IOP("Xti",   VDMOS_MOD_XTI,   IF_REAL,   "Body diode saturation current temperature exponent"),
+ IOP("xti",   VDMOS_MOD_XTI,   IF_REAL,   "Body diode saturation current temperature exponent"),
  IOP("is",    VDMOS_MOD_IS,    IF_REAL,   "Body diode saturation current"),
  IOP("vj",    VDMOS_MOD_VJ,    IF_REAL,   "Body diode junction potential"),
 
@@ -108,16 +117,30 @@ IFparm VDMOSmPTable[] = { /* model parameters */
  IOPA("cgdmax", VDMOS_MOD_CGDMAX, IF_REAL, "Maximum non-linear G-D capacitance"),
  IOPA("a",      VDMOS_MOD_A,      IF_REAL, "Non-linear Cgd capacitance parameter"),
  IOPA("cgs",    VDMOS_MOD_CGS,    IF_REAL, "Gate-source capacitance"),
+
+/* self heating */
+ IOP("rthjc",  VDMOS_MOD_RTHJC,  IF_REAL, "Self-heating thermal resistance"),
+ IOP("rthca",  VDMOS_MOD_RTHCA,  IF_REAL, "Self-heating thermal resistance"),
+ IOP("cthj",   VDMOS_MOD_CTHJ,  IF_REAL, "Self-heating thermal capacitance"),
+ 
+/* soa check */
+ IOP("vgs_max", VDMOS_MOD_VGS_MAX, IF_REAL, "maximum voltage G-S branch"),
+ IOP("vgd_max", VDMOS_MOD_VGD_MAX, IF_REAL, "maximum voltage G-D branch"),
+ IOP("vds_max", VDMOS_MOD_VDS_MAX, IF_REAL, "maximum voltage D-S branch"),
+ IOP("vgsr_max", VDMOS_MOD_VGSR_MAX, IF_REAL, "maximum voltage G-S branch"),
+ IOP("vgdr_max", VDMOS_MOD_VGDR_MAX, IF_REAL, "maximum voltage G-D branch"),
 };
 
 char *VDMOSnames[] = {
     "Drain",
     "Gate",
-    "Source"
+    "Source",
+    "Temp",
+    "Tcase"
 };
 
-int	VDMOSnSize = NUMELEMS(VDMOSnames);
-int	VDMOSpTSize = NUMELEMS(VDMOSpTable);
-int	VDMOSmPTSize = NUMELEMS(VDMOSmPTable);
-int	VDMOSiSize = sizeof(VDMOSinstance);
-int	VDMOSmSize = sizeof(VDMOSmodel);
+int VDMOSnSize = NUMELEMS(VDMOSnames);
+int VDMOSpTSize = NUMELEMS(VDMOSpTable);
+int VDMOSmPTSize = NUMELEMS(VDMOSmPTable);
+int VDMOSiSize = sizeof(VDMOSinstance);
+int VDMOSmSize = sizeof(VDMOSmodel);

src/spicelib/devices/vdmos/vdmosacld.c

@@ -2,6 +2,7 @@
 Copyright 1990 Regents of the University of California.  All rights reserved.
 Author: 1985 Thomas L. Quarles
 Modified: 2000 AlansFixes
+VDMOS: 2018 Holger Vogt, 2020 Dietmar Warning
 **********/
 /*
  */
@@ -20,15 +21,21 @@ VDMOSacLoad(GENmodel *inModel, CKTcircuit *ckt)
     VDMOSinstance *here;
     int xnrm;
     int xrev;
-    double xgs;
-    double xgd;
-    double capgs;
-    double capgd;
+    double xgs, xcgT;
+    double xgd, xcdT;
+    double capgs, cgT;
+    double capgd, cdT;
+    double cTt, gTtt, gTtg, gTtdp, gTtsp;
+    double GmT;
+    double xcsT, xcTt;
+
+    register int selfheat;
 
     for( ; model != NULL; model = VDMOSnextModel(model)) {
         for(here = VDMOSinstances(model); here!= NULL;
                 here = VDMOSnextInstance(here)) {
-        
+
+            selfheat = (here->VDMOStnodeoutGiven) && (model->VDMOSrthjc != 0.0);
             if (here->VDMOSmode < 0) {
                 xnrm=0;
                 xrev=1;
@@ -36,6 +43,27 @@ VDMOSacLoad(GENmodel *inModel, CKTcircuit *ckt)
                 xnrm=1;
                 xrev=0;
             }
+
+            if (here->VDMOSmode >= 0) {
+                GmT = model->VDMOStype * here->VDMOSgmT;
+                cgT  = model->VDMOStype * here->VDMOScgT;
+                cdT  = model->VDMOStype * here->VDMOScdT;
+                cTt = model->VDMOScthj;
+                gTtg  = here->VDMOSgtempg;
+                gTtdp = here->VDMOSgtempd;
+                gTtt  = here->VDMOSgtempT;
+                gTtsp = - (gTtg + gTtdp);
+            } else {
+                GmT = -model->VDMOStype * here->VDMOSgmT;
+                cgT  = -model->VDMOStype * here->VDMOScgT;
+                cdT  = -model->VDMOStype * here->VDMOScdT;
+                cTt = -model->VDMOScthj;
+                gTtg  = -here->VDMOSgtempg;
+                gTtdp = -here->VDMOSgtempd;
+                gTtt  = -here->VDMOSgtempT;
+                gTtsp = gTtg + gTtdp;
+            }
+
             /*
              *     VDMOS cap model parameters
              */
@@ -46,7 +74,12 @@ VDMOSacLoad(GENmodel *inModel, CKTcircuit *ckt)
             xgs = capgs * ckt->CKTomega;
             xgd = capgd * ckt->CKTomega;
 
-            /* bulk diode */
+            xcgT = cgT * ckt->CKTomega;
+            xcdT = cdT * ckt->CKTomega;
+            xcsT = -(cgT + cdT) * ckt->CKTomega;
+            xcTt = cTt * ckt->CKTomega;
+
+            /* body diode */
             double gspr, geq, xceq;
             gspr = here->VDIOtConductance;
             geq = *(ckt->CKTstate0 + here->VDIOconduct);
@@ -62,6 +95,7 @@ VDMOSacLoad(GENmodel *inModel, CKTcircuit *ckt)
             *(here->VDMOSGPspPtr +1) -= xgs;
             *(here->VDMOSDPgpPtr +1) -= xgd;
             *(here->VDMOSSPgpPtr +1) -= xgs;
+
             *(here->VDMOSDdPtr) += here->VDMOSdrainConductance;
             *(here->VDMOSSsPtr) += here->VDMOSsourceConductance;
             *(here->VDMOSDPdpPtr) += here->VDMOSdrainConductance+
@@ -78,10 +112,10 @@ VDMOSacLoad(GENmodel *inModel, CKTcircuit *ckt)
             *(here->VDMOSSPdpPtr) -= here->VDMOSgds+xrev*(here->VDMOSgm);
             /* gate resistor */
             *(here->VDMOSGgPtr) += (here->VDMOSgateConductance);
-            *(here->VDMOSGPgpPtr) += (here->VDMOSgateConductance)/* + ?? FIXME */;
+            *(here->VDMOSGPgpPtr) += (here->VDMOSgateConductance);
             *(here->VDMOSGgpPtr) -= here->VDMOSgateConductance;
             *(here->VDMOSGPgPtr) -= here->VDMOSgateConductance;
-            /* bulk diode */
+            /* body diode */
             *(here->VDMOSSsPtr) += gspr;
             *(here->VDMOSDdPtr) += geq;
             *(here->VDMOSDdPtr +1) += xceq;
@@ -93,6 +127,29 @@ VDMOSacLoad(GENmodel *inModel, CKTcircuit *ckt)
             *(here->VDIORPsPtr) -= gspr;
             *(here->VDIORPdPtr) -= geq;
             *(here->VDIORPdPtr +1) -= xceq;
+            if (selfheat)
+            {
+               *(here->VDMOSDPtempPtr)       += GmT;
+               *(here->VDMOSSPtempPtr)       += -GmT;
+
+               *(here->VDMOSTemptempPtr)     += gTtt + 1/model->VDMOSrthjc;
+               *(here->VDMOSTempgpPtr)       += gTtg;
+               *(here->VDMOSTempdpPtr)       += gTtdp;
+               *(here->VDMOSTempspPtr)       += gTtsp;
+               *(here->VDMOSTemptcasePtr)    += -1/model->VDMOSrthjc;
+               *(here->VDMOSTcasetempPtr)    += -1/model->VDMOSrthjc;
+               *(here->VDMOSTcasetcasePtr)   +=  1/model->VDMOSrthjc + 1/model->VDMOSrthca;
+               *(here->VDMOSTptpPtr)         +=  1/model->VDMOSrthca;
+               *(here->VDMOSTptcasePtr)      += -1/model->VDMOSrthca;
+               *(here->VDMOSTcasetpPtr)      += -1/model->VDMOSrthca;
+               *(here->VDMOSCktTtpPtr)       +=  1.0;
+               *(here->VDMOSTpcktTPtr)       +=  1.0;
+
+               *(here->VDMOSTemptempPtr + 1) += xcTt;
+               *(here->VDMOSDPtempPtr + 1)   += xcdT;
+               *(here->VDMOSSPtempPtr + 1)   += xcsT;
+               *(here->VDMOSGPtempPtr + 1)   += xcgT;
+            }
         }
     }
     return(OK);

src/spicelib/devices/vdmos/vdmosask.c

@@ -43,19 +43,22 @@ VDMOSask(CKTcircuit *ckt, GENinstance *inst, int which, IFvalue *value,
             return(OK);
         case VDMOS_L:
             value->rValue = here->VDMOSl;
-                return(OK);
+            return(OK);
         case VDMOS_W:
             value->rValue = here->VDMOSw;
-                return(OK);
+            return(OK);
         case VDMOS_OFF:
-            value->rValue = here->VDMOSoff;
-                return(OK);
+            value->iValue = here->VDMOSoff;
+            return(OK);
+        case VDMOS_TNODEOUT:
+            value->iValue = here->VDMOStnodeout;
+            return(OK);
         case VDMOS_IC_VDS:
             value->rValue = here->VDMOSicVDS;
-                return(OK);
+            return(OK);
         case VDMOS_IC_VGS:
             value->rValue = here->VDMOSicVGS;
-                return(OK);
+            return(OK);
         case VDMOS_DNODE:
             value->iValue = here->VDMOSdNode;
             return(OK);
@@ -65,6 +68,9 @@ VDMOSask(CKTcircuit *ckt, GENinstance *inst, int which, IFvalue *value,
         case VDMOS_SNODE:
             value->iValue = here->VDMOSsNode;
             return(OK);
+        case VDMOS_TNODE:
+            value->iValue = here->VDMOStempNode;
+            return(OK);
         case VDMOS_SNODEPRIME:
             value->iValue = here->VDMOSsNodePrime;
             return(OK);
@@ -89,15 +95,6 @@ VDMOSask(CKTcircuit *ckt, GENinstance *inst, int which, IFvalue *value,
         case VDMOS_VON:
             value->rValue = here->VDMOSvon;
             return(OK);
-        case VDMOS_VDSAT:
-            value->rValue = here->VDMOSvdsat;
-            return(OK);
-        case VDMOS_SOURCEVCRIT:
-            value->rValue = here->VDMOSsourceVcrit;
-            return(OK);
-        case VDMOS_DRAINVCRIT:
-            value->rValue = here->VDMOSdrainVcrit;
-            return(OK);
         case VDMOS_CD:
             value->rValue = here->VDMOScd;
             return(OK);

src/spicelib/devices/vdmos/vdmosdefs.h

@@ -2,6 +2,7 @@
 Copyright 1990 Regents of the University of California.  All rights reserved.
 Author: 1985 Thomas L. Quarles
 Modified: 2000 AlansFixes
+VDMOS: 2018 Holger Vogt, 2020 Dietmar Warning
 **********/
 
 #ifndef VDMOS
@@ -41,10 +42,15 @@ typedef struct sVDMOSinstance {
     const int VDMOSdNode;  /* number of the gate node of the mosfet */
     const int VDMOSgNode;  /* number of the gate node of the mosfet */
     const int VDMOSsNode;  /* number of the source node of the mosfet */
+    const int VDMOStempNode;  /* number of the temperature node of the mosfet */
+    const int VDMOStcaseNode;  /* number of the 2nd temperature node of the mosfet */
     int VDMOSdNodePrime; /* number of the internal drain node of the mosfet */
     int VDMOSsNodePrime; /* number of the internal source node of the mosfet */
     int VDMOSgNodePrime; /* number of the internal gate node of the mosfet */
-    int VDIOposPrimeNode; /* number of the internal node of the bulk diode */
+    int VDMOStNodePrime; /* number of the internal temp node between voltage source and Rthca */
+    int VDIOposPrimeNode; /* number of the internal node of the body diode */
+
+    int VDMOSvcktTbranch; /* equation number of branch equation added for cktTemp source */ 
 
     double VDMOSm;   /* parallel device multiplier */
 
@@ -52,42 +58,37 @@ typedef struct sVDMOSinstance {
     double VDMOSw;   /* the width of the channel region */
     double VDMOSsourceConductance;   /*conductance of source(or 0):set in setup*/
     double VDMOSdrainConductance;    /*conductance of drain(or 0):set in setup*/
+    double VDMOSdrainResistance;    /*resistance of drain(or 0): set in temp*/
+    double VDMOSqsResistance;    /*resistance of drain: set in temp*/
     double VDMOSgateConductance;    /*conductance of gate(or 0):set in setup*/
     double VDMOSdsConductance;    /*conductance of drain to source:set in setup*/
     double VDMOStemp;    /* operating temperature of this instance */
     double VDMOSdtemp;   /* operating temperature of the instance relative to circuit temperature*/
+    int    VDMOStnodeout;  /* flag indicate self heating on */
 
     double VDMOStTransconductance;   /* temperature corrected transconductance*/
     double VDMOStPhi;                /* temperature corrected Phi */
-    double VDMOStVto;                /* temperature corrected Vto */
-    double VDMOStSatCur;             /* temperature corrected saturation Cur. */
+    double VDMOStVth;                /* temperature corrected Vth */
 
     double VDMOSicVDS;   /* initial condition D-S voltage */
     double VDMOSicVGS;   /* initial condition G-S voltage */
     double VDMOSvon;
     double VDMOSvdsat;
-    double VDMOSsourceVcrit; /* Vcrit for pos. vds */
-    double VDMOSdrainVcrit;  /* Vcrit for pos. vds */
     double VDMOScd;
     double VDMOSgm;
     double VDMOSgds;
-    double VDMOSf2d;
-    double VDMOSf3d;
-    double VDMOSf4d;
-    double VDMOSf2s;
-    double VDMOSf3s;
-    double VDMOSf4s;
 
     double VDIOcap;
-    double VDIOtSatCur; /* temperature corrected saturation Cur. density*/
+    double VDIOtSatCur;     /* temperature corrected saturation Cur. density*/
     double VDIOinitCond;
     double VDIOtVcrit;
+    double VDIOconductance;
     double VDIOtConductance;
     double VDIOtBrkdwnV;
     double VDIOtJctCap;
     double VDIOtDepCap;     /* temperature adjusted transition point in */
-                             /* the cureve matching Fc * Vj */
-    double VDIOtJctPot;    /* temperature corrected Bulk potential */
+                            /* the curve matching Fc * Vj */
+    double VDIOtJctPot;     /* temperature corrected junction potential */
     double VDIOtGradingCoeff;
 
     double VDIOtTransitTime;
@@ -95,6 +96,16 @@ typedef struct sVDMOSinstance {
     double VDIOtF2;
     double VDIOtF3;
 
+    double VDMOSTempSH;      /* for portability of SH temp to noise analysis */
+
+    double VDMOSgmT;
+    double VDMOSgtempg;
+    double VDMOSgtempd;
+    double VDMOSgtempT;
+    double VDMOScgT;
+    double VDMOScdT;
+    double VDMOScth;         /* current alias power */
+
 /*
  * naming convention:
  * x = vgs
@@ -135,21 +146,20 @@ typedef struct sVDMOSinstance {
 
     int VDMOSmode;       /* device mode : 1 = normal, -1 = inverse */
 
-
     unsigned VDMOSoff:1;  /* non-zero to indicate device is off for dc analysis*/
     unsigned VDMOStempGiven :1;  /* instance temperature specified */
     unsigned VDMOSdtempGiven :1;  /* instance delta temperature specified */
     unsigned VDMOSmGiven :1;
     unsigned VDMOSlGiven :1;
     unsigned VDMOSwGiven :1;
-    unsigned VDMOSdNodePrimeSet  :1;
-    unsigned VDMOSsNodePrimeSet  :1;
+    unsigned VDMOSdNodePrimeSet :1;
+    unsigned VDMOSsNodePrimeSet :1;
     unsigned VDMOSicVDSGiven :1;
     unsigned VDMOSicVGSGiven :1;
-    unsigned VDMOSvonGiven   :1;
+    unsigned VDMOStnodeoutGiven : 1; /* flag indicate self heating on */
+    unsigned VDMOSvonGiven : 1;
     unsigned VDMOSvdsatGiven :1;
-    unsigned VDMOSmodeGiven  :1;
-
+    unsigned VDMOSmodeGiven :1;
 
     double *VDMOSDdPtr;      /* pointer to sparse matrix element at
                                      * (Drain node,drain node) */
@@ -201,7 +211,7 @@ typedef struct sVDMOSinstance {
                              * (source node, drain node) */
     double *VDMOSSdPtr;    /* pointer to sparse matrix element at
                              * (drain node, source node) */
-    /* bulk diode */
+    /* body diode */
     double *VDIORPdPtr;    /* pointer to sparse matrix element at
                              * (diode prime node, drain node) */
     double *VDIODrpPtr;    /* pointer to sparse matrix element at
@@ -212,32 +222,57 @@ typedef struct sVDMOSinstance {
                             * (source node, diode prime node) */
     double *VDIORPsPtr;    /* pointer to sparse matrix element at
                             * (diode prime node, source node) */
+    /* self heating */
+    double *VDMOSTemptempPtr;
+    double *VDMOSTempdpPtr;
+    double *VDMOSTempspPtr;
+    double *VDMOSTempgpPtr;
+    double *VDMOSGPtempPtr;
+    double *VDMOSDPtempPtr;
+    double *VDMOSSPtempPtr;
+
+    double *VDMOSTcasetcasePtr; /* for Rthjc */
+    double *VDMOSTcasetempPtr;
+    double *VDMOSTemptcasePtr;
+
+    double *VDMOSTptpPtr;       /* for Rthca */
+    double *VDMOSTptcasePtr;
+    double *VDMOSTcasetpPtr;
+
+    double *VDMOSCktTcktTPtr;   /* for VcktTemp */
+    double *VDMOSCktTtpPtr;
+    double *VDMOSTpcktTPtr;
 
 } VDMOSinstance ;
 
 #define VDMOSvgs VDMOSstates+ 0   /* gate-source voltage */
 #define VDMOSvds VDMOSstates+ 1   /* drain-source voltage */
+#define VDMOSdeltemp VDMOSstates+ 2
 
-#define VDMOScapgs VDMOSstates+2  /* gate-source capacitor value */
-#define VDMOSqgs VDMOSstates+ 3   /* gate-source capacitor charge */
-#define VDMOScqgs VDMOSstates+ 4  /* gate-source capacitor current */
+#define VDMOScapgs VDMOSstates+3  /* gate-source capacitor value */
+#define VDMOSqgs VDMOSstates+ 4   /* gate-source capacitor charge */
+#define VDMOScqgs VDMOSstates+ 5  /* gate-source capacitor current */
 
-#define VDMOScapgd VDMOSstates+ 5 /* gate-drain capacitor value */
-#define VDMOSqgd VDMOSstates+ 6   /* gate-drain capacitor charge */
-#define VDMOScqgd VDMOSstates+ 7  /* gate-drain capacitor current */
+#define VDMOScapgd VDMOSstates+ 6 /* gate-drain capacitor value */
+#define VDMOSqgd VDMOSstates+ 7   /* gate-drain capacitor charge */
+#define VDMOScqgd VDMOSstates+ 8  /* gate-drain capacitor current */
 
-#define VDIOvoltage VDMOSstates+ 8
-#define VDIOcurrent VDMOSstates+ 9
-#define VDIOconduct VDMOSstates+ 10
-#define VDIOcapCharge VDMOSstates+ 11
-#define VDIOcapCurrent VDMOSstates+ 12
+#define VDIOvoltage VDMOSstates+ 9
+#define VDIOcurrent VDMOSstates+ 10
+#define VDIOconduct VDMOSstates+ 11
+#define VDIOcapCharge VDMOSstates+ 12
+#define VDIOcapCurrent VDMOSstates+ 13
 
-#define VDMOSnumStates 13
+#define VDMOScapth VDMOSstates+ 14  /* thermal capacitor value */
+#define VDMOSqth VDMOSstates+ 15  /* thermal capacitor charge */
+#define VDMOScqth VDMOSstates+ 16 /* thermal capacitor current */
+
+#define VDMOSnumStates 17
 
 
 /* per model data */
 
-    /* NOTE:  parameters marked 'input - use xxxx' are paramters for
+    /* NOTE:  parameters marked 'input - use xxxx' are parameters for
      * which a temperature correction is applied in VDMOStemp, thus
      * the VDMOSxxxx value in the per-instance structure should be used
      * instead in all calculations 
@@ -262,7 +297,7 @@ typedef struct sVDMOSmodel {       /* model structure for a resistor */
     double VDMOSqsVoltage;
     double VDMOStransconductance;    /* input - use tTransconductance */
     double VDMOSoxideCapFactor;
-    double VDMOSvt0; /* input - use tVto */
+    double VDMOSvth0; /* input - use tVth */
     double VDMOSphi; /* input - use tPhi */
     double VDMOSlambda;
     double VDMOStheta;
@@ -273,34 +308,44 @@ typedef struct sVDMOSmodel {       /* model structure for a resistor */
     double VDMOScgdmax;
     double VDMOSa;
     double VDMOScgs;
-    double VDMOSsubsl;
     double VDMOSsubshift;
     double VDMOSksubthres;
     double VDMOSmtr;
 
-    /* bulk diode */
+    /* body diode */
     double VDIOjunctionCap;   /* input - use tCj */
-    double VDIOjunctionPot;    /* input - use tBulkPot */
+    double VDIOjunctionPot;    /* input - use tJctPot */
     double VDIOdepletionCapCoeff;
     double VDIOjctSatCur;   /* input - use tSatCur */
-    double VDMOSDbv;
-    double VDMOSDibv;
+    double VDMOSbv;
+    double VDMOSibv;
     double VDIObrkdEmissionCoeff;
     double VDIOresistance;
-    double VDIOresistTemp1;
-    double VDIOresistTemp2;
-    double VDIOconductance;
     double VDMOSrds;
-    double VDMOSDn;
+    double VDMOSn;
     double VDIOtransitTime;
     double VDIOtranTimeTemp1;
     double VDIOtranTimeTemp2;
-    double VDMOSDeg;
-    double VDMOSDxti;
+    double VDMOSeg;
+    double VDMOSxti;
     double VDIOgradCoeff;
     double VDIOgradCoeffTemp1;
     double VDIOgradCoeffTemp2;
 
+    double VDMOSrthjc;
+    double VDMOSrthca;
+    double VDMOScthj;
+    double VDMOSmu;
+    double VDMOStexp0;
+    double VDMOStexp1;
+    double VDMOStcvth;
+
+    double VDMOSvgsMax;
+    double VDMOSvgdMax;
+    double VDMOSvdsMax;
+    double VDMOSvgsrMax;
+    double VDMOSvgdrMax;
+
     unsigned VDMOStypeGiven  :1;
     unsigned VDIOjctSatCurGiven :1;
     unsigned VDMOSdrainResistanceGiven   :1;
@@ -310,7 +355,7 @@ typedef struct sVDMOSmodel {       /* model structure for a resistor */
     unsigned VDMOSqsVoltageGiven    :1;
     unsigned VDMOSqsGiven    :1;
     unsigned VDMOStransconductanceGiven  :1;
-    unsigned VDMOSvt0Given   :1;
+    unsigned VDMOSvth0Given   :1;
     unsigned VDIOgradCoeffGiven    :1;
     unsigned VDIOdepletionCapCoeffGiven :1;
     unsigned VDMOSphiGiven   :1;
@@ -324,23 +369,35 @@ typedef struct sVDMOSmodel {       /* model structure for a resistor */
     unsigned VDMOScgdmaxGiven   :1;
     unsigned VDMOScgsGiven   :1;
     unsigned VDMOSaGiven   :1;
-    unsigned VDMOSsubslGiven   :1;
     unsigned VDMOSsubshiftGiven   :1;
     unsigned VDMOSksubthresGiven :1;
     unsigned VDMOSmtrGiven   :1;
 
-    unsigned VDMOSDbvGiven   :1;
-    unsigned VDMOSDibvGiven   :1;
+    unsigned VDMOSbvGiven   :1;
+    unsigned VDMOSibvGiven   :1;
     unsigned VDIOjunctionCapGiven :1;
     unsigned VDIOjunctionPotGiven :1;
     unsigned VDIObrkdEmissionCoeffGiven :1;
     unsigned VDIOresistanceGiven :1;
     unsigned VDMOSrdsGiven   :1;
-    unsigned VDMOSDnGiven   :1;
+    unsigned VDMOSnGiven   :1;
     unsigned VDIOtransitTimeGiven :1;
-    unsigned VDMOSDegGiven   :1;
-    unsigned VDMOSDxtiGiven   :1;
-
+    unsigned VDMOSegGiven   :1;
+    unsigned VDMOSxtiGiven   :1;
+
+    unsigned VDMOSrthjcGiven :1;
+    unsigned VDMOSrthcaGiven :1;
+    unsigned VDMOScthjGiven :1;
+    unsigned VDMOSmuGiven :1;
+    unsigned VDMOStexp0Given :1;
+    unsigned VDMOStexp1Given :1;
+    unsigned VDMOStcvthGiven :1;
+
+    unsigned VDMOSvgsMaxGiven  :1;
+    unsigned VDMOSvgdMaxGiven  :1;
+    unsigned VDMOSvdsMaxGiven  :1;
+    unsigned VDMOSvgsrMaxGiven  :1;
+    unsigned VDMOSvgdrMaxGiven  :1;
 } VDMOSmodel;
 
 #ifndef NMOS
@@ -362,11 +419,12 @@ enum {
     VDMOS_TEMP,
     VDMOS_M,
     VDMOS_DTEMP,
+    VDMOS_TNODEOUT,
 };
 
-/* model paramerers */
+/* model parameters */
 enum {
-    VDMOS_MOD_VTO = 101,
+    VDMOS_MOD_VTH = 101,
     VDMOS_MOD_KP,
     VDMOS_MOD_PHI,
     VDMOS_MOD_LAMBDA,
@@ -394,7 +452,6 @@ enum {
     VDMOS_MOD_CGS,
     VDMOS_MOD_RB,
     VDMOS_MOD_MTRIODE,
-    VDMOS_MOD_SUBSLOPE,
     VDMOS_MOD_SUBSHIFT,
     VDMOS_MOD_KSUBTHRES,
     VDMOS_MOD_BV,
@@ -405,6 +462,18 @@ enum {
     VDMOS_MOD_TT,
     VDMOS_MOD_EG,
     VDMOS_MOD_XTI,
+    VDMOS_MOD_RTHJC,
+    VDMOS_MOD_RTHCA,
+    VDMOS_MOD_CTHJ,
+    VDMOS_MOD_MU,
+    VDMOS_MOD_TEXP0,
+    VDMOS_MOD_TEXP1,
+    VDMOS_MOD_TCVTH,
+    VDMOS_MOD_VGS_MAX,
+    VDMOS_MOD_VGD_MAX,
+    VDMOS_MOD_VDS_MAX,
+    VDMOS_MOD_VGSR_MAX,
+    VDMOS_MOD_VGDR_MAX,
 };
 
 /* device questions */
@@ -415,14 +484,13 @@ enum {
     VDMOS_DNODE,
     VDMOS_GNODE,
     VDMOS_SNODE,
+    VDMOS_TNODE,
+    VDMOS_TCASE,
     VDMOS_DNODEPRIME,
     VDMOS_SNODEPRIME,
     VDMOS_SOURCECONDUCT,
     VDMOS_DRAINCONDUCT,
     VDMOS_VON,
-    VDMOS_VDSAT,
-    VDMOS_SOURCEVCRIT,
-    VDMOS_DRAINVCRIT,
     VDMOS_CD,
     VDMOS_GM,
     VDMOS_GDS,

src/spicelib/devices/vdmos/vdmosdset.c

@@ -21,8 +21,7 @@ VDMOSdSetup(GENmodel *inModel, CKTcircuit *ckt)
     VDMOSmodel *model = (VDMOSmodel *) inModel;
     VDMOSinstance *here;
     double Beta;
-    double DrainSatCur;
-    double SourceSatCur;
+    double OxideCap;
     double gm;
     double gds;
     double vgst;
@@ -52,11 +51,10 @@ VDMOSdSetup(GENmodel *inModel, CKTcircuit *ckt)
 
             vt = CONSTKoverQ * here->VDMOStemp;
 
-            DrainSatCur = here->VDMOSm * here->VDMOStSatCur;
-            SourceSatCur = here->VDMOSm * here->VDMOStSatCur;
+            Beta = here->VDMOStTransconductance;
 
-            Beta = here->VDMOStTransconductance * here->VDMOSm *
-                    here->VDMOSw/here->VDMOSl;
+            OxideCap = model->VDMOSoxideCapFactor * here->VDMOSl * 
+                    here->VDMOSm * here->VDMOSw;
 
             vgs = model->VDMOStype * ( 
                 *(ckt->CKTrhsOld+here->VDMOSgNode) -
@@ -84,7 +82,7 @@ VDMOSdSetup(GENmodel *inModel, CKTcircuit *ckt)
             /*
              *     this block of code evaluates the drain current and its 
              *     derivatives using the shichman-hodges model and the 
-             *     charges associated with the gate, channel and bulk for 
+             *     charges associated with the gate and channel for 
              *     mosfets
              *
              */
@@ -93,59 +91,50 @@ VDMOSdSetup(GENmodel *inModel, CKTcircuit *ckt)
              * it is obvious that they can be made global 
              */
             {
-                double betap;
-
-                vgst=(here->VDMOSmode==1?vgs:vgd);
-                vdsat=MAX(vgst,0);
-
-                if (vgst <= 0) {
-                    /*
-                     *     cutoff region
-                     */
-                /* cdrain = 0 */
-                    gm=0;
-                    gds=0;
-                    gm2=gds2=0;
-                    gmds=0;
-                    gm3=gds3=0;
-                    gm2ds=gmds2=0;
-                } else {
-                    /*
-                     *     saturation region
-                     */
-
-                    betap=Beta*(1+model->VDMOSlambda*(vds*here->VDMOSmode));
-                    /* cdrain = betap * vgst * vgst * 0.5; */
-                    if (vgst <= (vds*here->VDMOSmode)){
-                        gm=betap*vgst;
-                        gds=model->VDMOSlambda*Beta*vgst*vgst*.5;
+                double von = here->VDMOStVth * model->VDMOStype;
+                vgst = (here->VDMOSmode == 1 ? vgs : vgd) - von;
+                vdsat = MAX(vgst, 0);
+                double slope = model->VDMOSksubthres;
+                double lambda = model->VDMOSlambda;
+                double theta = model->VDMOStheta;
+                double shift = model->VDMOSsubshift;
+                double mtr = model->VDMOSmtr;
+
+                /* scale vds with mtr (except with lambda) */
+                double vdss = vds*mtr*here->VDMOSmode;
+                double t0 = 1 + lambda*vds;
+                double t1 = 1 + theta*vgs;
+                double betap = Beta*t0/t1;
+                double dbetapdvgs = -Beta*theta*t0/(t1*t1);
+                double dbetapdvds = Beta*lambda/t1;
+
+                double t2 = exp((vgst-shift)/slope);
+                vgst = slope * log(1 + t2);
+                double dvgstdvgs = t2/(t2+1);
+
+                if (vgst <= vdss) {
+                    /* saturation region */
+                    gm = betap*vgst*dvgstdvgs + 0.5*dbetapdvgs*vgst*vgst;
+                    gds = .5*dbetapdvds*vgst*vgst;
                         gm2 = betap;
                         gds2 = 0;
-                        gmds = vgst*model->VDMOSlambda*Beta;
+                        gmds = Beta*lambda*vgst;
                         gm3 = 0;
                         gds3 = 0;
-                        gm2ds = Beta * model->VDMOSlambda;
+                        gm2ds = Beta*lambda;
                         gmds2 = 0;
-
-                    } else {
-                        /*
-                         *     linear region
-                         */
-                        /* cdrain = betap * vds * (vgst - vds/2); */
-                        gm=betap*(vds*here->VDMOSmode);
-                        gds= Beta * model->VDMOSlambda*(vgst*
-                        vds*here->VDMOSmode - vds*vds*0.5) +
-                        betap*(vgst - vds*here->VDMOSmode);
+                }
+                else {
+                    /* linear region */
+                    gm = betap*vdss*dvgstdvgs + vdss*dbetapdvgs*(vgst-.5*vdss);
+                    gds = vdss*dbetapdvds*(vgst-.5*vdss) + betap*mtr*(vgst-.5*vdss) - .5*vdss*betap*mtr;
                         gm2 = 0;
-                        gds2 = 2*Beta * model->VDMOSlambda*(vgst - 
-                            vds*here->VDMOSmode) - betap;
-                        gmds = Beta * model->VDMOSlambda* vds *
-                            here->VDMOSmode + betap;
+                        gds2 = 2*Beta * lambda*(vgst - vds*here->VDMOSmode) - betap;
+                        gmds = Beta * lambda * vds * here->VDMOSmode + betap;
                         gm3=0;
-                        gds3 = -Beta*model->VDMOSlambda*3.;
+                        gds3 = -Beta*lambda*3.;
                         gm2ds=0;
-                        gmds2 = 2*model->VDMOSlambda*Beta;
-                    }
+                        gmds2 = 2*lambda*Beta;
                 }
                 /*
                  *     finished
@@ -156,15 +145,6 @@ VDMOSdSetup(GENmodel *inModel, CKTcircuit *ckt)
             /*
              *  COMPUTE EQUIVALENT DRAIN CURRENT SOURCE
              */
-                /* 
-                 * now we do the hard part of the bulk-drain and bulk-source
-                 * diode - we evaluate the non-linear capacitance and
-                 * charge
-                 *
-                 * the basic equations are not hard, but the implementation
-                 * is somewhat long in an attempt to avoid log/exponential
-                 * evaluations
-                 */
             /*
              *     meyer's capacitor model
              */
@@ -180,7 +160,6 @@ VDMOSdSetup(GENmodel *inModel, CKTcircuit *ckt)
              * expressions for the charge are available
              */
 
-     
             {
 
                 double phi;
@@ -191,7 +170,7 @@ VDMOSdSetup(GENmodel *inModel, CKTcircuit *ckt)
                 /* von, vgst and vdsat have already been adjusted for 
                     possible source-drain interchange */
                 phi = here->VDMOStPhi;
-                cox = 0;/*FIXME: can we do disto without knowing the oxide thickness?*/
+                cox = OxideCap; /*FIXME: using a guess for the oxide thickness of 1e-07 in vdmostemp.c */
                 lcapgs2=lcapgs3=lcapgd2=lcapgd3=0;
                 if (vgst <= 0) {
                     lcapgs2 = cox/(3*phi);

src/spicelib/devices/vdmos/vdmosext.h

@@ -20,3 +20,4 @@ extern int VDMOSconvTest(GENmodel*,CKTcircuit*);
 extern int VDMOSdisto(int,GENmodel*,CKTcircuit*);
 extern int VDMOSnoise(int,int,GENmodel*,CKTcircuit*,Ndata*,double*);
 extern int VDMOSdSetup(GENmodel*,CKTcircuit*);
+extern int VDMOSsoaCheck(CKTcircuit *, GENmodel *);

src/spicelib/devices/vdmos/vdmosinit.c

@@ -58,7 +58,7 @@ SPICEdev VDMOSinfo = {
     .DEVsenTrunc = NULL,
     .DEVdisto = VDMOSdisto,
     .DEVnoise = VDMOSnoise,
-    .DEVsoaCheck = NULL,
+    .DEVsoaCheck = VDMOSsoaCheck,
     .DEVinstSize = &VDMOSiSize,
     .DEVmodSize = &VDMOSmSize,
 

src/spicelib/devices/vdmos/vdmosload.c

@@ -2,7 +2,7 @@
 Copyright 1990 Regents of the University of California.  All rights reserved.
 Author: 1985 Thomas L. Quarles
 Modified: 2000 AlansFixes
-VDMOS: 2018 Holger Vogt
+VDMOS: 2018 Holger Vogt, 2020 Dietmar Warning
 **********/
 
 #include "ngspice/ngspice.h"
@@ -14,9 +14,35 @@ VDMOS: 2018 Holger Vogt
 #include "ngspice/sperror.h"
 #include "ngspice/suffix.h"
 
+/* VDMOSlimitlog(deltemp, deltemp_old, LIM_TOL, check)
+ * Logarithmic damping the per-iteration change of deltemp beyond LIM_TOL.
+ */
 static double
-cweakinv2(double sl, double shift, double vgst, double vds, double lambda, double beta, double vt, double mtr, double theta);
-
+VDMOSlimitlog(
+    double deltemp,
+    double deltemp_old,
+    double LIM_TOL,
+    int *check)
+{
+    *check = 0;
+    if (isnan (deltemp) || isnan (deltemp_old))
+    {
+        fprintf(stderr, "Alberto says:  YOU TURKEY!  The limiting function received NaN.\n");
+        fprintf(stderr, "New prediction returns to 0.0!\n");
+        deltemp = 0.0;
+        *check = 1;
+    }
+    /* Logarithmic damping of deltemp beyond LIM_TOL */
+    if (deltemp > deltemp_old + LIM_TOL) {
+        deltemp = deltemp_old + LIM_TOL + log10((deltemp-deltemp_old)/LIM_TOL);
+        *check = 1;
+    }
+    else if (deltemp < deltemp_old - LIM_TOL) {
+        deltemp = deltemp_old - LIM_TOL - log10((deltemp_old-deltemp)/LIM_TOL);
+        *check = 1;
+    }
+    return deltemp;
+}
 
 int
 VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
@@ -27,8 +53,6 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
     VDMOSmodel *model = (VDMOSmodel *)inModel;
     VDMOSinstance *here;
     double Beta;
-    double DrainSatCur;
-    double SourceSatCur;
     double arg;
     double cdhat;
     double cdrain;
@@ -55,14 +79,18 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
 #ifndef PREDICTOR
     double xfact = 0.0;
 #endif
-    int xnrm;
-    int xrev;
+    int xnrm, xrev;
     double capgs = 0.0;   /* total gate-source capacitance */
     double capgd = 0.0;   /* total gate-drain capacitance */
-    int Check;
+    double capth = 0.0;   /* total thermal capacitance */
+    int Check_mos, Check_diode;
     int error;
 
-    double CGBdummy;
+    register int selfheat;
+    double rd0T, rd1T, dBeta_dT, drd0T_dT, drd1T_dT, dIds_dT;
+    double deldelTemp, delTemp, delTemp1, Temp, Vds, Vgs;
+    double ceqqth=0.0;
+    double GmT, gTtg, gTtdp, gTtt, gTtsp, gcTt=0.0;
 
     /*  loop through all the VDMOS device models */
     for (; model != NULL; model = VDMOSnextModel(model)) {
@@ -76,34 +104,60 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
         for (here = VDMOSinstances(model); here != NULL;
                 here = VDMOSnextInstance(here)) {
 
+            selfheat = (here->VDMOStnodeoutGiven) && (model->VDMOSrthjc != 0.0);
+            if (selfheat)
+                Check_mos = 1;
+            else
+                Check_mos = 0;
+
             vt = CONSTKoverQ * here->VDMOStemp;
-            Check = 1;
 
             /* first, we compute a few useful values - these could be
              * pre-computed, but for historical reasons are still done
              * here.  They may be moved at the expense of instance size
              */
 
-            DrainSatCur = here->VDMOSm * here->VDMOStSatCur;
-            SourceSatCur = here->VDMOSm * here->VDMOStSatCur;
-            Beta = here->VDMOStTransconductance * here->VDMOSm *
-                   here->VDMOSw / here->VDMOSl;
+            delTemp = 0.0;
+            if ((ckt->CKTmode & MODEINITSMSIG)) {
+                vgs = *(ckt->CKTstate0 + here->VDMOSvgs);
+                vds = *(ckt->CKTstate0 + here->VDMOSvds);
+                delTemp = *(ckt->CKTstate0 + here->VDMOSdeltemp);
+            } else if ((ckt->CKTmode & MODEINITTRAN)) {
+                vgs = *(ckt->CKTstate1 + here->VDMOSvgs);
+                vds = *(ckt->CKTstate1 + here->VDMOSvds);
+                delTemp = *(ckt->CKTstate1 + here->VDMOSdeltemp);
+            } else if ((ckt->CKTmode & MODEINITJCT) && !here->VDMOSoff) {
+                /* ok - not one of the simple cases, so we have to
+                 * look at all of the possibilities for why we were
+                 * called.  We still just initialize the two voltages
+                 */
+                vds = model->VDMOStype * here->VDMOSicVDS;
+                vgs = model->VDMOStype * here->VDMOSicVGS;
+                delTemp = 0.0;
+                if ((vds == 0.0) && (vgs == 0.0) &&
+                        ((ckt->CKTmode & (MODETRAN | MODEAC|MODEDCOP |
+                                          MODEDCTRANCURVE)) || (!(ckt->CKTmode & MODEUIC))))
+                {
+                    vgs = model->VDMOStype * model->VDMOSvth0 + 0.1;
+                    vds = 0.0;
+                }
+            } else if ((ckt->CKTmode & (MODEINITJCT | MODEINITFIX)) && (here->VDMOSoff)) {
+                delTemp = vgs = vds = 0.0;
 
             /*
              * ok - now to do the start-up operations
              *
-             * we must get values for vbs, vds, and vgs from somewhere
+             * we must get values for vds and vgs from somewhere
              * so we either predict them or recover them from last iteration
              * These are the two most common cases - either a prediction
              * step or the general iteration step and they
              * share some code, so we put them first - others later on
              */
-
-            if ((ckt->CKTmode & (MODEINITFLOAT | MODEINITPRED | MODEINITSMSIG
-                                 | MODEINITTRAN)) ||
-                    ((ckt->CKTmode & MODEINITFIX) && (!here->VDMOSoff))) {
+            }
+            else
+            {
 #ifndef PREDICTOR
-                if (ckt->CKTmode & (MODEINITPRED | MODEINITTRAN)) {
+                if (ckt->CKTmode & MODEINITPRED) {
 
                     /* predictor step */
 
@@ -116,7 +170,13 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
                         *(ckt->CKTstate1 + here->VDMOSvds);
                     vds = (1 + xfact)* (*(ckt->CKTstate1 + here->VDMOSvds))
                           - (xfact * (*(ckt->CKTstate2 + here->VDMOSvds)));
-                } else {
+                    *(ckt->CKTstate0 + here->VDMOSdeltemp) =
+                        *(ckt->CKTstate1 + here->VDMOSdeltemp);
+                    delTemp = (1 + xfact)* (*(ckt->CKTstate1 + here->VDMOSdeltemp))
+                          - (xfact * (*(ckt->CKTstate2 + here->VDMOSdeltemp)));
+                }
+                else
+                {
 #endif /* PREDICTOR */
 
                     /* general iteration */
@@ -127,6 +187,10 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
                     vds = model->VDMOStype * (
                               *(ckt->CKTrhsOld + here->VDMOSdNodePrime) -
                               *(ckt->CKTrhsOld + here->VDMOSsNodePrime));
+                    if (selfheat)
+                        delTemp = *(ckt->CKTrhsOld + here->VDMOStempNode);
+                    else
+                        delTemp = 0.0;
 #ifndef PREDICTOR
                 }
 #endif /* PREDICTOR */
@@ -140,61 +204,73 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
                 delvds = vds - *(ckt->CKTstate0 + here->VDMOSvds);
                 delvgd = vgd - vgdo;
 
+                deldelTemp = delTemp - *(ckt->CKTstate0 + here->VDMOSdeltemp);
+
                 /* these are needed for convergence testing */
 
                 if (here->VDMOSmode >= 0) {
                     cdhat =
-                            here->VDMOScd
-                          + here->VDMOSgm * delvgs
-                          + here->VDMOSgds * delvds;
+                        here->VDMOScd
+                      + here->VDMOSgm * delvgs
+                      + here->VDMOSgds * delvds
+                      + here->VDMOSgmT * deldelTemp;
                 } else {
                     cdhat =
-                            here->VDMOScd
-                          - here->VDMOSgm * delvgd
-                          + here->VDMOSgds * delvds;
+                        here->VDMOScd
+                      - here->VDMOSgm * delvgd
+                      + here->VDMOSgds * delvds
+                      + here->VDMOSgmT * deldelTemp;
                 }
 
 #ifndef NOBYPASS
                 /* now lets see if we can bypass (ugh) */
-                if ((!(ckt->CKTmode &
-                        (MODEINITPRED | MODEINITTRAN | MODEINITSMSIG))) &&
+                if ((!(ckt->CKTmode & MODEINITPRED)) &&
                         (ckt->CKTbypass) &&
                         (fabs(delvgs) < (ckt->CKTreltol *
                                          MAX(fabs(vgs),
                                              fabs(*(ckt->CKTstate0 +
                                                     here->VDMOSvgs))) +
                                          ckt->CKTvoltTol)) &&
-                        (fabs(delvds) < (ckt->CKTreltol *
-                                         MAX(fabs(vds),
-                                             fabs(*(ckt->CKTstate0 +
-                                                    here->VDMOSvds))) +
-                                         ckt->CKTvoltTol)) &&
-                        (fabs(cdhat - here->VDMOScd) < (ckt->CKTreltol *
-                                                        MAX(fabs(cdhat),
-                                                                fabs(here->VDMOScd)) +
-                                                        ckt->CKTabstol))) {
-                    /* bypass code */
-                    /* nothing interesting has changed since last
-                     * iteration on this device, so we just
-                     * copy all the values computed last iteration out
-                     * and keep going
-                     */
-                    vgs = *(ckt->CKTstate0 + here->VDMOSvgs);
-                    vds = *(ckt->CKTstate0 + here->VDMOSvds);
-                    vgd = vgs - vds;
-                    cdrain = here->VDMOSmode * (here->VDMOScd);
-                    if (ckt->CKTmode & (MODETRAN | MODETRANOP)) {
-                        capgs = (*(ckt->CKTstate0 + here->VDMOScapgs) +
-                                 *(ckt->CKTstate1 + here->VDMOScapgs));
-                        capgd = (*(ckt->CKTstate0 + here->VDMOScapgd) +
-                                 *(ckt->CKTstate1 + here->VDMOScapgd));
-
-                    }
-                    goto bypass;
-                }
+                            (fabs(delvds) < (ckt->CKTreltol *
+                                             MAX(fabs(vds),
+                                                 fabs(*(ckt->CKTstate0 +
+                                                        here->VDMOSvds))) +
+                                             ckt->CKTvoltTol)) &&
+                                (fabs(cdhat - here->VDMOScd) < (ckt->CKTreltol *
+                                                                MAX(fabs(cdhat),
+                                                                        fabs(here->VDMOScd)) +
+                                                                ckt->CKTabstol)) &&
+                                    ((here->VDMOStempNode == 0) ||
+                                            (fabs(deldelTemp) < (ckt->CKTreltol * MAX(fabs(delTemp),
+                                                                                      fabs(*(ckt->CKTstate0+here->VDMOSdeltemp)))
+                                                                 + ckt->CKTvoltTol*1e4))))
+                                        {
+                                            /* bypass code */
+                                            /* nothing interesting has changed since last
+                                             * iteration on this device, so we just
+                                             * copy all the values computed last iteration out
+                                             * and keep going
+                                             */
+                                            vgs = *(ckt->CKTstate0 + here->VDMOSvgs);
+                                            vds = *(ckt->CKTstate0 + here->VDMOSvds);
+                                            vgd = vgs - vds;
+                                            delTemp = *(ckt->CKTstate0 + here->VDMOSdeltemp);
+                                            /*  calculate Vds for temperature conductance calculation
+                                                in bypass (used later when filling Temp node matrix)  */
+                                            Vds = here->VDMOSmode > 0 ? vds : -vds;
+                                            cdrain = here->VDMOSmode * (here->VDMOScd);
+                                            if (ckt->CKTmode & (MODETRAN | MODETRANOP)) {
+                                                capgs = (*(ckt->CKTstate0 + here->VDMOScapgs) +
+                                                         *(ckt->CKTstate1 + here->VDMOScapgs));
+                                                capgd = (*(ckt->CKTstate0 + here->VDMOScapgd) +
+                                                         *(ckt->CKTstate1 + here->VDMOScapgd));
+                                                capth = (*(ckt->CKTstate0 + here->VDMOScapth) +
+                                                         *(ckt->CKTstate1 + here->VDMOScapth));
+                                            }
+                                            goto bypass;
+                                        }
 #endif /*NOBYPASS*/
 
-
                 /* ok - bypass is out, do it the hard way */
 
                 von = model->VDMOStype * here->VDMOSvon;
@@ -202,7 +278,7 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
 #ifndef NODELIMITING
                 /*
                  * limiting
-                 *  we want to keep device voltages from changing
+                 * we want to keep device voltages from changing
                  * so fast that the exponentials churn out overflows
                  * and similar rudeness
                  */
@@ -222,32 +298,42 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
                     }
                     vgs = vgd + vds;
                 }
+                if (selfheat)
+                    delTemp = VDMOSlimitlog(delTemp,
+                          *(ckt->CKTstate0 + here->VDMOSdeltemp),100,&Check_mos);
+                else
+                    delTemp = 0.0;
 #endif /*NODELIMITING*/
 
+            }
 
-            } else {
-
-                /* ok - not one of the simple cases, so we have to
-                 * look at all of the possibilities for why we were
-                 * called.  We still just initialize the three voltages
-                 */
-
-                if ((ckt->CKTmode & MODEINITJCT) && !here->VDMOSoff) {
-                    vds = model->VDMOStype * here->VDMOSicVDS;
-                    vgs = model->VDMOStype * here->VDMOSicVGS;
-                    if ((vds == 0) && (vgs == 0) &&
-                            ((ckt->CKTmode &
-                              (MODETRAN | MODEDCOP | MODEDCTRANCURVE)) ||
-                             (!(ckt->CKTmode & MODEUIC)))) {
-                        vgs = model->VDMOStype * here->VDMOStVto;
-                        vds = 0;
-                    }
-                } else {
-                    vgs = vds = 0;
+            Temp = delTemp + here->VDMOStemp;
+            here->VDMOSTempSH = Temp; /* added for portability of SH Temp for noise analysis */
+
+            /*  Calculate temperature dependent values for self-heating effect  */
+            if (selfheat) {
+                double TempRatio = Temp / here->VDMOStemp;
+                Beta = here->VDMOStTransconductance * pow(TempRatio,model->VDMOSmu);
+                dBeta_dT = here->VDMOStTransconductance * model->VDMOSmu / (here->VDMOStemp * pow(TempRatio,1-model->VDMOSmu));
+                rd0T =  here->VDMOSdrainResistance * pow(TempRatio, model->VDMOStexp0);
+                drd0T_dT = rd0T * model->VDMOStexp0 / Temp;
+                rd1T = 0.0;
+                drd1T_dT = 0.0;
+                if (model->VDMOSqsGiven) {
+                    rd1T = here->VDMOSqsResistance * pow(TempRatio, model->VDMOStexp1);
+                    drd1T_dT = rd1T * model->VDMOStexp1 / Temp;
                 }
+            } else {
+                Beta = here->VDMOStTransconductance;
+                dBeta_dT = 0.0;
+                rd0T = here->VDMOSdrainResistance;
+                drd0T_dT = 0.0;
+                rd1T = 0.0;
+                if (model->VDMOSqsGiven)
+                    rd1T = here->VDMOSqsResistance;
+                drd1T_dT = 0.0;
             }
 
-
             /*
              * now all the preliminaries are over - we can start doing the
              * real work
@@ -261,125 +347,63 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
             if (vds >= 0) {
                 /* normal mode */
                 here->VDMOSmode = 1;
+                Vds = vds;
+                Vgs = vgs;
             } else {
                 /* inverse mode */
                 here->VDMOSmode = -1;
+                Vds = -vds;
+                Vgs = vgd;
             }
 
             {
                 /*
                  *     this block of code evaluates the drain current and its
                  *     derivatives using the shichman-hodges model and the
-                 *     charges associated with the gate, channel and bulk for
+                 *     charges associated with the gate and channel for
                  *     mosfets
                  *
                  */
 
-                /* the following 2 variables are local to this code block until
-                 * it is obvious that they can be made global
-                 */
-                double betap;
-                double vgst;
-
-                von = (model->VDMOSvt0*model->VDMOStype);
-                vgst = (here->VDMOSmode == 1 ? vgs : vgd) - von;
+                von = here->VDMOStVth * model->VDMOStype;
+                double vgst = (here->VDMOSmode == 1 ? vgs : vgd) - von;
                 vdsat = MAX(vgst, 0);
-                if (model->VDMOSksubthresGiven) {
-                /* Alternative simple weak inversion model, according to https://www.anasoft.co.uk/MOS1Model.htm
+                /* Simple weak inversion model, according to https://www.anasoft.co.uk/MOS1Model.htm
                  * Scale the voltage overdrive vgst logarithmically in weak inversion.
                  * Best fits LTSPICE curves with shift=0
-                 * Drain current including subthreshold current */
-
-                    double slope = model->VDMOSksubthres;
-                    double lambda = model->VDMOSlambda;
-                    double theta = model->VDMOStheta;
-                    double shift = model->VDMOSsubshift;
-                    double mtr = model->VDMOSmtr;
-
-                    /* scale vds with mtr (except with lambda) */
-                    double vdss = vds*mtr*here->VDMOSmode;
-                    double t0 = 1 + lambda*vds;
-                    double t1 = 1 + theta*vgs;
-                    betap = Beta*t0/t1;
-                    double dbetapdvgs = -Beta*theta*t0/(t1*t1);
-                    double dbetapdvds = Beta*lambda/t1;
-
-                    double t2 = exp((vgst-shift)/slope);
-                    vgst = slope * log(1 + t2);
-                    double dvgstdvgs = t2/(t2+1);
-
-                    if (vgst <= vdss) {
-                        /* saturation region */
-                        cdrain = betap*vgst*vgst*.5;
-                        here->VDMOSgm = betap*vgst*dvgstdvgs + 0.5*dbetapdvgs*vgst*vgst;
-                        here->VDMOSgds = .5*dbetapdvds*vgst*vgst;
-                    }
-                    else {
-                        /* linear region */
-                        cdrain = betap * vdss * (vgst - .5 * vdss);
-                        here->VDMOSgm = betap*vdss*dvgstdvgs + vdss*dbetapdvgs*(vgst-.5*vdss);
-                        here->VDMOSgds = vdss*dbetapdvds*(vgst-.5*vdss) + betap*mtr*(vgst-.5*vdss) - .5*vdss*betap*mtr;
-                    }
+                 */
+                double slope = model->VDMOSksubthres;
+                double lambda = model->VDMOSlambda;
+                double theta = model->VDMOStheta;
+                double shift = model->VDMOSsubshift;
+                double mtr = model->VDMOSmtr;
+
+                /* scale vds with mtr (except with lambda) */
+                double vdss = vds*mtr*here->VDMOSmode;
+                double t0 = 1 + lambda*vds;
+                double t1 = 1 + theta*vgs;
+                double betap = Beta*t0/t1;
+                double dbetapdvgs = -Beta*theta*t0/(t1*t1);
+                double dbetapdvds = Beta*lambda/t1;
+                double dbetapdT = dBeta_dT*t0/t1;
+
+                double t2 = exp((vgst-shift)/slope);
+                vgst = slope * log(1 + t2);
+                double dvgstdvgs = t2/(t2+1);
+
+                if (vgst <= vdss) {
+                    /* saturation region */
+                    cdrain = betap * vgst*vgst * .5;
+                    here->VDMOSgm = betap*vgst*dvgstdvgs + 0.5*dbetapdvgs*vgst*vgst;
+                    here->VDMOSgds = .5*dbetapdvds*vgst*vgst;
+                    dIds_dT = dbetapdT * vgst*vgst * .5;
                 }
-                else if (model->VDMOSsubslGiven) {
-                    /* numerical differentiation for gd and gm with a delta of 2 mV */
-                    double vdsm = vds * here->VDMOSmode;
-                    double delta = 0.001;
-                    cdrain = cweakinv2(model->VDMOSsubsl, model->VDMOSsubshift, vgst, vdsm, model->VDMOSlambda,
-                        Beta, vt, model->VDMOSmtr, model->VDMOStheta);
-                    /* gd */
-                    double vds1 = vdsm + delta;
-                    double cdrp = cweakinv2(model->VDMOSsubsl, model->VDMOSsubshift, vgst, vds1, model->VDMOSlambda,
-                        Beta, vt, model->VDMOSmtr, model->VDMOStheta);
-                    vds1 = vdsm - delta;
-                    double cdrm = cweakinv2(model->VDMOSsubsl, model->VDMOSsubshift, vgst, vds1, model->VDMOSlambda,
-                        Beta, vt, model->VDMOSmtr, model->VDMOStheta);
-                    here->VDMOSgds = (cdrp - cdrm) / (2. * delta);
-                    /* gm */
-                    double vgst1 = vgst + delta;
-                    cdrp = cweakinv2(model->VDMOSsubsl, model->VDMOSsubshift, vgst1, vdsm, model->VDMOSlambda,
-                        Beta, vt, model->VDMOSmtr, model->VDMOStheta);
-                    vgst1 = vgst - delta;
-                    cdrm = cweakinv2(model->VDMOSsubsl, model->VDMOSsubshift, vgst1, vdsm, model->VDMOSlambda,
-                        Beta, vt, model->VDMOSmtr, model->VDMOStheta);
-                    here->VDMOSgm = (cdrp - cdrm) / (2. * delta);
-                } else {
-                    double onfg, fgate, Betam, dfgdvg;
-                    onfg = 1.0+model->VDMOStheta*vgst;
-                    fgate = 1.0/onfg;
-                    Betam = Beta * fgate;
-                    dfgdvg = -model->VDMOStheta*fgate*fgate;
-                    if (vgst <= 0) {
-                        /*
-                         *     cutoff region
-                         */
-                        cdrain = 0;
-                        here->VDMOSgm = 0;
-                        here->VDMOSgds = 0;
-                    } else {
-                        /* scale vds with mtr */
-                        double mtr = model->VDMOSmtr;
-                        betap = Betam*(1 + model->VDMOSlambda*(vds*here->VDMOSmode));
-                        if (vgst <= (vds * here->VDMOSmode) * mtr) {
-                            /*
-                             *     saturation region
-                             */
-                            cdrain = betap*vgst*vgst*.5;
-                            here->VDMOSgm = betap*vgst * fgate + dfgdvg * cdrain;
-                            here->VDMOSgds = model->VDMOSlambda*Betam*vgst*vgst*.5;
-                        } else {
-                            /*
-                             *     linear region
-                             */
-                            cdrain = betap * (vds * here->VDMOSmode) * mtr *
-                                     (vgst - .5 * (vds*here->VDMOSmode) * mtr);
-                            here->VDMOSgm = betap * (vds * here->VDMOSmode) * mtr  * fgate + dfgdvg * cdrain;
-                            here->VDMOSgds = betap * (vgst - (vds * here->VDMOSmode) * mtr) +
-                                             model->VDMOSlambda * Betam *
-                                             (vds * here->VDMOSmode) * mtr *
-                                             (vgst - .5 * (vds * here->VDMOSmode) * mtr);
-                        }
-                    }
+                else {
+                    /* linear region */
+                    cdrain = betap * vdss * (vgst - .5 * vdss);
+                    here->VDMOSgm = betap*vdss*dvgstdvgs + vdss*dbetapdvgs*(vgst-.5*vdss);
+                    here->VDMOSgds = vdss*dbetapdvds*(vgst-.5*vdss) + betap*mtr*(vgst-.5*vdss) - .5*vdss*betap*mtr;
+                    dIds_dT = dbetapdT * vdss * (vgst - .5 * vdss);
                 }
             }
 
@@ -398,7 +422,47 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
 
             *(ckt->CKTstate0 + here->VDMOSvgs) = vgs;
             *(ckt->CKTstate0 + here->VDMOSvds) = vds;
+            *(ckt->CKTstate0 + here->VDMOSdeltemp) = delTemp;
 
+            /* quasi saturation
+             * according to Vincenzo d'Alessandro's Quasi-Saturation Model, simplified:
+             V. D'Alessandro, F. Frisina, N. Rinaldi: A New SPICE Model of VDMOS Transistors
+             Including Thermal and Quasi-saturation Effects, 9th European Conference on Power
+             Electronics and applications (EPE), Graz, Austria, August 2001, pp. P.1 − P.10.
+             */
+            if (model->VDMOSqsGiven && (here->VDMOSmode == 1)) {
+                double vdsn = model->VDMOStype * (
+                    *(ckt->CKTrhsOld + here->VDMOSdNode) -
+                    *(ckt->CKTrhsOld + here->VDMOSsNode));
+                double rd = rd0T + rd1T * (vdsn / (vdsn + fabs(model->VDMOSqsVoltage)));
+                if (rd > 0)
+                    here->VDMOSdrainConductance = 1 / rd + ckt->CKTgmin;
+                else
+                    here->VDMOSdrainConductance = 1 / rd0T;
+            } else {
+                if (rd0T > 0)
+                    here->VDMOSdrainConductance = 1 / rd0T;
+            }
+
+            if (selfheat) {
+                GmT = dIds_dT;
+                here->VDMOSgmT = GmT;
+            } else {
+                GmT = 0.0;
+                here->VDMOSgmT = 0.0;
+            }
+
+            if (selfheat) {
+                /*  note that sign is switched because power flows out
+                    of device into the temperature node. */
+                here->VDMOSgtempg = -model->VDMOStype*here->VDMOSgm * Vds;
+                here->VDMOSgtempT = -GmT * Vds;
+                here->VDMOSgtempd = -model->VDMOStype* (here->VDMOSgds * Vds + cdrain);
+                here->VDMOScth = - cdrain * Vds
+                                 - 1/here->VDMOSdrainConductance * cdrain*cdrain
+                                 - model->VDMOStype * (here->VDMOSgtempg * Vgs + here->VDMOSgtempd * Vds)
+                                 - here->VDMOSgtempT * delTemp;
+            }
 
             /*
              * vdmos capacitor model
@@ -406,7 +470,7 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
             if (ckt->CKTmode & (MODETRAN | MODETRANOP | MODEINITSMSIG)) {
                 /*
                  * calculate gate - drain, gate - source capacitors
-                 * drain-source capacitor is evaluated with the bulk diode below
+                 * drain-source capacitor is evaluated with the body diode below
                  */
                 /*
                  * this just evaluates at the current time,
@@ -417,23 +481,24 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
                  */
                 DevCapVDMOS(vgd, cgdmin, cgdmax, a, cgs,
                             (ckt->CKTstate0 + here->VDMOScapgs),
-                            (ckt->CKTstate0 + here->VDMOScapgd),
-                             &CGBdummy);
+                            (ckt->CKTstate0 + here->VDMOScapgd));
+                *(ckt->CKTstate0 + here->VDMOScapth) = model->VDMOScthj; /* always constant */
 
                 vgs1 = *(ckt->CKTstate1 + here->VDMOSvgs);
                 vgd1 = vgs1 - *(ckt->CKTstate1 + here->VDMOSvds);
+                delTemp1 = *(ckt->CKTstate1 + here->VDMOSdeltemp);
                 if (ckt->CKTmode & (MODETRANOP | MODEINITSMSIG)) {
                     capgs = 2 * *(ckt->CKTstate0 + here->VDMOScapgs);
                     capgd = 2 * *(ckt->CKTstate0 + here->VDMOScapgd);
+                    capth = 2 * *(ckt->CKTstate0 + here->VDMOScapth);
                 } else {
                     capgs = (*(ckt->CKTstate0 + here->VDMOScapgs) +
                              *(ckt->CKTstate1 + here->VDMOScapgs));
                     capgd = (*(ckt->CKTstate0 + here->VDMOScapgd) +
                              *(ckt->CKTstate1 + here->VDMOScapgd));
+                    capth = (*(ckt->CKTstate0 + here->VDMOScapth) +
+                             *(ckt->CKTstate1 + here->VDMOScapth));
                 }
-                /*
-
-                */
 
 #ifndef PREDICTOR
                 if (ckt->CKTmode & (MODEINITPRED | MODEINITTRAN)) {
@@ -443,6 +508,9 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
                     *(ckt->CKTstate0 + here->VDMOSqgd) =
                         (1 + xfact) * *(ckt->CKTstate1 + here->VDMOSqgd)
                         - xfact * *(ckt->CKTstate2 + here->VDMOSqgd);
+                    *(ckt->CKTstate0 + here->VDMOSqth) =
+                        (1 + xfact) * *(ckt->CKTstate1 + here->VDMOSqth)
+                        - xfact * *(ckt->CKTstate2 + here->VDMOSqth);
                 } else {
 #endif /*PREDICTOR*/
                     if (ckt->CKTmode & MODETRAN) {
@@ -450,19 +518,21 @@ VDMOSload(GENmodel *inModel, CKTcircuit *ckt)
                                                              *(ckt->CKTstate1 + here->VDMOSqgs);
                         *(ckt->CKTstate0 + here->VDMOSqgd) = (vgd - vgd1)*capgd +
                                                              *(ckt->CKTstate1 + here->VDMOSqgd);
+                        *(ckt->CKTstate0 + here->VDMOSqth) = (delTemp-delTemp1)*capth +
+                                                             *(ckt->CKTstate1 + here->VDMOSqth);
                     } else {
                         /* TRANOP only */
                         *(ckt->CKTstate0 + here->VDMOSqgs) = vgs*capgs;
                         *(ckt->CKTstate0 + here->VDMOSqgd) = vgd*capgd;
+                        *(ckt->CKTstate0 + here->VDMOSqth) = delTemp*capth;
                     }
 #ifndef PREDICTOR
                 }
 #endif /*PREDICTOR*/
             }
 #ifndef NOBYPASS
-bypass :
+bypass:
 #endif
-
             if ((ckt->CKTmode & (MODEINITTRAN)) ||
                     (!(ckt->CKTmode & (MODETRAN)))) {
                 /*
@@ -473,12 +543,15 @@ bypass :
                 ceqgs = 0;
                 gcgd = 0;
                 ceqgd = 0;
+                gcTt = 0.0;
+                ceqqth = 0.0;
             } else {
                 if (capgs == 0) *(ckt->CKTstate0 + here->VDMOScqgs) = 0;
                 if (capgd == 0) *(ckt->CKTstate0 + here->VDMOScqgd) = 0;
+                if (capth == 0) *(ckt->CKTstate0 + here->VDMOScqth) = 0;
                 /*
                  *    calculate equivalent conductances and currents for
-                 *    meyer"s capacitors
+                 *    vdmos capacitors
                  */
                 error = NIintegrate(ckt, &gcgs, &ceqgs, capgs, here->VDMOSqgs);
                 if (error) return(error);
@@ -488,51 +561,67 @@ bypass :
                         *(ckt->CKTstate0 + here->VDMOSqgs);
                 ceqgd = ceqgd - gcgd*vgd + ckt->CKTag[0] *
                         *(ckt->CKTstate0 + here->VDMOSqgd);
+                if (selfheat)
+                {
+                    error = NIintegrate(ckt, &gcTt, &ceqqth, 0.0, here->VDMOSqth);
+                    gcTt = model->VDMOScthj * ckt->CKTag[0];
+                    ceqqth = *(ckt->CKTstate0 + here->VDMOScqth) - gcTt * delTemp;
+                }
             }
 
             /*
              *  load current vector
              */
+
+            if (selfheat) {
+                if (here->VDMOSmode >= 0) {
+                    GmT = model->VDMOStype * here->VDMOSgmT;
+                    gTtg  = here->VDMOSgtempg;
+                    gTtdp = here->VDMOSgtempd;
+                    gTtt  = here->VDMOSgtempT;
+                    gTtsp = - (gTtg + gTtdp);
+                } else {
+                    GmT = -model->VDMOStype * here->VDMOSgmT;
+                    gTtg  = here->VDMOSgtempg;
+                    gTtsp = here->VDMOSgtempd;
+                    gTtt  = here->VDMOSgtempT;
+                    gTtdp = - (gTtg + gTtsp);
+                }
+            } else {
+                GmT = 0.0;
+                gTtg  = 0.0;
+                gTtdp = 0.0;
+                gTtt  = 0.0;
+                gTtsp = 0.0;
+            }
+
             if (here->VDMOSmode >= 0) {
                 xnrm = 1;
                 xrev = 0;
-                cdreq = model->VDMOStype*(cdrain - here->VDMOSgds*vds -
-                                          here->VDMOSgm*vgs);
+                cdreq = model->VDMOStype*(cdrain - here->VDMOSgds*vds
+                                                 - here->VDMOSgm*vgs)
+                        - GmT * delTemp;
             } else {
                 xnrm = 0;
                 xrev = 1;
-                cdreq = -(model->VDMOStype)*(cdrain - here->VDMOSgds*(-vds) -
-                                             here->VDMOSgm*vgd);
+                cdreq = -(model->VDMOStype)*(cdrain - here->VDMOSgds*(-vds)
+                                                    - here->VDMOSgm*vgd)
+                        - GmT * delTemp;
             }
-            *(ckt->CKTrhs + here->VDMOSgNodePrime) -=
-                (model->VDMOStype * (ceqgs + ceqgd));
-            *(ckt->CKTrhs + here->VDMOSdNodePrime) +=
-                (-cdreq + model->VDMOStype * ceqgd);
-            *(ckt->CKTrhs + here->VDMOSsNodePrime) +=
-                cdreq + model->VDMOStype * ceqgs;
 
-
-            /* quasi saturation
-             * according to Vincenzo d'Alessandro's Quasi-Saturation Model, simplified:
-             V. D'Alessandro, F. Frisina, N. Rinaldi: A New SPICE Model of VDMOS Transistors
-             Including Thermal and Quasi-saturation Effects, 9th European Conference on Power
-             Electronics and applications (EPE), Graz, Austria, August 2001, pp. P.1 − P.10.
-             */
-            if (model->VDMOSqsGiven && (here->VDMOSmode == 1)) {
-                double vdsn = model->VDMOStype * (
-                    *(ckt->CKTrhsOld + here->VDMOSdNode) -
-                    *(ckt->CKTrhsOld + here->VDMOSsNode));
-                double rd = model->VDMOSdrainResistance + model->VDMOSqsResistance *
-                    (vdsn / (vdsn + fabs(model->VDMOSqsVoltage)));
-                here->VDMOSdrainConductance = 1 / rd;
+            *(ckt->CKTrhs + here->VDMOSgNodePrime) -= (model->VDMOStype * (ceqgs + ceqgd));
+            *(ckt->CKTrhs + here->VDMOSdNodePrime) += (-cdreq + model->VDMOStype * ceqgd);
+            *(ckt->CKTrhs + here->VDMOSsNodePrime) +=   cdreq + model->VDMOStype * ceqgs;
+            if (selfheat) {
+                *(ckt->CKTrhs + here->VDMOStempNode) -= here->VDMOScth + ceqqth; /* dissipated power + Cthj current */
+                *(ckt->CKTrhs + here->VDMOSvcktTbranch) = ckt->CKTtemp-CONSTCtoK; /* ckt temperature */
             }
 
-
             /*
              *  load y matrix
              */
             *(here->VDMOSDdPtr) += (here->VDMOSdrainConductance + here->VDMOSdsConductance);
-            *(here->VDMOSGgPtr) += (here->VDMOSgateConductance); //((gcgd + gcgs + gcgb));
+            *(here->VDMOSGgPtr) += (here->VDMOSgateConductance);
             *(here->VDMOSSsPtr) += (here->VDMOSsourceConductance + here->VDMOSdsConductance);
             *(here->VDMOSDPdpPtr) +=
                 (here->VDMOSdrainConductance + here->VDMOSgds +
@@ -560,26 +649,40 @@ bypass :
             *(here->VDMOSDsPtr) += (-here->VDMOSdsConductance);
             *(here->VDMOSSdPtr) += (-here->VDMOSdsConductance);
 
+            if (selfheat)
+            {
+                (*(here->VDMOSDPtempPtr)     += GmT);
+                (*(here->VDMOSSPtempPtr)     += -GmT);
+                (*(here->VDMOSGPtempPtr)     += 0.0);
+                (*(here->VDMOSTemptempPtr)   += gTtt + 1/model->VDMOSrthjc + gcTt);
+                (*(here->VDMOSTempgpPtr)     += gTtg);
+                (*(here->VDMOSTempdpPtr)     += gTtdp);
+                (*(here->VDMOSTempspPtr)     += gTtsp);
+                (*(here->VDMOSTemptcasePtr)  += -1/model->VDMOSrthjc);
+                (*(here->VDMOSTcasetempPtr)  += -1/model->VDMOSrthjc);
+                (*(here->VDMOSTcasetcasePtr) +=  1/model->VDMOSrthjc + 1/model->VDMOSrthca);
+                (*(here->VDMOSTptpPtr)       +=  1/model->VDMOSrthca);
+                (*(here->VDMOSTptcasePtr)    += -1/model->VDMOSrthca);
+                (*(here->VDMOSTcasetpPtr)    += -1/model->VDMOSrthca);
+                (*(here->VDMOSCktTtpPtr)     +=  1.0);
+                (*(here->VDMOSTpcktTPtr)     +=  1.0);
+            }
 
-            /* bulk diode model
+            /* body diode model
              * Delivers reverse conduction and forward breakdown
              * of VDMOS transistor
              */
 
             double vd;      /* current diode voltage */
-            double vdtemp;
             double vte;
             double vtebrk, vbrknp;
             double cd, cdb, csat, cdeq;
-            double czero;
-            double czof2;
             double capd;
             double gd, gdb, gspr;
             double delvd;   /* change in diode voltage temporary */
-            double diffcharge, deplcharge, diffcap, deplcap;
-            double evd, evrev;
+            double evrev;
 #ifndef NOBYPASS
-            double tol;     /* temporary for tolerence calculations */
+            double tol;     /* temporary for tolerance calculations */
 #endif
 
             cd = 0.0;
@@ -588,11 +691,11 @@ bypass :
             gdb = 0.0;
             csat = here->VDIOtSatCur;
             gspr = here->VDIOtConductance;
-            vte = model->VDMOSDn * vt;
+            vte = model->VDMOSn * vt;
             vtebrk = model->VDIObrkdEmissionCoeff * vt;
             vbrknp = here->VDIOtBrkdwnV;
 
-            Check = 1;
+            Check_diode = 1;
             if (ckt->CKTmode & MODEINITSMSIG) {
                 vd = *(ckt->CKTstate0 + here->VDIOvoltage);
             } else if (ckt->CKTmode & MODEINITTRAN) {
@@ -646,29 +749,31 @@ bypass :
                 /*
                 *   limit new junction voltage
                 */
-                if ((model->VDMOSDbvGiven) &&
+                if ((model->VDMOSbvGiven) &&
                         (vd < MIN(0, -vbrknp + 10 * vtebrk))) {
+                    double vdtemp;
                     vdtemp = -(vd + vbrknp);
                     vdtemp = DEVpnjlim(vdtemp,
                                        -(*(ckt->CKTstate0 + here->VDIOvoltage) +
                                          vbrknp), vtebrk,
-                                       here->VDIOtVcrit, &Check);
+                                       here->VDIOtVcrit, &Check_diode);
                     vd = -(vdtemp + vbrknp);
                 } else {
                     vd = DEVpnjlim(vd, *(ckt->CKTstate0 + here->VDIOvoltage),
-                                   vte, here->VDIOtVcrit, &Check);
+                                   vte, here->VDIOtVcrit, &Check_diode);
                 }
             }
             /*
             *   compute dc current and derivatives
             */
             if (vd >= -3 * vte) {                 /* bottom current forward */
+                double evd;
 
                 evd = exp(vd / vte);
                 cdb = csat*(evd - 1);
                 gdb = csat*evd / vte;
 
-            } else if ((!(model->VDMOSDbvGiven)) ||
+            } else if ((!(model->VDMOSbvGiven)) ||
                        vd >= -vbrknp) { /* reverse */
 
                 arg = 3 * vte / (vd*CONSTe);
@@ -696,6 +801,7 @@ bypass :
                 /*
                 *   charge storage elements
                 */
+                double czero, czof2, diffcharge, deplcharge, diffcap, deplcap;
                 czero = here->VDIOtJctCap;
                 if (vd < here->VDIOtDepCap) {
                     arg = 1 - vd / here->VDIOtJctPot;
@@ -751,7 +857,7 @@ bypass :
             *   check convergence
             */
 
-            if (Check == 1) {
+            if ((Check_mos == 1) || (Check_diode == 1)) {
                 ckt->CKTnoncon++;
                 ckt->CKTtroubleElt = (GENinstance *)here;
             }
@@ -761,7 +867,7 @@ bypass :
             *(ckt->CKTstate0 + here->VDIOconduct) = gd;
 
 #ifndef NOBYPASS
-load :
+load:
 #endif
             /*
             *   load current vector
@@ -789,47 +895,3 @@ load :
     }
     return(OK);
 }
-
-
-/* scaling function, sine function interpolating between 0 and 1
- * nf2: empirical setting of sine 'speed'
- */
-
-static double
-scalef(double nf2, double vgst)
-{
-    double vgstsin = vgst / nf2;
-    if (vgstsin > 1)
-        return 1;
-    else if (vgstsin < -1)
-        return 0;
-    else
-        return 0.5 * sin(vgstsin * M_PI / 2) + 0.5;
-}
-
-
-/* Calculate D/S current including weak inversion.
- * Uses a single function covering weak-moderate-stong inversion, as well
- * as linear and saturation regions, with an interpolation method according to
- * Tvividis, McAndrew: "Operation and Modeling of the MOS Transistor", Oxford, 2011, p. 209.
- * A single parameter n sets the slope of the weak inversion current. The weak inversion
- * current is independent from vds, as in long channel devices.
- * The following modification has been added for VDMOS compatibility:
- * n and lambda are depending on vgst with a sine function interpolating between 0 and 1.
- */
-
-static double
-cweakinv2(double slope, double shift, double vgst, double vds, double lambda, double beta, double vt, double mtr, double theta)
-{
-    double betam = beta / (1.0+theta*vgst);
-    vgst += shift * (1 - scalef(0.5, vgst));
-    double n = slope / 2.3 / 0.0256; /* Tsividis, p. 208 */
-    double n1 = n + (1 - n) * scalef(0.7, vgst); /* n < n1 < 1 */
-    double first = log(1 + exp(vgst / (2 * n1 * vt)));
-    double second = log(1 + exp((vgst - vds * mtr * n1) / (2 * n1 * vt)));
-    double cds =
-        betam * n1 * 2 * vt * vt * (1 + scalef(1, vgst) * lambda * vds) *
-        (first * first - second * second);
-    return cds;
-}
-

src/spicelib/devices/vdmos/vdmosmask.c

@@ -24,8 +24,8 @@ VDMOSmAsk(CKTcircuit *ckt, GENmodel *inst, int which, IFvalue *value)
         case VDMOS_MOD_TNOM:
             value->rValue = model->VDMOStnom-CONSTCtoK;
             return(OK);
-        case VDMOS_MOD_VTO:
-            value->rValue = model->VDMOSvt0;
+        case VDMOS_MOD_VTH:
+            value->rValue = model->VDMOSvth0;
             return(OK);
         case VDMOS_MOD_KP:
             value->rValue = model->VDMOStransconductance;
@@ -57,9 +57,6 @@ VDMOSmAsk(CKTcircuit *ckt, GENmodel *inst, int which, IFvalue *value)
         case VDMOS_MOD_MTRIODE:
             value->rValue = model->VDMOSmtr;
             return(OK);
-        case VDMOS_MOD_SUBSLOPE:
-            value->rValue = model->VDMOSsubsl;
-            return(OK);
         case VDMOS_MOD_SUBSHIFT:
             value->rValue = model->VDMOSsubshift;
             return(OK);
@@ -93,7 +90,7 @@ VDMOSmAsk(CKTcircuit *ckt, GENmodel *inst, int which, IFvalue *value)
             value->rValue = model->VDIOjctSatCur;
             return(OK);
         case VDMOS_MOD_N:
-            value->rValue = model->VDMOSDn;
+            value->rValue = model->VDMOSn;
             return(OK);
         case VDMOS_MOD_VJ:
             value->rValue = model->VDIOjunctionPot;
@@ -105,10 +102,10 @@ VDMOSmAsk(CKTcircuit *ckt, GENmodel *inst, int which, IFvalue *value)
             value->rValue = model->VDIOgradCoeff;
             return(OK);
         case VDMOS_MOD_BV:
-            value->rValue = model->VDMOSDbv;
+            value->rValue = model->VDMOSbv;
             return(OK);
         case VDMOS_MOD_IBV:
-            value->rValue = model->VDMOSDibv;
+            value->rValue = model->VDMOSibv;
             return(OK);
         case VDMOS_MOD_NBV:
             value->rValue = model->VDIObrkdEmissionCoeff;
@@ -123,10 +120,46 @@ VDMOSmAsk(CKTcircuit *ckt, GENmodel *inst, int which, IFvalue *value)
             value->rValue = model->VDIOtransitTime;
             return(OK);
         case VDMOS_MOD_EG:
-            value->rValue = model->VDMOSDeg;
+            value->rValue = model->VDMOSeg;
             return(OK);
         case VDMOS_MOD_XTI:
-            value->rValue = model->VDMOSDxti;
+            value->rValue = model->VDMOSxti;
+            return(OK);
+        case VDMOS_MOD_RTHJC:
+            value->rValue = model->VDMOSrthjc; 
+            return(OK);
+        case VDMOS_MOD_RTHCA:
+            value->rValue = model->VDMOSrthca; 
+            return(OK);
+        case VDMOS_MOD_CTHJ:
+            value->rValue = model->VDMOScthj; 
+            return(OK);
+        case VDMOS_MOD_MU:
+            value->rValue = model->VDMOSmu; 
+            return(OK);
+        case VDMOS_MOD_TEXP0:
+            value->rValue = model->VDMOStexp0; 
+            return(OK);
+        case VDMOS_MOD_TEXP1:
+            value->rValue = model->VDMOStexp1; 
+            return(OK);
+        case VDMOS_MOD_TCVTH:
+            value->rValue = model->VDMOStcvth;
+            return(OK);
+        case VDMOS_MOD_VGS_MAX:
+            value->rValue = model->VDMOSvgsMax;
+            return(OK);
+        case VDMOS_MOD_VGD_MAX:
+            value->rValue = model->VDMOSvgdMax;
+            return(OK);
+        case VDMOS_MOD_VDS_MAX:
+            value->rValue = model->VDMOSvdsMax;
+            return(OK);
+        case VDMOS_MOD_VGSR_MAX:
+            value->rValue = model->VDMOSvgsrMax;
+            return(OK);
+        case VDMOS_MOD_VGDR_MAX:
+            value->rValue = model->VDMOSvgdrMax;
             return(OK);
         default:
             return(E_BADPARM);

src/spicelib/devices/vdmos/vdmosmpar.c

@@ -1,6 +1,7 @@
 /**********
 Copyright 1990 Regents of the University of California.  All rights reserved.
 Author: 1985 Thomas L. Quarles
+VDMOS: 2018 Holger Vogt, 2020 Dietmar Warning
 **********/
 
 #include "ngspice/ngspice.h"
@@ -19,9 +20,9 @@ VDMOSmParam(int param, IFvalue *value, GENmodel *inModel)
             model->VDMOStnom = value->rValue + CONSTCtoK;
             model->VDMOStnomGiven = TRUE;
             break;
-        case VDMOS_MOD_VTO:
-            model->VDMOSvt0 = value->rValue;
-            model->VDMOSvt0Given = TRUE;
+        case VDMOS_MOD_VTH:
+            model->VDMOSvth0 = value->rValue;
+            model->VDMOSvth0Given = TRUE;
             break;
         case VDMOS_MOD_KP:
             model->VDMOStransconductance = value->rValue;
@@ -54,20 +55,14 @@ VDMOSmParam(int param, IFvalue *value, GENmodel *inModel)
         case VDMOS_MOD_RQ:
             model->VDMOSqsResistance = value->rValue;
             model->VDMOSqsResistanceGiven = TRUE;
-            if (model->VDMOSqsVoltageGiven)
-                model->VDMOSqsGiven = TRUE;
             break;
         case VDMOS_MOD_VQ:
             model->VDMOSqsVoltage = value->rValue;
             model->VDMOSqsVoltageGiven = TRUE;
-            if (model->VDMOSqsResistanceGiven)
-                model->VDMOSqsGiven = TRUE;
             break;
         case VDMOS_MOD_RB:
             model->VDIOresistance = value->rValue;
             model->VDIOresistanceGiven = TRUE;
-            model->VDIOresistTemp1 = 0;
-            model->VDIOresistTemp2 = 0;
             break;
         case VDMOS_MOD_IS:
             model->VDIOjctSatCur = value->rValue;
@@ -137,10 +132,6 @@ VDMOSmParam(int param, IFvalue *value, GENmodel *inModel)
             model->VDMOSmtr = value->rValue;
             model->VDMOSmtrGiven = TRUE;
             break;
-        case VDMOS_MOD_SUBSLOPE:
-            model->VDMOSsubsl = value->rValue;
-            model->VDMOSsubslGiven = TRUE;
-            break;
         case VDMOS_MOD_SUBSHIFT:
             model->VDMOSsubshift = value->rValue;
             model->VDMOSsubshiftGiven = TRUE;
@@ -150,12 +141,12 @@ VDMOSmParam(int param, IFvalue *value, GENmodel *inModel)
             model->VDMOSksubthresGiven = TRUE;
             break;
         case VDMOS_MOD_BV:
-            model->VDMOSDbv = value->rValue;
-            model->VDMOSDbvGiven = TRUE;
+            model->VDMOSbv = value->rValue;
+            model->VDMOSbvGiven = TRUE;
             break;
         case VDMOS_MOD_IBV:
-            model->VDMOSDibv = value->rValue;
-            model->VDMOSDibvGiven = TRUE;
+            model->VDMOSibv = value->rValue;
+            model->VDMOSibvGiven = TRUE;
             break;
         case VDMOS_MOD_NBV:
             model->VDIObrkdEmissionCoeff = value->rValue;
@@ -166,8 +157,8 @@ VDMOSmParam(int param, IFvalue *value, GENmodel *inModel)
             model->VDMOSrdsGiven = TRUE;
             break;
         case VDMOS_MOD_N:
-            model->VDMOSDn = value->rValue;
-            model->VDMOSDnGiven = TRUE;
+            model->VDMOSn = value->rValue;
+            model->VDMOSnGiven = TRUE;
             break;
         case VDMOS_MOD_TT:
             model->VDIOtransitTime = value->rValue;
@@ -176,12 +167,60 @@ VDMOSmParam(int param, IFvalue *value, GENmodel *inModel)
             model->VDIOtranTimeTemp2 = 0;
             break;
         case VDMOS_MOD_EG:
-            model->VDMOSDeg = value->rValue;
-            model->VDMOSDegGiven = TRUE;
+            model->VDMOSeg = value->rValue;
+            model->VDMOSegGiven = TRUE;
             break;
         case VDMOS_MOD_XTI:
-            model->VDMOSDxti = value->rValue;
-            model->VDMOSDxtiGiven = TRUE;
+            model->VDMOSxti = value->rValue;
+            model->VDMOSxtiGiven = TRUE;
+            break;
+        case  VDMOS_MOD_RTHJC:
+            model->VDMOSrthjc = value->rValue;
+            model->VDMOSrthjcGiven = TRUE;
+            break;
+        case  VDMOS_MOD_RTHCA:
+            model->VDMOSrthca = value->rValue;
+            model->VDMOSrthcaGiven = TRUE;
+            break;
+        case  VDMOS_MOD_CTHJ:
+            model->VDMOScthj = value->rValue;
+            model->VDMOScthjGiven = TRUE;
+            break;
+        case  VDMOS_MOD_MU:
+            model->VDMOSmu = value->rValue;
+            model->VDMOSmuGiven = TRUE;
+            break;
+        case  VDMOS_MOD_TEXP0:
+            model->VDMOStexp0 = value->rValue;
+            model->VDMOStexp0Given = TRUE;
+            break;
+        case  VDMOS_MOD_TEXP1:
+            model->VDMOStexp1 = value->rValue;
+            model->VDMOStexp1Given = TRUE;
+            break;
+        case VDMOS_MOD_TCVTH:
+            model->VDMOStcvth = value->rValue;
+            model->VDMOStcvthGiven = TRUE;
+            break;
+        case VDMOS_MOD_VGS_MAX:
+            model->VDMOSvgsMax = value->rValue;
+            model->VDMOSvgsMaxGiven = TRUE;
+            break;
+        case VDMOS_MOD_VGD_MAX:
+            model->VDMOSvgdMax = value->rValue;
+            model->VDMOSvgdMaxGiven = TRUE;
+            break;
+        case VDMOS_MOD_VDS_MAX:
+            model->VDMOSvdsMax = value->rValue;
+            model->VDMOSvdsMaxGiven = TRUE;
+            break;
+        case VDMOS_MOD_VGSR_MAX:
+            model->VDMOSvgsrMax = value->rValue;
+            model->VDMOSvgsrMaxGiven = TRUE;
+            break;
+        case VDMOS_MOD_VGDR_MAX:
+            model->VDMOSvgdrMax = value->rValue;
+            model->VDMOSvgdrMaxGiven = TRUE;
             break;
         default:
             return(E_BADPARM);

src/spicelib/devices/vdmos/vdmosnoi.c

@@ -36,6 +36,7 @@ VDMOSnoise (int mode, int operation, GENmodel *genmodel, CKTcircuit *ckt,
     double noizDens[VDMOSNSRCS];
     double lnNdens[VDMOSNSRCS];
     int i;
+    double tempRatioSH;
 
     /* define the names of the noise sources */
 
@@ -90,17 +91,21 @@ VDMOSnoise (int mode, int operation, GENmodel *genmodel, CKTcircuit *ckt,
         switch (mode) {
 
         case N_DENS:
+            if ((inst->VDMOStnodeoutGiven) && (model->VDMOSrthjc != 0.0))  
+                tempRatioSH = inst->VDMOSTempSH / ckt->CKTtemp;
+            else
+                tempRatioSH = 1.0;
             NevalSrc(&noizDens[VDMOSRDNOIZ],&lnNdens[VDMOSRDNOIZ],
                  ckt,THERMNOISE,inst->VDMOSdNodePrime,inst->VDMOSdNode,
-                 inst->VDMOSdrainConductance);
+                 inst->VDMOSdrainConductance * tempRatioSH);
 
             NevalSrc(&noizDens[VDMOSRSNOIZ],&lnNdens[VDMOSRSNOIZ],
                  ckt,THERMNOISE,inst->VDMOSsNodePrime,inst->VDMOSsNode,
-                 inst->VDMOSsourceConductance);
+                 inst->VDMOSsourceConductance * tempRatioSH);
 
             NevalSrc(&noizDens[VDMOSIDNOIZ],&lnNdens[VDMOSIDNOIZ],
                  ckt,THERMNOISE,inst->VDMOSdNodePrime,inst->VDMOSsNodePrime,
-                                 (2.0/3.0 * fabs(inst->VDMOSgm)));
+                                 (2.0/3.0 * fabs(inst->VDMOSgm)) * tempRatioSH);
 
             NevalSrc(&noizDens[VDMOSFLNOIZ], NULL, ckt,
                  N_GAIN,inst->VDMOSdNodePrime, inst->VDMOSsNodePrime,

src/spicelib/devices/vdmos/vdmospar.c

@@ -2,6 +2,7 @@
 Copyright 1990 Regents of the University of California.  All rights reserved.
 Author: 1985 Thomas L. Quarles
 Modified: 2000 AlansFixes
+VDMOS: 2018 Holger Vogt, 2020 Dietmar Warning
 **********/
 /*
  */
@@ -60,6 +61,10 @@ VDMOSparam(int param, IFvalue *value, GENinstance *inst, IFvalue *select)
             here->VDMOSicVGS = value->rValue;
             here->VDMOSicVGSGiven = TRUE;
             break;
+        case VDMOS_TNODEOUT:
+            here->VDMOStnodeout = (value->iValue != 0);
+            here->VDMOStnodeoutGiven = TRUE;
+            break;
         case VDMOS_IC:
             switch(value->v.numValue){
                 case 2:

src/spicelib/devices/vdmos/vdmospzld.c

@@ -44,8 +44,13 @@ VDMOSpzLoad(GENmodel *inModel, CKTcircuit *ckt, SPcomplex *s)
             capgd = ( 2* *(ckt->CKTstate0+here->VDMOScapgd));
             xgs = capgs;
             xgd = capgd;
-            /*printf("vdmos: xgs=%g, xgd=%g\n",
-                    xgs,xgd);*/
+
+            /* body diode */
+            double gspr, geq, xceq;
+            gspr = here->VDIOtConductance;
+            geq = *(ckt->CKTstate0 + here->VDIOconduct);
+            xceq = *(ckt->CKTstate0 + here->VDIOcapCurrent);
+
             /*
              *    load matrix
              */
@@ -67,21 +72,34 @@ VDMOSpzLoad(GENmodel *inModel, CKTcircuit *ckt, SPcomplex *s)
             *(here->VDMOSDdPtr) += here->VDMOSdrainConductance;
             *(here->VDMOSSsPtr) += here->VDMOSsourceConductance;
             *(here->VDMOSDPdpPtr) += here->VDMOSdrainConductance+
-                    here->VDMOSgds+
-                    xrev*(here->VDMOSgm);
+                    here->VDMOSgds+xrev*(here->VDMOSgm);
             *(here->VDMOSSPspPtr) += here->VDMOSsourceConductance+
-                    here->VDMOSgds+
-                    xnrm*(here->VDMOSgm);
+                    here->VDMOSgds+xnrm*(here->VDMOSgm);
             *(here->VDMOSDdpPtr) -= here->VDMOSdrainConductance;
             *(here->VDMOSSspPtr) -= here->VDMOSsourceConductance;
             *(here->VDMOSDPdPtr) -= here->VDMOSdrainConductance;
-            *(here->VDMOSDPgPtr) += (xnrm-xrev)*here->VDMOSgm;
-            *(here->VDMOSDPspPtr) -= here->VDMOSgds+
-                    xnrm*(here->VDMOSgm);
-            *(here->VDMOSSPgPtr) -= (xnrm-xrev)*here->VDMOSgm;
+            *(here->VDMOSDPgpPtr) += (xnrm-xrev)*here->VDMOSgm;
+            *(here->VDMOSDPspPtr) -= here->VDMOSgds+xnrm*(here->VDMOSgm);
+            *(here->VDMOSSPgpPtr) -= (xnrm-xrev)*here->VDMOSgm;
             *(here->VDMOSSPsPtr) -= here->VDMOSsourceConductance;
-            *(here->VDMOSSPdpPtr) -= here->VDMOSgds+
-                    xrev*(here->VDMOSgm);
+            *(here->VDMOSSPdpPtr) -= here->VDMOSgds+xrev*(here->VDMOSgm);
+            /* gate resistor */
+            *(here->VDMOSGgPtr) += (here->VDMOSgateConductance);
+            *(here->VDMOSGPgpPtr) += (here->VDMOSgateConductance);
+            *(here->VDMOSGgpPtr) -= here->VDMOSgateConductance;
+            *(here->VDMOSGPgPtr) -= here->VDMOSgateConductance;
+            /* body diode */
+            *(here->VDMOSSsPtr) += gspr;
+            *(here->VDMOSDdPtr) += geq + xceq * s->real;
+            *(here->VDMOSDdPtr +1 ) += xceq * s->imag;
+            *(here->VDIORPrpPtr) += geq + gspr + xceq * s->real;
+            *(here->VDIORPrpPtr +1) += xceq * s->imag;
+            *(here->VDIOSrpPtr) -= gspr;
+            *(here->VDIODrpPtr) -= geq + xceq * s->real;
+            *(here->VDIODrpPtr +1) -= xceq * s->imag;
+            *(here->VDIORPsPtr) -= gspr;
+            *(here->VDIORPdPtr) -= geq + xceq * s->real;
+            *(here->VDIORPdPtr +1 ) -= xceq * s->imag;
 
         }
     }

src/spicelib/devices/vdmos/vdmosset.c

@@ -2,6 +2,7 @@
 Copyright 1990 Regents of the University of California.  All rights reserved.
 Author: 1985 Thomas L. Quarles
 Modified: 2000 AlansFixes
+VDMOS: 2018 Holger Vogt, 2020 Dietmar Warning
 **********/
 
 /* load the VDMOS device structure with those pointers needed later
@@ -27,87 +28,136 @@ VDMOSsetup(SMPmatrix *matrix, GENmodel *inModel, CKTcircuit *ckt,
     /*  loop through all the VDMOS device models */
     for (; model != NULL; model = VDMOSnextModel(model)) {
 
-        if (!model->VDMOStypeGiven) {
+        if (!model->VDMOStypeGiven)
             model->VDMOStype = NMOS;
-        }
-        if (!model->VDIOjctSatCurGiven) {
+
+        if (!model->VDIOjctSatCurGiven)
             model->VDIOjctSatCur = 1e-14;
-        }
-        if (!model->VDMOStransconductanceGiven) {
+
+        if (!model->VDMOStransconductanceGiven)
             model->VDMOStransconductance = 1;
-        }
-        if (!model->VDMOSvt0Given) {
-            model->VDMOSvt0 = 0;
-        }
-        if (!model->VDIOjunctionPotGiven) {
+
+        if (!model->VDMOSvth0Given)
+            model->VDMOSvth0 = 0;
+
+        if (!model->VDIOjunctionPotGiven)
             model->VDIOjunctionPot = .8;
-        }
-        if (!model->VDIOgradCoeffGiven) {
+
+        if (!model->VDIOgradCoeffGiven)
             model->VDIOgradCoeff = .5;
-        }
-        if (!model->VDIOdepletionCapCoeffGiven) {
+
+        if (!model->VDIOdepletionCapCoeffGiven)
             model->VDIOdepletionCapCoeff = .5;
-        }
-        if (!model->VDMOSphiGiven) {
+
+        if (!model->VDMOSphiGiven)
             model->VDMOSphi = .6;
-        }
-        if (!model->VDMOSlambdaGiven) {
+
+        if (!model->VDMOSlambdaGiven)
             model->VDMOSlambda = 0;
-        }
-        if (!model->VDMOSthetaGiven) {
+
+        if (!model->VDMOSthetaGiven)
             model->VDMOStheta = 0;
-        }
-        if (!model->VDMOSfNcoefGiven) {
+
+        if (!model->VDMOSfNcoefGiven)
             model->VDMOSfNcoef = 0;
-        }
-        if (!model->VDMOSfNexpGiven) {
+
+        if (!model->VDMOSfNexpGiven)
             model->VDMOSfNexp = 1;
-        }
-        if (!model->VDMOScgdminGiven) {
+
+        if (!model->VDMOScgdminGiven)
             model->VDMOScgdmin = 0;
-        }
-        if (!model->VDMOScgdmaxGiven) {
+
+        if (!model->VDMOScgdmaxGiven)
             model->VDMOScgdmax = 0;
-        }
-        if (!model->VDMOScgsGiven) {
+
+        if (!model->VDMOScgsGiven)
             model->VDMOScgs = 0;
-        }
-        if (!model->VDMOSaGiven) {
+
+        if (!model->VDMOSaGiven)
             model->VDMOSa = 1.;
-        }
-        if (!model->VDMOSsubslGiven) {
-            model->VDMOSsubsl = 0;
-        }
-        if (!model->VDMOSsubshiftGiven) {
+
+        if (!model->VDMOSsubshiftGiven)
             model->VDMOSsubshift = 0;
-        }
-        if (!model->VDMOSksubthresGiven) {
-            model->VDMOSksubthres = 0;
-        }
-        if (!model->VDMOSmtrGiven) {
+
+        if (!model->VDMOSksubthresGiven)
+            model->VDMOSksubthres = 0.1;
+
+        if (!model->VDMOSmtrGiven)
             model->VDMOSmtr = 1.;
-        }
-        if (!model->VDMOSDbvGiven) {
-            model->VDMOSDbv = 1.0e30;
-        }
-        if (!model->VDMOSDibvGiven) {
-            model->VDMOSDibv = 1.0e-10;
-        }
-        if (!model->VDIObrkdEmissionCoeffGiven) {
+
+        if (!model->VDMOSbvGiven)
+            model->VDMOSbv = 1.0e30;
+
+        if (!model->VDMOSibvGiven)
+            model->VDMOSibv = 1.0e-10;
+
+        if (!model->VDIObrkdEmissionCoeffGiven)
             model->VDIObrkdEmissionCoeff = 1.;
-        }
-        if (!model->VDMOSrdsGiven) {
-            model->VDMOSrds = 1.0e30;
-        }
-        if (!model->VDMOSDnGiven) {
-            model->VDMOSDn = 1.;
-        }
-        if (!model->VDIOtransitTimeGiven) {
+
+        if (!model->VDMOSdrainResistanceGiven)
+            model->VDMOSdrainResistance = 1.0e-03;
+
+        if (!model->VDMOSsourceResistanceGiven)
+            model->VDMOSsourceResistance = 1.0e-03;
+
+        if (!model->VDMOSgateResistanceGiven)
+            model->VDMOSgateResistance = 1.0e-03;
+
+        if (!model->VDMOSrdsGiven)
+            model->VDMOSrds = 1.0e+15;
+
+        if (!model->VDIOresistanceGiven)
+            model->VDIOresistance = 1.0e-03;
+
+        if (!model->VDMOSnGiven)
+            model->VDMOSn = 1.;
+
+        if (!model->VDIOtransitTimeGiven)
             model->VDIOtransitTime = 0.;
-        }
-        if (!model->VDMOSDegGiven) {
-            model->VDMOSDeg = 1.11;
-        }
+
+        if (!model->VDMOSegGiven)
+            model->VDMOSeg = 1.11;
+
+        if (!model->VDMOSrthjcGiven)
+            model->VDMOSrthjc = 0;
+
+        if (!model->VDMOSrthcaGiven)
+            model->VDMOSrthca = 100;
+
+        if (!model->VDMOScthjGiven)
+            model->VDMOScthj = 10e-06;
+
+        if (!model->VDMOSmuGiven)
+            model->VDMOSmu = -1.5;
+
+        if (!model->VDMOStcvthGiven)
+            model->VDMOStcvth = 0.0;
+
+        if (!model->VDMOStexp0Given)
+            model->VDMOStexp0 = 1.5;
+
+        if (!model->VDMOStexp1Given)
+            model->VDMOStexp1 = 0.3;
+
+        if (!model->VDMOSvgsMaxGiven)
+            model->VDMOSvgsMax = 1e99;
+
+        if (!model->VDMOSvgdMaxGiven)
+            model->VDMOSvgdMax = 1e99;
+
+        if (!model->VDMOSvdsMaxGiven)
+            model->VDMOSvdsMax = 1e99;
+
+        if (!model->VDMOSvgsrMaxGiven)
+            model->VDMOSvgsrMax = 1e99;
+
+        if (!model->VDMOSvgdrMaxGiven)
+            model->VDMOSvgdrMax = 1e99;
+
+        if ((model->VDMOSqsResistanceGiven) && (model->VDMOSqsVoltageGiven))
+            model->VDMOSqsGiven = 1;
+        else
+            model->VDMOSqsGiven = 0;
 
         /* loop through all the instances of the model */
         for (here = VDMOSinstances(model); here != NULL;
@@ -129,6 +179,45 @@ VDMOSsetup(SMPmatrix *matrix, GENmodel *inModel, CKTcircuit *ckt,
             if (!here->VDMOSvonGiven) {
                 here->VDMOSvon = 0;
             }
+            if(!here->VDMOSmGiven) {
+                here->VDMOSm = 1;
+            }
+            if(!here->VDMOSlGiven) {
+                here->VDMOSl = 1;
+            }
+            if(!here->VDMOSwGiven) {
+                here->VDMOSw = 1;
+            }
+            if (model->VDMOSdrainResistance != 0) {
+                here->VDMOSdrainConductance = here->VDMOSm / model->VDMOSdrainResistance;
+            } else {
+                here->VDMOSdrainConductance = here->VDMOSm / 1.0e-03;
+            }
+            if (model->VDMOSsourceResistance != 0) {
+                here->VDMOSsourceConductance = here->VDMOSm / model->VDMOSsourceResistance;
+            } else {
+                here->VDMOSsourceConductance = here->VDMOSm / 1.0e-03;
+            }
+            if (model->VDMOSgateResistance != 0) {
+                here->VDMOSgateConductance = here->VDMOSm / model->VDMOSgateResistance;
+            } else {
+                here->VDMOSgateConductance = here->VDMOSm / 1.0e-03;
+            }
+            if (model->VDMOSrdsGiven) {
+                if (model->VDMOSrds != 0) {
+                    here->VDMOSdsConductance = here->VDMOSm / model->VDMOSrds;
+                } else {
+                    here->VDMOSdsConductance = 1e-15;
+                }
+            } else {
+                here->VDMOSdsConductance = 1e-15;
+            }
+            if (model->VDIOresistance != 0) {
+                here->VDIOconductance = here->VDMOSm / model->VDIOresistance;
+            } else {
+                here->VDIOconductance = here->VDMOSm / 1.0e-03;
+            }
+
             if (model->VDMOSdrainResistance != 0) {
                 if (here->VDMOSdNodePrime == 0) {
                     error = CKTmkVolt(ckt, &tmp, here->VDMOSname, "drain");
@@ -153,35 +242,34 @@ VDMOSsetup(SMPmatrix *matrix, GENmodel *inModel, CKTcircuit *ckt,
                 here->VDMOSdNodePrime = here->VDMOSdNode;
             }
 
-            if (model->VDMOSsourceResistance != 0) {
-                if (here->VDMOSsNodePrime == 0) {
-                    error = CKTmkVolt(ckt, &tmp, here->VDMOSname, "source");
+            if (model->VDMOSgateResistance != 0 ) {
+                if (here->VDMOSgNodePrime == 0) {
+                    error = CKTmkVolt(ckt, &tmp, here->VDMOSname, "gate");
                     if (error) return(error);
-                    here->VDMOSsNodePrime = tmp->number;
+                    here->VDMOSgNodePrime = tmp->number;
 
                     if (ckt->CKTcopyNodesets) {
                         CKTnode *tmpNode;
                         IFuid tmpName;
 
-                        if (CKTinst2Node(ckt, here, 3, &tmpNode, &tmpName) == OK) {
+                        if (CKTinst2Node(ckt, here, 2, &tmpNode, &tmpName) == OK) {
                             if (tmpNode->nsGiven) {
                                 tmp->nodeset = tmpNode->nodeset;
                                 tmp->nsGiven = tmpNode->nsGiven;
                             }
                         }
                     }
-
                 }
             }
             else {
-                here->VDMOSsNodePrime = here->VDMOSsNode;
+                here->VDMOSgNodePrime = here->VDMOSgNode;
             }
 
-            if (model->VDMOSgateResistance != 0 ) {
-                if (here->VDMOSgNodePrime == 0) {
-                    error = CKTmkVolt(ckt, &tmp, here->VDMOSname, "gate");
+            if (model->VDMOSsourceResistance != 0) {
+                if (here->VDMOSsNodePrime == 0) {
+                    error = CKTmkVolt(ckt, &tmp, here->VDMOSname, "source");
                     if (error) return(error);
-                    here->VDMOSgNodePrime = tmp->number;
+                    here->VDMOSsNodePrime = tmp->number;
 
                     if (ckt->CKTcopyNodesets) {
                         CKTnode *tmpNode;
@@ -194,15 +282,16 @@ VDMOSsetup(SMPmatrix *matrix, GENmodel *inModel, CKTcircuit *ckt,
                             }
                         }
                     }
+
                 }
             }
             else {
-                here->VDMOSgNodePrime = here->VDMOSgNode;
+                here->VDMOSsNodePrime = here->VDMOSsNode;
             }
 
             if (model->VDIOresistance != 0 ) {
                 if (here->VDIOposPrimeNode == 0) {
-                    error = CKTmkVolt(ckt, &tmp, here->VDMOSname, "bulk diode");
+                    error = CKTmkVolt(ckt, &tmp, here->VDMOSname, "body diode");
                     if (error) return(error);
                     here->VDIOposPrimeNode = tmp->number;
 
@@ -223,6 +312,18 @@ VDMOSsetup(SMPmatrix *matrix, GENmodel *inModel, CKTcircuit *ckt,
                 here->VDIOposPrimeNode = here->VDMOSsNode;
             }
 
+            if ((here->VDMOStnodeoutGiven) && (model->VDMOSrthjc!=0.0)) {
+                if(here->VDMOSvcktTbranch == 0) {
+                    error = CKTmkCur(ckt,&tmp,here->VDMOSname,"VcktTemp");
+                    if(error) return(error);
+                    here->VDMOSvcktTbranch = tmp->number;
+                } 
+                if (here->VDMOStNodePrime == 0) {
+                    error = CKTmkVolt(ckt, &tmp, here->VDMOSname, "thermal node");
+                    if (error) return(error);
+                    here->VDMOStNodePrime = tmp->number;
+                }
+            }
 
             /* macro to make elements with built in test for out of memory */
 #define TSTALLOC(ptr,first,second) \
@@ -230,6 +331,25 @@ do { if((here->ptr = SMPmakeElt(matrix, here->first, here->second)) == NULL){\
     return(E_NOMEM);\
 } } while(0)
 
+            if ((here->VDMOStnodeoutGiven) && (model->VDMOSrthjc!=0.0)) {
+                TSTALLOC(VDMOSTemptempPtr, VDMOStempNode, VDMOStempNode);
+                TSTALLOC(VDMOSTempdpPtr, VDMOStempNode, VDMOSdNodePrime);
+                TSTALLOC(VDMOSTempspPtr, VDMOStempNode, VDMOSsNodePrime);
+                TSTALLOC(VDMOSTempgpPtr, VDMOStempNode, VDMOSgNodePrime);
+                TSTALLOC(VDMOSGPtempPtr, VDMOSgNodePrime, VDMOStempNode);
+                TSTALLOC(VDMOSDPtempPtr, VDMOSdNodePrime, VDMOStempNode);
+                TSTALLOC(VDMOSSPtempPtr, VDMOSsNodePrime, VDMOStempNode);
+
+                TSTALLOC(VDMOSTcasetcasePtr, VDMOStcaseNode, VDMOStcaseNode); /* Rthjc between tj and tcase*/
+                TSTALLOC(VDMOSTcasetempPtr, VDMOStcaseNode, VDMOStempNode);
+                TSTALLOC(VDMOSTemptcasePtr, VDMOStempNode, VDMOStcaseNode);
+                TSTALLOC(VDMOSTptpPtr, VDMOStNodePrime, VDMOStNodePrime); /* Rthca between tcase and Vsrc */
+                TSTALLOC(VDMOSTptcasePtr, VDMOStNodePrime, VDMOStempNode);
+                TSTALLOC(VDMOSTcasetpPtr, VDMOStempNode, VDMOStNodePrime);
+                TSTALLOC(VDMOSCktTcktTPtr, VDMOSvcktTbranch, VDMOSvcktTbranch); /* Vsrc=cktTemp to gnd */
+                TSTALLOC(VDMOSCktTtpPtr, VDMOSvcktTbranch, VDMOStNodePrime);
+                TSTALLOC(VDMOSTpcktTPtr, VDMOStNodePrime, VDMOSvcktTbranch);
+            }
             TSTALLOC(VDMOSDdPtr, VDMOSdNode, VDMOSdNode);
             TSTALLOC(VDMOSGgPtr, VDMOSgNode, VDMOSgNode);
             TSTALLOC(VDMOSSsPtr, VDMOSsNode, VDMOSsNode);
@@ -294,6 +414,16 @@ VDMOSunsetup(GENmodel *inModel, CKTcircuit *ckt)
                 && here->VDIOposPrimeNode != here->VDMOSsNode)
                 CKTdltNNum(ckt, here->VDIOposPrimeNode);
             here->VDIOposPrimeNode = 0;
+
+            if ((here->VDMOStnodeoutGiven) && (model->VDMOSrthjc!=0.0)) {
+                if (here->VDMOStNodePrime > 0)
+                    CKTdltNNum(ckt, here->VDMOStNodePrime);
+                here->VDMOStNodePrime = 0;
+                if (here->VDMOSvcktTbranch > 0)
+                    CKTdltNNum(ckt, here->VDMOSvcktTbranch);
+                here->VDMOSvcktTbranch = 0;
+            }
+
         }
     }
     return OK;

src/spicelib/devices/vdmos/vdmossoachk.c

@@ -0,0 +1,142 @@
+/**********
+Copyright 2013 Dietmar Warning. All rights reserved.
+Author:   2013 Dietmar Warning
+**********/
+
+#include "ngspice/ngspice.h"
+#include "ngspice/cktdefs.h"
+#include "vdmosdefs.h"
+#include "ngspice/trandefs.h"
+#include "ngspice/sperror.h"
+#include "ngspice/suffix.h"
+#include "ngspice/cpdefs.h"
+
+
+int
+VDMOSsoaCheck(CKTcircuit *ckt, GENmodel *inModel)
+{
+    VDMOSmodel *model = (VDMOSmodel *) inModel;
+    VDMOSinstance *here;
+    double vgs, vgd, vds;    /* actual mos voltages */
+    int maxwarns;
+    static int warns_vgs = 0, warns_vgd = 0, warns_vds = 0;
+
+    if (!ckt) {
+        warns_vgs = 0;
+        warns_vgd = 0;
+        warns_vds = 0;
+        return OK;
+    }
+
+    maxwarns = ckt->CKTsoaMaxWarns;
+
+    for (; model; model = VDMOSnextModel(model)) {
+
+        for (here = VDMOSinstances(model); here; here = VDMOSnextInstance(here)) {
+
+            vgs = ckt->CKTrhsOld [here->VDMOSgNode] -
+                  ckt->CKTrhsOld [here->VDMOSsNodePrime];
+
+            vgd = ckt->CKTrhsOld [here->VDMOSgNode] -
+                  ckt->CKTrhsOld [here->VDMOSdNodePrime];
+
+            vds = ckt->CKTrhsOld [here->VDMOSdNodePrime] -
+                  ckt->CKTrhsOld [here->VDMOSsNodePrime];
+
+            if (!model->VDMOSvgsrMaxGiven) {
+                if (fabs(vgs) > model->VDMOSvgsMax)
+                    if (warns_vgs < maxwarns) {
+                        soa_printf(ckt, (GENinstance*) here,
+                                   "Vgs=%g has exceeded Vgs_max=%g\n",
+                                   vgs, model->VDMOSvgsMax);
+                        warns_vgs++;
+                    }
+            } else {
+                if (model->VDMOStype > 0) {
+                    if (vgs > model->VDMOSvgsMax)
+                        if (warns_vgs < maxwarns) {
+                            soa_printf(ckt, (GENinstance*) here,
+                                       "Vgs=%g has exceeded Vgs_max=%g\n",
+                                       vgs, model->VDMOSvgsMax);
+                            warns_vgs++;
+                        }
+                    if (-1*vgs > model->VDMOSvgsrMax)
+                        if (warns_vgs < maxwarns) {
+                            soa_printf(ckt, (GENinstance*) here,
+                                       "Vgs=%g has exceeded Vgsr_max=%g\n",
+                                       vgs, model->VDMOSvgsrMax);
+                            warns_vgs++;
+                        }
+                } else {
+                    if (vgs > model->VDMOSvgsrMax)
+                        if (warns_vgs < maxwarns) {
+                            soa_printf(ckt, (GENinstance*) here,
+                                       "Vgs=%g has exceeded Vgsr_max=%g\n",
+                                       vgs, model->VDMOSvgsrMax);
+                            warns_vgs++;
+                        }
+                    if (-1*vgs > model->VDMOSvgsMax)
+                        if (warns_vgs < maxwarns) {
+                            soa_printf(ckt, (GENinstance*) here,
+                                       "Vgs=%g has exceeded Vgs_max=%g\n",
+                                       vgs, model->VDMOSvgsMax);
+                            warns_vgs++;
+                        }
+                }
+            }
+
+            if (!model->VDMOSvgdrMaxGiven) {
+                if (fabs(vgd) > model->VDMOSvgdMax)
+                    if (warns_vgd < maxwarns) {
+                        soa_printf(ckt, (GENinstance*) here,
+                                   "Vgd=%g has exceeded Vgd_max=%g\n",
+                                   vgd, model->VDMOSvgdMax);
+                        warns_vgd++;
+                    }
+            } else {
+                if (model->VDMOStype > 0) {
+                    if (vgd > model->VDMOSvgdMax)
+                        if (warns_vgd < maxwarns) {
+                            soa_printf(ckt, (GENinstance*) here,
+                                       "Vgd=%g has exceeded Vgd_max=%g\n",
+                                       vgd, model->VDMOSvgdMax);
+                            warns_vgd++;
+                        }
+                    if (-1*vgd > model->VDMOSvgdrMax)
+                        if (warns_vgd < maxwarns) {
+                            soa_printf(ckt, (GENinstance*) here,
+                                       "Vgd=%g has exceeded Vgdr_max=%g\n",
+                                       vgd, model->VDMOSvgdrMax);
+                            warns_vgd++;
+                        }
+                } else {
+                    if (vgd > model->VDMOSvgdrMax)
+                        if (warns_vgd < maxwarns) {
+                            soa_printf(ckt, (GENinstance*) here,
+                                       "Vgd=%g has exceeded Vgdr_max=%g\n",
+                                       vgd, model->VDMOSvgdrMax);
+                            warns_vgd++;
+                        }
+                    if (-1*vgd > model->VDMOSvgdMax)
+                        if (warns_vgd < maxwarns) {
+                            soa_printf(ckt, (GENinstance*) here,
+                                       "Vgd=%g has exceeded Vgd_max=%g\n",
+                                       vgd, model->VDMOSvgdMax);
+                            warns_vgd++;
+                        }
+                }
+            }
+
+            if (fabs(vds) > model->VDMOSvdsMax)
+                if (warns_vds < maxwarns) {
+                    soa_printf(ckt, (GENinstance*) here,
+                               "Vds=%g has exceeded Vds_max=%g\n",
+                               vds, model->VDMOSvdsMax);
+                    warns_vds++;
+                }
+
+        }
+    }
+
+    return OK;
+}

src/spicelib/devices/vdmos/vdmostemp.c

@@ -2,6 +2,7 @@
 Copyright 1990 Regents of the University of California.  All rights reserved.
 Author: 1985 Thomas L. Quarles
 Modified: 2000 AlansFixes
+VDMOS: 2018 Holger Vogt, 2020 Dietmar Warning
 **********/
 
 #include "ngspice/ngspice.h"
@@ -21,7 +22,7 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
     double fact1,fact2;
     double kt,kt1;
     double arg1;
-    double ratio,ratio4;
+    double ratio;
     double phio;
     double pbfact1,pbfact;
     double vt,vtnom;
@@ -45,22 +46,21 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
         pbfact1 = -2*vtnom *(1.5*log(fact1)+CHARGE*arg1);
 
         /* now model parameter preprocessing */
-
         if (model->VDMOSphi <= 0.0) {
             SPfrontEnd->IFerrorf(ERR_FATAL,
                 "%s: Phi is not positive.", model->VDMOSmodName);
             return(E_BADPARM);
         }
 
-        model->VDMOSoxideCapFactor = 0;
+        model->VDMOSoxideCapFactor = 3.9 * 8.854214871e-12 / 1e-07; /* use default Tox of 100nm */
 
-        /* bulk diode model */
+        /* body diode model */
         /* limit activation energy to min of .1 */
-        if (model->VDMOSDeg<.1) {
+        if (model->VDMOSeg<.1) {
             SPfrontEnd->IFerrorf(ERR_WARNING,
-                "%s: bulk diode activation energy too small, limited to 0.1",
+                "%s: body diode activation energy too small, limited to 0.1",
                 model->VDMOSmodName);
-            model->VDMOSDeg = .1;
+            model->VDMOSeg = .1;
         }
         /* limit depletion cap coeff to max of .95 */
         if (model->VDIOdepletionCapCoeff>.95) {
@@ -72,12 +72,7 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
         /* set lower limit of saturation current */
         if (model->VDIOjctSatCur < ckt->CKTepsmin)
             model->VDIOjctSatCur = ckt->CKTepsmin;
-        if ((!model->VDIOresistanceGiven) || (model->VDIOresistance == 0)) {
-            model->VDIOconductance = 0.0;
-        }
-        else {
-            model->VDIOconductance = 1 / model->VDIOresistance;
-        }
+
         xfc = log(1 - model->VDIOdepletionCapCoeff);
 
         /* loop through all instances of the model */
@@ -86,15 +81,29 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
             double arg;     /* 1 - fc */
 
             /* perform the parameter defaulting */
-
             if(!here->VDMOSdtempGiven) {
                 here->VDMOSdtemp = 0.0;
             }
             if(!here->VDMOStempGiven) {
                 here->VDMOStemp = ckt->CKTtemp + here->VDMOSdtemp;
             }
-            vt = here->VDMOStemp * CONSTKoverQ;
+
+            double dt = here->VDMOStemp - model->VDMOStnom;
+
+            /* vdmos temperature model */
             ratio = here->VDMOStemp/model->VDMOStnom;
+            here->VDMOStTransconductance = model->VDMOStransconductance 
+                                           * here->VDMOSm * here->VDMOSw / here->VDMOSl 
+                                           * pow(ratio, model->VDMOSmu);
+
+            here->VDMOStVth = model->VDMOSvth0 - model->VDMOStype * model->VDMOStcvth * dt;
+
+            here->VDMOSdrainResistance =  model->VDMOSdrainResistance / here->VDMOSm * pow(ratio, model->VDMOStexp0);
+
+            if (model->VDMOSqsGiven)
+                here->VDMOSqsResistance = model->VDMOSqsResistance / here->VDMOSm * pow(ratio, model->VDMOStexp1);
+
+            vt = here->VDMOStemp * CONSTKoverQ;
             fact2 = here->VDMOStemp/REFTEMP;
             kt = here->VDMOStemp * CONSTboltz;
             egfet = 1.16-(7.02e-4*here->VDMOStemp*here->VDMOStemp)/
@@ -102,82 +111,15 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
             arg = -egfet/(kt+kt)+1.1150877/(CONSTboltz*(REFTEMP+REFTEMP));
             pbfact = -2*vt *(1.5*log(fact2)+CHARGE*arg);
 
-            if(!here->VDMOSmGiven) {
-                here->VDMOSm = 1;
-            }
-            if(!here->VDMOSlGiven) {
-                here->VDMOSl = 1;
-            }
-            if(!here->VDMOSwGiven) {
-                here->VDMOSw = 1;
-            }
-
-            ratio4 = ratio * sqrt(ratio);
-            here->VDMOStTransconductance = model->VDMOStransconductance / ratio4;
             phio = (model->VDMOSphi - pbfact1) / fact1;
-            here->VDMOStPhi = fact2 * phio + pbfact;
-            here->VDMOStVto = model->VDMOSvt0;
-
-            here->VDMOSf2d = 0;
-            here->VDMOSf3d = 0;
-            here->VDMOSf4d = 0;
+            here->VDMOStPhi = fact2 * phio + pbfact; /* needed for distortion analysis */
 
-            here->VDMOSf2s = 0;
-            here->VDMOSf3s = 0;
-            here->VDMOSf4s = 0;
-
-
-            if (model->VDMOSdrainResistanceGiven) {
-                if (model->VDMOSdrainResistance != 0) {
-                    here->VDMOSdrainConductance = here->VDMOSm /
-                        model->VDMOSdrainResistance;
-                }
-                else {
-                    here->VDMOSdrainConductance = 0;
-                }
-            } else {
-                here->VDMOSdrainConductance = 0;
-            }
-            if(model->VDMOSsourceResistanceGiven) {
-                if(model->VDMOSsourceResistance != 0) {
-                   here->VDMOSsourceConductance = here->VDMOSm /
-                                         model->VDMOSsourceResistance;
-                } else {
-                    here->VDMOSsourceConductance = 0;
-                }
-            } else {
-                here->VDMOSsourceConductance = 0;
-            }
-            if (model->VDMOSgateResistanceGiven) {
-                if (model->VDMOSgateResistance != 0) {
-                    here->VDMOSgateConductance = here->VDMOSm /
-                        model->VDMOSgateResistance;
-                } else {
-                    here->VDMOSgateConductance = 0;
-                }
-            } else {
-                here->VDMOSgateConductance = 0;
-            }
-            if (model->VDMOSrdsGiven) {
-                if (model->VDMOSrds != 0) {
-                    here->VDMOSdsConductance = here->VDMOSm /
-                        model->VDMOSrds;
-                }
-                else {
-                    here->VDMOSdsConductance = 0;
-                }
-            }
-            else {
-                here->VDMOSdsConductance = 0;
-            }
-
-            /* bulk diode model */
+            /* body diode temperature model */
             double pbo, gmaold;
             double gmanew, factor;
             double tBreakdownVoltage, vte, cbv;
-            double xbv, xcbv, tol, iter, dt;
+            double xbv, xcbv, tol, iter;
 
-            dt = here->VDMOStemp - model->VDMOStnom;
             /* Junction grading temperature adjust */
             factor = 1.0 + (model->VDIOgradCoeffTemp1 * dt)
                 + (model->VDIOgradCoeffTemp2 * dt * dt);
@@ -185,7 +127,7 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
 
             pbo = (model->VDIOjunctionPot - pbfact1) / fact1;
             gmaold = (model->VDIOjunctionPot - pbo) / pbo;
-            here->VDIOtJctCap = model->VDIOjunctionCap /
+            here->VDIOtJctCap = here->VDMOSm * model->VDIOjunctionCap /
                 (1 + here->VDIOtGradingCoeff*
                 (400e-6*(model->VDMOStnom - REFTEMP) - gmaold));
             here->VDIOtJctPot = pbfact + fact2*pbo;
@@ -193,10 +135,10 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
             here->VDIOtJctCap *= 1 + here->VDIOtGradingCoeff*
                 (400e-6*(here->VDMOStemp - REFTEMP) - gmanew);
 
-            here->VDIOtSatCur = model->VDIOjctSatCur * exp(
+            here->VDIOtSatCur = here->VDMOSm * model->VDIOjctSatCur * exp(
                 ((here->VDMOStemp / model->VDMOStnom) - 1) *
-                model->VDMOSDeg / (model->VDMOSDn*vt) +
-                model->VDMOSDxti / model->VDMOSDn *
+                model->VDMOSeg / (model->VDMOSn*vt) +
+                model->VDMOSxti / model->VDMOSn *
                 log(here->VDMOStemp / model->VDMOStnom));
 
             /* the defintion of f1, just recompute after temperature adjusting
@@ -209,14 +151,14 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
                 here->VDIOtJctPot;
 
             /* and Vcrit */
-            vte = model->VDMOSDn*vt;
+            vte = model->VDMOSn*vt;
 
             here->VDIOtVcrit = vte * log(vte / (CONSTroot2*here->VDIOtSatCur));
 
             /* limit junction potential to max of 1/FC */
-            if (here->VDIOtDepCap > 1.0) {
-                here->VDIOtJctPot = 1.0 / model->VDMOSDn;
-                here->VDIOtDepCap = model->VDMOSDn*here->VDIOtJctPot;
+            if (here->VDIOtDepCap > 2.5) {
+                here->VDIOtJctPot = 2.5 / model->VDMOSn;
+                here->VDIOtDepCap = model->VDMOSn*here->VDIOtJctPot;
                 SPfrontEnd->IFerrorf(ERR_WARNING,
                     "%s: junction potential VJ too large, limited to %f",
                     model->VDMOSmodName, here->VDIOtJctPot);
@@ -224,16 +166,16 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
 
             /* and now to compute the breakdown voltage, again, using
             * temperature adjusted basic parameters */
-            if (model->VDMOSDbvGiven) {
+            if (model->VDMOSbvGiven) {
                 /* tlev == 0 */
-                tBreakdownVoltage = fabs(model->VDMOSDbv);
+                tBreakdownVoltage = fabs(model->VDMOSbv);
 
-                cbv = model->VDMOSDibv;
+                cbv = model->VDMOSibv;
 
                 if (cbv < here->VDIOtSatCur * tBreakdownVoltage / vt) {
                     cbv = here->VDIOtSatCur * tBreakdownVoltage / vt;
 #ifdef TRACE
-                    SPfrontEnd->IFerrorf(ERR_WARNING, "%s: breakdown current increased to %g to resolve", here->DIOname, cbv);
+                    SPfrontEnd->IFerrorf(ERR_WARNING, "%s: breakdown current increased to %g to resolve", here->VDMOSname, cbv);
                     SPfrontEnd->IFerrorf(ERR_WARNING,
                         "incompatibility with specified saturation current");
 #endif
@@ -252,7 +194,7 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
                         if (fabs(xcbv - cbv) <= tol) goto matched;
                     }
 #ifdef TRACE
-                    SPfrontEnd->IFerrorf(ERR_WARNING, "%s: unable to match forward and reverse diode regions: bv = %g, ibv = %g", here->DIOname, xbv, xcbv);
+                    SPfrontEnd->IFerrorf(ERR_WARNING, "%s: unable to match forward and reverse diode regions: bv = %g, ibv = %g", here->VDMOSname, xbv, xcbv);
 #endif
                 }
             matched:
@@ -264,13 +206,8 @@ VDMOStemp(GENmodel *inModel, CKTcircuit *ckt)
                 + (model->VDIOtranTimeTemp2 * dt * dt);
             here->VDIOtTransitTime = model->VDIOtransitTime * factor;
 
-            /* Series resistance temperature adjust */
-            here->VDIOtConductance = model->VDIOconductance;
-            if (model->VDIOresistanceGiven && model->VDIOresistance != 0.0) {
-                factor = 1.0 + (model->VDIOresistTemp1) * dt
-                    + (model->VDIOresistTemp2 * dt * dt);
-                here->VDIOtConductance = model->VDIOconductance / factor;
-            }
+            /* Series resistance temperature adjust (not implemented yet) */
+            here->VDIOtConductance = here->VDIOconductance;
 
             here->VDIOtF2 = exp((1 + here->VDIOtGradingCoeff)*xfc);
             here->VDIOtF3 = 1 - model->VDIOdepletionCapCoeff*

src/spicelib/parser/inp2m.c

@@ -33,7 +33,7 @@ model_numnodes(int type)
 
     if (type == INPtypelook("VDMOS"))       /* 3 ; VDMOSnames */
     {
-        return 3;
+        return 5;
     }
 
     return 4;

visualc/sharedspice.vcxproj

@@ -2106,6 +2106,7 @@
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmospar.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmospzld.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmosset.c" />
+    <ClCompile Include="..\src\spicelib\devices\vdmos\vdmossoachk.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmostemp.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmostrun.c" />
     <ClCompile Include="..\src\spicelib\devices\vsrc\vsrc.c" />

visualc/vngspice-fftw.vcxproj

@@ -2547,6 +2547,7 @@ lib /machine:x64 /def:..\..\fftw-3.3-dll64\libfftw3-3.def /out:$(IntDir)libfftw3
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmospar.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmospzld.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmosset.c" />
+    <ClCompile Include="..\src\spicelib\devices\vdmos\vdmossoachk.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmostemp.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmostrun.c" />
     <ClCompile Include="..\src\spicelib\devices\vsrc\vsrc.c" />

visualc/vngspice.vcxproj

@@ -2555,6 +2555,7 @@
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmospar.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmospzld.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmosset.c" />
+    <ClCompile Include="..\src\spicelib\devices\vdmos\vdmossoachk.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmostemp.c" />
     <ClCompile Include="..\src\spicelib\devices\vdmos\vdmostrun.c" />
     <ClCompile Include="..\src\spicelib\devices\vsrc\vsrc.c" />