1 /**********
2 Copyright 1990 Regents of the University of California.  All rights reserved.
3 Author: 1985 Thomas L. Quarles
4 Modified: 2000 AlansFixes
5 Modified by Paolo Nenzi 2003 and Dietmar Warning 2012
6 **********/
7 
8 #include "ngspice/ngspice.h"
9 #include "ngspice/devdefs.h"
10 #include "ngspice/cktdefs.h"
11 #include "diodefs.h"
12 #include "ngspice/const.h"
13 #include "ngspice/trandefs.h"
14 #include "ngspice/sperror.h"
15 #include "ngspice/suffix.h"
16 
17 int
DIOload(GENmodel * inModel,CKTcircuit * ckt)18 DIOload(GENmodel *inModel, CKTcircuit *ckt)
19         /* actually load the current resistance value into the
20          * sparse matrix previously provided
21          */
22 {
23     DIOmodel *model = (DIOmodel*)inModel;
24     DIOinstance *here;
25     double arg;
26     double argsw;
27     double capd;
28     double cd, cdb, cdsw, cdb_dT, cdsw_dT;
29     double cdeq;
30     double cdhat;
31     double ceq;
32     double csat;    /* area-scaled saturation current */
33     double csatsw;  /* perimeter-scaled saturation current */
34     double czero;
35     double czof2;
36     double argSW;
37     double czeroSW;
38     double czof2SW;
39     double sargSW;
40     double sqrt_ikr;
41     double sqrt_ikf;
42     double ikf_area_m;
43     double ikr_area_m;
44 
45     double delvd;   /* change in diode voltage temporary */
46     double evd;
47     double evrev;
48     double gd, gdb, gdsw, gen_fac, gen_fac_vd;
49     double t1, evd_rec, cdb_rec, gdb_rec, cdb_rec_dT;
50     double geq;
51     double gspr;    /* area-scaled conductance */
52     double sarg;
53 #ifndef NOBYPASS
54     double tol;     /* temporary for tolerence calculations */
55 #endif
56     double vd;      /* current diode voltage */
57     double vdtemp;
58     double vt;      /* K t / Q */
59     double vte, vtesw, vtetun, vtebrk;
60     int Check_dio=0, Check_th;
61     int error;
62     int SenCond=0;    /* sensitivity condition */
63     double diffcharge, diffchargeSW, deplcharge, deplchargeSW, diffcap, diffcapSW, deplcap, deplcapSW;
64 
65     double deldelTemp, delTemp, Temp;
66     double ceqqth=0.0, Ith=0.0, gcTt=0.0, vrs=0.0;
67     double dIdio_dT, dIth_dVdio=0.0, dIrs_dT=0.0, dIth_dVrs=0.0, dIth_dT=0.0;
68     double argsw_dT, csat_dT, csatsw_dT;
69 
70     /*  loop through all the diode models */
71     for( ; model != NULL; model = DIOnextModel(model)) {
72 
73         /* loop through all the instances of the model */
74         for (here = DIOinstances(model); here != NULL ;
75                 here=DIOnextInstance(here)) {
76 
77             int selfheat = ((here->DIOtempNode > 0) && (here->DIOthermal) && (model->DIOrth0Given));
78 
79             /*
80              *     this routine loads diodes for dc and transient analyses.
81              */
82 
83             if (selfheat)
84                 Check_th = 1;
85             else
86                 Check_th = 0;
87 
88             if(ckt->CKTsenInfo){
89                 if((ckt->CKTsenInfo->SENstatus == PERTURBATION)
90                         && (here->DIOsenPertFlag == OFF))continue;
91                 SenCond = here->DIOsenPertFlag;
92 
93 #ifdef SENSDEBUG
94                 printf("DIOload \n");
95 #endif /* SENSDEBUG */
96 
97             }
98             cdsw = 0.0;
99             cdsw_dT = 0.0;
100             gdsw = 0.0;
101             delTemp = 0.0;
102             vt = CONSTKoverQ * here->DIOtemp;
103             vte = model->DIOemissionCoeff * vt;
104             vtebrk = model->DIObrkdEmissionCoeff * vt;
105             gspr = here->DIOtConductance;
106             /*
107              *   initialization
108              */
109 
110             if(SenCond){
111 
112 #ifdef SENSDEBUG
113                 printf("DIOsenPertFlag = ON \n");
114 #endif /* SENSDEBUG */
115 
116                 if((ckt->CKTsenInfo->SENmode == TRANSEN)&&
117                         (ckt->CKTmode & MODEINITTRAN)) {
118                     vd = *(ckt->CKTstate1 + here->DIOvoltage);
119                     delTemp = *(ckt->CKTstate1 + here->DIOdeltemp);
120                 } else{
121                     vd = *(ckt->CKTstate0 + here->DIOvoltage);
122                     delTemp = *(ckt->CKTstate0 + here->DIOdeltemp);
123                 }
124 
125 #ifdef SENSDEBUG
126                 printf("vd = %.7e \n",vd);
127 #endif /* SENSDEBUG */
128                 goto next1;
129             }
130 
131             Check_dio=1;
132             if(ckt->CKTmode & MODEINITSMSIG) {
133                 vd= *(ckt->CKTstate0 + here->DIOvoltage);
134                 delTemp = *(ckt->CKTstate0 + here->DIOdeltemp);
135             } else if (ckt->CKTmode & MODEINITTRAN) {
136                 vd= *(ckt->CKTstate1 + here->DIOvoltage);
137                 delTemp = *(ckt->CKTstate1 + here->DIOdeltemp);
138             } else if ( (ckt->CKTmode & MODEINITJCT) &&
139                     (ckt->CKTmode & MODETRANOP) && (ckt->CKTmode & MODEUIC) ) {
140                 vd=here->DIOinitCond;
141             } else if ( (ckt->CKTmode & MODEINITJCT) && here->DIOoff) {
142                 vd=0;
143                 delTemp = 0.0;
144             } else if ( ckt->CKTmode & MODEINITJCT) {
145                 vd=here->DIOtVcrit;
146                 delTemp = 0.0;
147             } else if ( ckt->CKTmode & MODEINITFIX && here->DIOoff) {
148                 vd=0;
149                 delTemp = 0.0;
150             } else {
151 #ifndef PREDICTOR
152                 if (ckt->CKTmode & MODEINITPRED) {
153                     *(ckt->CKTstate0 + here->DIOvoltage) =
154                             *(ckt->CKTstate1 + here->DIOvoltage);
155                     vd = DEVpred(ckt,here->DIOvoltage);
156                     *(ckt->CKTstate0 + here->DIOcurrent) =
157                             *(ckt->CKTstate1 + here->DIOcurrent);
158                     *(ckt->CKTstate0 + here->DIOconduct) =
159                             *(ckt->CKTstate1 + here->DIOconduct);
160                     *(ckt->CKTstate0 + here->DIOdeltemp) =
161                             *(ckt->CKTstate1 + here->DIOdeltemp);
162                     delTemp = DEVpred(ckt,here->DIOdeltemp);
163                     *(ckt->CKTstate0 + here->DIOdIdio_dT) =
164                             *(ckt->CKTstate1 + here->DIOdIdio_dT);
165                     *(ckt->CKTstate0+here->DIOqth) =
166                             *(ckt->CKTstate1+here->DIOqth);
167                 } else {
168 #endif /* PREDICTOR */
169                     vd = *(ckt->CKTrhsOld+here->DIOposPrimeNode)-
170                             *(ckt->CKTrhsOld + here->DIOnegNode);
171                     if (selfheat)
172                         delTemp = *(ckt->CKTrhsOld + here->DIOtempNode);
173                     else
174                         delTemp = 0.0;
175                     *(ckt->CKTstate0+here->DIOqth) = model->DIOcth0 * delTemp;
176                     if((ckt->CKTmode & MODEINITTRAN)) {
177                         *(ckt->CKTstate1+here->DIOqth) =
178                             *(ckt->CKTstate0+here->DIOqth);
179                     }
180 #ifndef PREDICTOR
181                 }
182 #endif /* PREDICTOR */
183                 delvd=vd- *(ckt->CKTstate0 + here->DIOvoltage);
184                 deldelTemp = delTemp - *(ckt->CKTstate0 + here->DIOdeltemp);
185                 cdhat= *(ckt->CKTstate0 + here->DIOcurrent) +
186                         *(ckt->CKTstate0 + here->DIOconduct) * delvd +
187                         *(ckt->CKTstate0 + here->DIOdIdio_dT) * deldelTemp;
188                 /*
189                  *   bypass if solution has not changed
190                  */
191 #ifndef NOBYPASS
192                 if ((!(ckt->CKTmode & MODEINITPRED)) && (ckt->CKTbypass)) {
193                     tol=ckt->CKTvoltTol + ckt->CKTreltol*
194                         MAX(fabs(vd),fabs(*(ckt->CKTstate0 +here->DIOvoltage)));
195                     if (fabs(delvd) < tol){
196                         tol=ckt->CKTreltol* MAX(fabs(cdhat),
197                                 fabs(*(ckt->CKTstate0 + here->DIOcurrent)))+
198                                 ckt->CKTabstol;
199                         if (fabs(cdhat- *(ckt->CKTstate0 + here->DIOcurrent))
200                                 < tol) {
201                             if ((here->DIOtempNode == 0) ||
202                                 (fabs(deldelTemp) < (ckt->CKTreltol * MAX(fabs(delTemp),
203                                       fabs(*(ckt->CKTstate0+here->DIOdeltemp)))+
204                                       ckt->CKTvoltTol*1e4))) {
205                                 vd= *(ckt->CKTstate0 + here->DIOvoltage);
206                                 cd= *(ckt->CKTstate0 + here->DIOcurrent);
207                                 gd= *(ckt->CKTstate0 + here->DIOconduct);
208                                 delTemp = *(ckt->CKTstate0 + here->DIOdeltemp);
209                                 dIdio_dT= *(ckt->CKTstate0 + here->DIOdIdio_dT);
210                                 goto load;
211                             }
212                         }
213                     }
214                 }
215 #endif /* NOBYPASS */
216                 /*
217                  *   limit new junction voltage
218                  */
219                 if ( (model->DIObreakdownVoltageGiven) &&
220                         (vd < MIN(0,-here->DIOtBrkdwnV+10*vtebrk))) {
221                     vdtemp = -(vd+here->DIOtBrkdwnV);
222                     vdtemp = DEVpnjlim(vdtemp,
223                             -(*(ckt->CKTstate0 + here->DIOvoltage) +
224                             here->DIOtBrkdwnV),vtebrk,
225                             here->DIOtVcrit,&Check_dio);
226                     vd = -(vdtemp+here->DIOtBrkdwnV);
227                 } else {
228                     vd = DEVpnjlim(vd,*(ckt->CKTstate0 + here->DIOvoltage),
229                             vte,here->DIOtVcrit,&Check_dio);
230                 }
231                 if (selfheat)
232                     delTemp = DEVlimitlog(delTemp,
233                         *(ckt->CKTstate0 + here->DIOdeltemp), 100, &Check_th);
234                 else
235                     delTemp = 0.0;
236             }
237             /*
238              *   compute dc current and derivitives
239              */
240 next1:
241             if (selfheat) {
242                 Temp = here->DIOtemp + delTemp;
243                 DIOtempUpdate(model, here, Temp, ckt);
244             } else {
245                 Temp = here->DIOtemp;
246             }
247 
248             csat = here->DIOtSatCur;
249             csat_dT = here->DIOtSatCur_dT;
250             csatsw = here->DIOtSatSWCur;
251             csatsw_dT = here->DIOtSatSWCur_dT;
252             gspr = here->DIOtConductance;
253             vt = CONSTKoverQ * Temp;
254             vte = model->DIOemissionCoeff * vt;
255             vtebrk = model->DIObrkdEmissionCoeff * vt;
256 
257             if (model->DIOsatSWCurGiven) {              /* sidewall current */
258 
259                 if (model->DIOswEmissionCoeffGiven) {   /* current with own characteristic */
260 
261                     vtesw = model->DIOswEmissionCoeff * vt;
262 
263                     if (vd >= -3*vtesw) {               /* forward */
264 
265                         evd = exp(vd/vtesw);
266                         cdsw = csatsw*(evd-1);
267                         gdsw = csatsw*evd/vtesw;
268                         cdsw_dT = csatsw_dT * (evd - 1) - csatsw * vd * evd / (vtesw * Temp);
269 
270                     } else if((!(model->DIObreakdownVoltageGiven)) ||
271                             vd >= -here->DIOtBrkdwnV) { /* reverse */
272 
273                         argsw = 3*vtesw/(vd*CONSTe);
274                         argsw = argsw * argsw * argsw;
275                         argsw_dT = 3 * argsw / Temp;
276                         cdsw = -csatsw*(1+argsw);
277                         gdsw = csatsw*3*argsw/vd;
278                         cdsw_dT = -csatsw_dT - (csatsw_dT*argsw + csatsw*argsw_dT);
279 
280                     } else {                            /* breakdown */
281                         double evrev_dT;
282 
283                         evrev = exp(-(here->DIOtBrkdwnV+vd)/vtebrk);
284                         evrev_dT = (here->DIOtBrkdwnV+vd)*evrev/(vtebrk*Temp);
285                         cdsw = -csatsw*evrev;
286                         gdsw = csatsw*evrev/vtebrk;
287                         cdsw_dT = -(csatsw_dT*evrev + csatsw*evrev_dT);
288 
289                     }
290 
291                 } else { /* merge saturation currents and use same characteristic as bottom diode */
292 
293                     csat = csat + csatsw;
294                     csat_dT = csat_dT + csatsw_dT;
295                     cdsw_dT = 0.0;
296 
297                 }
298 
299             }
300 
301             /*
302              *   temperature dependent diode saturation current and derivative
303              */
304 
305             if (vd >= -3*vte) {                 /* bottom current forward */
306 
307                 evd = exp(vd/vte);
308                 cdb = csat*(evd-1);
309                 gdb = csat*evd/vte;
310                 cdb_dT = csat_dT * (evd - 1) - csat * vd * evd / (vte * Temp);
311                 if (model->DIOrecSatCurGiven) { /* recombination current */
312                     double vterec = model->DIOrecEmissionCoeff*vt;
313                     evd_rec = exp(vd/(vterec));
314                     cdb_rec = here->DIOtRecSatCur*(evd_rec-1);
315                     gdb_rec = here->DIOtRecSatCur*evd_rec/vterec;
316                     cdb_rec_dT = here->DIOtRecSatCur_dT * (evd_rec - 1)
317                                 -here->DIOtRecSatCur * vd * evd_rec / (vterec*Temp);
318                     t1 = pow((1-vd/here->DIOtJctPot), 2) + 0.005;
319                     gen_fac = pow(t1, here->DIOtGradingCoeff/2);
320                     gen_fac_vd = -here->DIOtGradingCoeff * (1-vd/here->DIOtJctPot)
321                                                          * pow(t1, (here->DIOtGradingCoeff/2-1));
322                     cdb_rec = cdb_rec * gen_fac;
323                     gdb_rec = gdb_rec * gen_fac + cdb_rec * gen_fac_vd;
324                     cdb = cdb + cdb_rec;
325                     gdb = gdb + gdb_rec;
326                     cdb_dT = cdb_dT + cdb_rec_dT*gen_fac;
327                 }
328 
329             } else if((!(model->DIObreakdownVoltageGiven)) ||
330                     vd >= -here->DIOtBrkdwnV) { /* reverse */
331 
332                 double darg_dT;
333 
334                 arg = 3*vte/(vd*CONSTe);
335                 arg = arg * arg * arg;
336                 darg_dT = 3 * arg / Temp;
337                 cdb = -csat*(1+arg);
338                 gdb = csat*3*arg/vd;
339                 cdb_dT = -csat_dT - (csat_dT*arg + csat*darg_dT);
340 
341             } else {                            /* breakdown */
342                 double evrev_dT;
343 
344                 evrev = exp(-(here->DIOtBrkdwnV+vd)/vtebrk);
345                 evrev_dT = (here->DIOtBrkdwnV+vd)*evrev/(vtebrk*Temp);
346                 cdb = -csat*evrev;
347                 gdb = csat*evrev/vtebrk;
348                 cdb_dT = -(csat_dT*evrev + csat*evrev_dT);
349 
350             }
351 
352             if (model->DIOtunSatSWCurGiven) {    /* tunnel sidewall current */
353 
354                 vtetun = model->DIOtunEmissionCoeff * vt;
355                 evd = exp(-vd/vtetun);
356 
357                 cdsw = cdsw - here->DIOtTunSatSWCur * (evd - 1);
358                 gdsw = gdsw + here->DIOtTunSatSWCur * evd / vtetun;
359                 cdsw_dT = cdsw_dT - here->DIOtTunSatSWCur_dT * (evd - 1)
360                                   - here->DIOtTunSatSWCur * vd * evd / (vtetun * Temp);
361 
362             }
363 
364             if (model->DIOtunSatCurGiven) {      /* tunnel bottom current */
365 
366                 vtetun = model->DIOtunEmissionCoeff * vt;
367                 evd = exp(-vd/vtetun);
368 
369                 cdb = cdb - here->DIOtTunSatCur * (evd - 1);
370                 gdb = gdb + here->DIOtTunSatCur * evd / vtetun;
371                 cdb_dT = cdb_dT - here->DIOtTunSatCur_dT * (evd - 1)
372                                 - here->DIOtTunSatCur * vd * evd / (vtetun * Temp);
373 
374             }
375 
376             cd = cdb + cdsw;
377             gd = gdb + gdsw;
378             dIdio_dT = cdb_dT + cdsw_dT;
379 
380             if (vd >= -3*vte) { /* limit forward */
381 
382                 if( (model->DIOforwardKneeCurrent > 0.0) && (cd > 1.0e-18) ) {
383                     ikf_area_m = here->DIOforwardKneeCurrent;
384                     sqrt_ikf = sqrt(cd/ikf_area_m);
385                     gd = ((1+sqrt_ikf)*gd - cd*gd/(2*sqrt_ikf*ikf_area_m))/(1+2*sqrt_ikf + cd/ikf_area_m) + ckt->CKTgmin;
386                     cd = cd/(1+sqrt_ikf) + ckt->CKTgmin*vd;
387                 } else {
388                     gd = gd + ckt->CKTgmin;
389                     cd = cd + ckt->CKTgmin*vd;
390                 }
391 
392             } else {            /* limit reverse */
393 
394                 if( (model->DIOreverseKneeCurrent > 0.0) && (cd < -1.0e-18) ) {
395                     ikr_area_m = here->DIOreverseKneeCurrent;
396                     sqrt_ikr = sqrt(cd/(-ikr_area_m));
397                     gd = ((1+sqrt_ikr)*gd + cd*gd/(2*sqrt_ikr*ikr_area_m))/(1+2*sqrt_ikr - cd/ikr_area_m) + ckt->CKTgmin;
398                     cd = cd/(1+sqrt_ikr) + ckt->CKTgmin*vd;
399                 } else {
400                     gd = gd + ckt->CKTgmin;
401                     cd = cd + ckt->CKTgmin*vd;
402                 }
403 
404             }
405 
406             if ((ckt->CKTmode & (MODEDCTRANCURVE | MODETRAN | MODEAC | MODEINITSMSIG)) ||
407                      ((ckt->CKTmode & MODETRANOP) && (ckt->CKTmode & MODEUIC))) {
408               /*
409                *   charge storage elements
410                */
411                 czero=here->DIOtJctCap;
412                 if (vd < here->DIOtDepCap){
413                     arg=1-vd/here->DIOtJctPot;
414                     sarg=exp(-here->DIOtGradingCoeff*log(arg));
415                     deplcharge = here->DIOtJctPot*czero*(1-arg*sarg)/(1-here->DIOtGradingCoeff);
416                     deplcap = czero*sarg;
417                 } else {
418                     czof2=czero/here->DIOtF2;
419                     deplcharge = czero*here->DIOtF1+czof2*(here->DIOtF3*(vd-here->DIOtDepCap)+
420                                  (here->DIOtGradingCoeff/(here->DIOtJctPot+here->DIOtJctPot))*(vd*vd-here->DIOtDepCap*here->DIOtDepCap));
421                     deplcap = czof2*(here->DIOtF3+here->DIOtGradingCoeff*vd/here->DIOtJctPot);
422                 }
423                 czeroSW=here->DIOtJctSWCap;
424                 if (vd < here->DIOtDepSWCap){
425                     argSW=1-vd/here->DIOtJctSWPot;
426                     sargSW=exp(-model->DIOgradingSWCoeff*log(argSW));
427                     deplchargeSW = here->DIOtJctSWPot*czeroSW*(1-argSW*sargSW)/(1-model->DIOgradingSWCoeff);
428                     deplcapSW = czeroSW*sargSW;
429                 } else {
430                     czof2SW=czeroSW/here->DIOtF2SW;
431                     deplchargeSW = czeroSW*here->DIOtF1+czof2SW*(here->DIOtF3SW*(vd-here->DIOtDepSWCap)+
432                                    (model->DIOgradingSWCoeff/(here->DIOtJctSWPot+here->DIOtJctSWPot))*(vd*vd-here->DIOtDepSWCap*here->DIOtDepSWCap));
433                     deplcapSW = czof2SW*(here->DIOtF3SW+model->DIOgradingSWCoeff*vd/here->DIOtJctSWPot);
434                 }
435 
436                 diffcharge = here->DIOtTransitTime*cdb;
437                 diffchargeSW = here->DIOtTransitTime*cdsw;
438                 *(ckt->CKTstate0 + here->DIOcapCharge) =
439                         diffcharge + diffchargeSW + deplcharge + deplchargeSW;
440 
441                 diffcap = here->DIOtTransitTime*gdb;
442                 diffcapSW = here->DIOtTransitTime*gdsw;
443 
444                 capd = diffcap + diffcapSW + deplcap + deplcapSW + here->DIOcmetal + here->DIOcpoly;
445 
446                 here->DIOcap = capd;
447 
448                 /*
449                  *   store small-signal parameters
450                  */
451                 if( (!(ckt->CKTmode & MODETRANOP)) ||
452                         (!(ckt->CKTmode & MODEUIC)) ) {
453                     if (ckt->CKTmode & MODEINITSMSIG){
454                         *(ckt->CKTstate0 + here->DIOcapCurrent) = capd;
455 
456                         if(SenCond){
457                             *(ckt->CKTstate0 + here->DIOcurrent) = cd;
458                             *(ckt->CKTstate0 + here->DIOconduct) = gd;
459                             *(ckt->CKTstate0 + here->DIOdIdio_dT) = dIdio_dT;
460 #ifdef SENSDEBUG
461                             printf("storing small signal parameters\n");
462                             printf("cd = %.7e,vd = %.7e\n",cd,vd);
463                             printf("capd = %.7e ,gd = %.7e \n",capd,gd);
464 #endif /* SENSDEBUG */
465                         }
466                         continue;
467                     }
468 
469                     /*
470                      *   transient analysis
471                      */
472                     if(SenCond && (ckt->CKTsenInfo->SENmode == TRANSEN)){
473                         *(ckt->CKTstate0 + here->DIOcurrent) = cd;
474 #ifdef SENSDEBUG
475                         printf("storing parameters for transient sensitivity\n"
476                                 );
477                         printf("qd = %.7e, capd = %.7e,cd = %.7e\n",
478                                 *(ckt->CKTstate0 + here->DIOcapCharge),capd,cd);
479 #endif /* SENSDEBUG */
480                         continue;
481                     }
482 
483                     if (ckt->CKTmode & MODEINITTRAN) {
484                         *(ckt->CKTstate1 + here->DIOcapCharge) =
485                                 *(ckt->CKTstate0 + here->DIOcapCharge);
486                     }
487                     error = NIintegrate(ckt,&geq,&ceq,capd,here->DIOcapCharge);
488                     if(error) return(error);
489                     gd=gd+geq;
490                     cd=cd+*(ckt->CKTstate0 + here->DIOcapCurrent);
491                     if (ckt->CKTmode & MODEINITTRAN) {
492                         *(ckt->CKTstate1 + here->DIOcapCurrent) =
493                                 *(ckt->CKTstate0 + here->DIOcapCurrent);
494                     }
495                     if (selfheat)
496                     {
497                         error = NIintegrate(ckt, &gcTt, &ceqqth, model->DIOcth0, here->DIOqth);
498                         if (error) return(error);
499                         if (ckt->CKTmode & MODEINITTRAN) {
500                             *(ckt->CKTstate1 + here->DIOcqth) =
501                                     *(ckt->CKTstate0 + here->DIOcqth);
502                         }
503                     }
504                 }
505             }
506 
507             if(SenCond) goto next2;
508 
509             /*
510              *   check convergence
511              */
512             if ( (!(ckt->CKTmode & MODEINITFIX)) || (!(here->DIOoff))  ) {
513                 if ((Check_th == 1) || (Check_dio == 1)) {
514                     ckt->CKTnoncon++;
515                     ckt->CKTtroubleElt = (GENinstance *) here;
516                 }
517             }
518 next2:      *(ckt->CKTstate0 + here->DIOvoltage) = vd;
519             *(ckt->CKTstate0 + here->DIOcurrent) = cd;
520             *(ckt->CKTstate0 + here->DIOconduct) = gd;
521             *(ckt->CKTstate0 + here->DIOdeltemp) = delTemp;
522             *(ckt->CKTstate0 + here->DIOdIdio_dT) = dIdio_dT;
523 
524             if(SenCond)  continue;
525 
526 #ifndef NOBYPASS
527             load:
528 #endif
529             if (selfheat) {
530                 double dIrs_dVrs, dIrs_dgspr, dIth_dIrs;
531                 vrs = *(ckt->CKTrhsOld + here->DIOposNode) - *(ckt->CKTrhsOld + here->DIOposPrimeNode);
532                 Ith = vd*cd + vrs*vrs*gspr; /* Diode dissipated power */
533                 dIrs_dVrs = gspr;
534                 dIrs_dgspr = vrs;
535                 dIrs_dT = dIrs_dgspr * here->DIOtConductance_dT;
536                 dIth_dVrs = vrs*gspr;
537                 dIth_dIrs = vrs;
538                 dIth_dVrs = dIth_dVrs + dIth_dIrs*dIrs_dVrs;
539                 dIth_dT = dIth_dIrs*dIrs_dT + dIdio_dT*vd;
540                 dIth_dVdio = cd + vd*gd;
541                 here->DIOdIth_dVrs = dIth_dVrs;
542                 here->DIOdIth_dVdio = dIth_dVdio;
543                 here->DIOdIth_dT = dIth_dT;
544                 here->DIOgcTt = gcTt;
545                 here->DIOdIrs_dT = dIrs_dT;
546             }
547             /*
548              *   load current vector
549              */
550             cdeq=cd-gd*vd;
551             *(ckt->CKTrhs + here->DIOnegNode) += cdeq;
552             *(ckt->CKTrhs + here->DIOposPrimeNode) -= cdeq;
553             if (selfheat) {
554                 *(ckt->CKTrhs + here->DIOposNode)      +=  dIrs_dT*delTemp;
555                 *(ckt->CKTrhs + here->DIOposPrimeNode) +=  dIdio_dT*delTemp - dIrs_dT*delTemp;
556                 *(ckt->CKTrhs + here->DIOnegNode)      += -dIdio_dT*delTemp;
557                 *(ckt->CKTrhs + here->DIOtempNode)     +=  Ith - dIth_dVdio*vd - dIth_dVrs*vrs - dIth_dT*delTemp - ceqqth;
558             }
559             /*
560              *   load matrix
561              */
562             *(here->DIOposPosPtr) += gspr;
563             *(here->DIOnegNegPtr) += gd;
564             *(here->DIOposPrimePosPrimePtr) += (gd + gspr);
565             *(here->DIOposPosPrimePtr) -= gspr;
566             *(here->DIOnegPosPrimePtr) -= gd;
567             *(here->DIOposPrimePosPtr) -= gspr;
568             *(here->DIOposPrimeNegPtr) -= gd;
569             if (selfheat) {
570                 (*(here->DIOtempPosPtr)      += -dIth_dVrs);
571                 (*(here->DIOtempPosPrimePtr) += -dIth_dVdio + dIth_dVrs);
572                 (*(here->DIOtempNegPtr)      +=  dIth_dVdio);
573                 (*(here->DIOtempTempPtr)     += -dIth_dT + 1/model->DIOrth0 + gcTt);
574                 (*(here->DIOposTempPtr)      +=  dIrs_dT);
575                 (*(here->DIOposPrimeTempPtr) +=  dIdio_dT - dIrs_dT);
576                 (*(here->DIOnegTempPtr)      += -dIdio_dT);
577             }
578         }
579     }
580     return(OK);
581 }
582