Xyce  6.1
N_DEV_VDMOS.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------------
2 // Copyright Notice
3 //
4 // Copyright 2002 Sandia Corporation. Under the terms
5 // of Contract DE-AC04-94AL85000 with Sandia Corporation, the U.S.
6 // Government retains certain rights in this software.
7 //
8 // Xyce(TM) Parallel Electrical Simulator
9 // Copyright (C) 2002-2015 Sandia Corporation
10 //
11 // This program is free software: you can redistribute it and/or modify
12 // it under the terms of the GNU General Public License as published by
13 // the Free Software Foundation, either version 3 of the License, or
14 // (at your option) any later version.
15 //
16 // This program is distributed in the hope that it will be useful,
17 // but WITHOUT ANY WARRANTY; without even the implied warranty of
18 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 // GNU General Public License for more details.
20 //
21 // You should have received a copy of the GNU General Public License
22 // along with this program. If not, see <http://www.gnu.org/licenses/>.
23 //-----------------------------------------------------------------------------
24 
25 //-------------------------------------------------------------------------
26 // Filename : $RCSfile: N_DEV_VDMOS.C,v $
27 //
28 // Purpose : Implement the Vertical double-diffused power MOSFET
29 //
30 // Special Notes :
31 //
32 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
33 //
34 // Creation Date : 02/28/00
35 //
36 // Revision Information:
37 // ---------------------
38 //
39 // Revision Number: $Revision: 1.172.2.1 $
40 //
41 // Revision Date : $Date: 2015/04/02 18:20:11 $
42 //
43 // Current Owner : $Author: tvrusso $
44 //-------------------------------------------------------------------------
45 #include <Xyce_config.h>
46 
47 // ---------- Standard Includes ----------
48 #include <N_UTL_Math.h>
49 #include <climits>
50 
51 // ---------- Xyce Includes ----------
52 #include <N_DEV_Const.h>
53 #include <N_DEV_DeviceOptions.h>
54 #include <N_DEV_ExternData.h>
55 #include <N_DEV_MOSFET1.h>
56 #include <N_DEV_MatrixLoadData.h>
57 #include <N_DEV_Message.h>
58 #include <N_DEV_SolverState.h>
59 #include <N_DEV_VDMOS.h>
60 #include <N_ERH_ErrorMgr.h>
61 #include <N_LAS_Matrix.h>
62 #include <N_LAS_Vector.h>
63 #include <N_UTL_FeatureTest.h>
64 #include <N_UTL_MachDepParams.h>
65 
66 namespace Xyce {
67 namespace Device {
68 
69 
70 namespace VDMOS {
71 
72 
74 {
75 // Set up map for normal (double) param variables:
76  p.addPar ("L",0.0,&VDMOS::Instance::l)
77  .setUnit(U_METER)
78  .setCategory(CAT_GEOMETRY)
79  .setDescription("Channel length");
80 
81  p.addPar ("W",0.0,&VDMOS::Instance::w)
82  .setUnit(U_METER)
83  .setCategory(CAT_GEOMETRY)
84  .setDescription("Channel width");
85 
87  .setUnit(U_METER2)
88  .setCategory(CAT_GEOMETRY)
89  .setDescription("Drain diffusion area");
90 
92  .setUnit(U_METER2)
93  .setCategory(CAT_GEOMETRY)
94  .setDescription("Source diffusion area");
95 
97  .setUnit(U_SQUARES)
98  .setCategory(CAT_GEOMETRY)
99  .setDescription("Multiplier for RSH to yield parasitic resistance of drain");
100 
102  .setUnit(U_SQUARES)
103  .setCategory(CAT_GEOMETRY)
104  .setDescription("Multiplier for RSH to yield parasitic resistance of source");
105 
107  .setUnit(U_METER)
108  .setCategory(CAT_GEOMETRY)
109  .setDescription("Drain diffusion perimeter");
110 
112  .setUnit(U_METER)
113  .setCategory(CAT_GEOMETRY)
114  .setDescription("Source diffusion perimeter");
115 
117  .setUnit(U_NONE)
118  .setCategory(CAT_CONTROL)
119  .setDescription("Multiplier for M devices connected in parallel");
120 
121  p.addPar ("TEMP",0.0,&VDMOS::Instance::temp)
122  .setExpressionAccess(ParameterType::TIME_DEP)
123  .setUnit(STANDARD)
124  .setCategory(CAT_NONE)
125  .setDescription("Device temperature");
126 }
127 
129 {
130  p.addPar ("L0",0.0,&VDMOS::Model::l0)
131  .setUnit(U_METER)
132  .setCategory(CAT_GEOMETRY)
133  .setDescription("Gate length of nominal device");
134 
135  p.addPar ("W0",0.0,&VDMOS::Model::w0)
136  .setUnit(U_METER)
137  .setCategory(CAT_GEOMETRY)
138  .setDescription("Gate width of nominal device");
139 
140  p.addPar ("VTO",0.0,&VDMOS::Model::vt0)
141  .setUnit(U_VOLT)
142  .setCategory(CAT_VOLT)
143  .setDescription("Zero-bias threshold voltage");
144 
145  p.addPar ("PHI",0.6,&VDMOS::Model::phi)
146  .setUnit(U_VOLT)
147  .setCategory(CAT_PROCESS)
148  .setDescription("Surface potential");
149 
151  .setExpressionAccess(ParameterType::MIN_RES)
152  .setUnit(U_OHM)
153  .setCategory(CAT_RES)
154  .setDescription("Drain ohmic resistance");
155 
157  .setUnit(U_OHM)
158  .setCategory(CAT_RES)
159  .setDescription("Gate ohmic resistance");
160 
162  .setExpressionAccess(ParameterType::MIN_RES)
163  .setUnit(U_OHM)
164  .setCategory(CAT_RES)
165  .setDescription("Source ohmic resistance");
166 
167  p.addPar ("CBD",0.0,&VDMOS::Model::capBD)
168  .setGivenMember(&VDMOS::Model::capBDGiven)
169  .setExpressionAccess(ParameterType::MIN_CAP)
170  .setUnit(U_FARAD)
171  .setCategory(CAT_CAP)
172  .setDescription("Zero-bias bulk-drain p-n capacitance");
173 
174  p.addPar ("CBS",0.0,&VDMOS::Model::capBS)
175  .setGivenMember(&VDMOS::Model::capBSGiven)
176  .setExpressionAccess(ParameterType::MIN_CAP)
177  .setUnit(U_FARAD)
178  .setCategory(CAT_CAP)
179  .setDescription("Zero-bias bulk-source p-n capacitance");
180 
181  p.addPar ("TS",0.0,&VDMOS::Model::timeScale)
182  .setUnit(U_NONE)
184  .setDescription("");
185 
186  p.addPar ("IS",1e-14,&VDMOS::Model::jctSatCur)
187  .setUnit(U_AMP)
188  .setCategory(CAT_CURRENT)
189  .setDescription("Bulk p-n saturation current");
190 
192  .setUnit(U_VOLT)
193  .setCategory(CAT_VOLT)
194  .setDescription("Bulk p-n bottom potential");
195 
197  .setUnit(U_FARADMM1)
198  .setCategory(CAT_CAP)
199  .setDescription("Gate-source overlap capacitance/channel width");
200 
202  .setUnit(U_FARADMM1)
203  .setCategory(CAT_CAP)
204  .setDescription("Gate-drain overlap capacitance/channel width");
205 
207  .setUnit(U_FARADMM1)
208  .setCategory(CAT_CAP)
209  .setDescription("Gate-bulk overlap capacitance/channel length");
210 
212  .setUnit(U_OHM)
213  .setCategory(CAT_RES)
214  .setDescription("Drain,source diffusion sheet resistance");
215 
216  p.addPar ("CJ",0.0,&VDMOS::Model::bulkCapFactor)
217  .setGivenMember(&VDMOS::Model::bulkCapFactorGiven)
218  .setUnit(U_FARADMM2)
219  .setCategory(CAT_CAP)
220  .setDescription("Bulk p-n zero-bias bottom capacitance/area");
221 
223  .setUnit(U_NONE)
224  .setCategory(CAT_DOPING)
225  .setDescription("Bulk p-n bottom grading coefficient");
226 
228  .setGivenMember(&VDMOS::Model::sideWallCapFactorGiven)
229  .setUnit(U_FARADMM2)
230  .setCategory(CAT_CAP)
231  .setDescription("Bulk p-n zero-bias sidewall capacitance/area");
232 
234  .setUnit(U_NONE)
235  .setCategory(CAT_DOPING)
236  .setDescription("Bulk p-n sidewall grading coefficient");
237 
239  .setUnit(U_AMPMM2)
240  .setCategory(CAT_PROCESS)
241  .setDescription("Bulk p-n saturation current density");
242 
243  p.addPar ("TOX",1e-7,&VDMOS::Model::oxideThickness)
244  .setUnit(U_METER)
245  .setCategory(CAT_GEOMETRY)
246  .setDescription("Gate oxide thickness");
247 
248  p.addPar ("LD",0.0,&VDMOS::Model::latDiff)
249  .setUnit(U_METER)
250  .setCategory(CAT_DOPING)
251  .setDescription("Lateral diffusion length");
252 
253  p.addPar ("UO",280., &VDMOS::Model::surfaceMobility)
254  .setUnit(U_CMM2VM1SM1)
256  .setDescription("Surface mobility");
257 
258  p.addPar ("U0",280., &VDMOS::Model::surfaceMobility0)
259  .setUnit(U_CMM2VM1SM1)
260  .setCategory(CAT_PROCESS)
261  .setDescription("Surface mobility");
262 
264  .setUnit(U_NONE)
265  .setCategory(CAT_CAP)
266  .setDescription("Coefficient for forward-bias depletion capacitance formula");
267 
268  p.addPar ("NSUB",0.0,&VDMOS::Model::substrateDoping)
269  .setUnit(U_CMM3)
270  .setCategory(CAT_DOPING)
271  .setDescription("Substrate doping density");
272 
274  .setUnit(U_CMM2)
275  .setCategory(CAT_PROCESS)
276  .setDescription("Surface state density");
277 
278  p.addPar ("TNOM",0.0,&VDMOS::Model::tnom)
279  .setUnit(STANDARD)
280  .setCategory(CAT_NONE)
281  .setDescription("Parameter measurement temperature");
282 
283  p.addPar ("VMAX",4e+4,&VDMOS::Model::maxDriftVel)
284  .setUnit(U_MSM1)
285  .setCategory(CAT_PROCESS)
286  .setDescription("Maximum drift velocity for carriers");
287 
288  p.addPar ("XJ",0.0,&VDMOS::Model::junctionDepth)
289  .setUnit(U_METER)
290  .setCategory(CAT_GEOMETRY)
291  .setDescription("Metallurgical junction depth");
292 
293  p.addPar ("LAMBDA",0.048,&VDMOS::Model::lambda)
294  .setUnit(U_VOLTM1)
295  .setCategory(CAT_PROCESS)
296  .setDescription("Output conductance parameter");
297 
298  p.addPar ("DELTA",5.0,&VDMOS::Model::delta)
299  .setUnit(U_NONE)
300  .setCategory(CAT_NONE)
301  .setDescription("Transition width parameter");
302 
303  p.addPar ("ETA",1.32,&VDMOS::Model::eta)
304  .setUnit(U_NONE)
305  .setCategory(CAT_NONE)
306  .setDescription("Subthreshold ideality factor");
307 
308  p.addPar ("M",4.0,&VDMOS::Model::m)
309  .setUnit(U_NONE)
310  .setCategory(CAT_CONTROL)
311  .setDescription("Knee shape parameter");
312 
313  p.addPar ("MC",3.0,&VDMOS::Model::mc)
314  .setUnit(U_NONE)
316  .setDescription("");
317 
318  p.addPar ("SIGMA0",0.048,&VDMOS::Model::sigma0)
319  .setUnit(U_NONE)
320  .setCategory(CAT_NONE)
321  .setDescription("DIBL parameter");
322 
323  p.addPar ("VSIGMAT",1.7,&VDMOS::Model::vsigmat)
324  .setUnit(U_VOLT)
325  .setCategory(CAT_VOLT)
326  .setDescription("DIBL parameter");
327 
328  p.addPar ("VSIGMA",0.2,&VDMOS::Model::vsigma)
329  .setUnit(U_VOLT)
330  .setCategory(CAT_VOLT)
331  .setDescription("DIBL parameter");
332 
333  p.addPar ("THETA",0.0,&VDMOS::Model::theta)
334  .setUnit(U_MVM1)
335  .setCategory(CAT_PROCESS)
336  .setDescription("Mobility degradation parameter");
337 
338  p.addPar ("GAMMAS0",0.5,&VDMOS::Model::gammas0)
339  .setUnit(U_VOLTMH)
340  .setCategory(CAT_NONE)
341  .setDescription("Body effect constant in front of square root term");
342 
343  p.addPar ("GAMMAL0",0.0,&VDMOS::Model::gammal0)
344  .setUnit(U_NONE)
345  .setCategory(CAT_NONE)
346  .setDescription("Body effect constant in front of linear term");
347 
348  p.addPar ("LGAMMAS",0.0,&VDMOS::Model::lgammas)
349  .setUnit(U_VOLTMH)
350  .setCategory(CAT_NONE)
351  .setDescription("Sensitivity of gS on device length");
352 
353  p.addPar ("WGAMMAS",0.0,&VDMOS::Model::wgammas)
354  .setUnit(U_VOLTMH)
355  .setCategory(CAT_NONE)
356  .setDescription("Sensitivity of gS on device width");
357 
358  p.addPar ("LGAMMAL",0.0,&VDMOS::Model::lgammal_)
359  .setUnit(U_NONE)
360  .setCategory(CAT_NONE)
361  .setDescription("Sensitivity of gL on device length");
362 
363  p.addPar ("WGAMMAL",0.0,&VDMOS::Model::wgammal)
364  .setUnit(U_NONE)
365  .setCategory(CAT_NONE)
366  .setDescription("Sensitivity of gL on device width");
367 
368  p.addPar ("K",0.0,&VDMOS::Model::k)
369  .setUnit(U_NONE)
371  .setDescription("");
372 
373  p.addPar ("KVT",0.0,&VDMOS::Model::kvt)
374  .setUnit(U_NONE)
376  .setDescription("");
377 
378  p.addPar ("MDTEMP",0.0,&VDMOS::Model::mdtemp)
379  .setUnit(U_NONE)
381  .setDescription("");
382 
383  p.addPar ("KVS",0.0,&VDMOS::Model::kvs)
384  .setUnit(U_NONE)
386  .setDescription("");
387 
388  p.addPar ("TVS",0.0,&VDMOS::Model::tvs)
389  .setUnit(U_NONE)
391  .setDescription("");
392 
393  p.addPar ("MTH",0.0,&VDMOS::Model::mth)
394  .setUnit(U_NONE)
396  .setDescription("");
397 
398  p.addPar ("ARTD",0.0,&VDMOS::Model::artd)
399  .setUnit(U_NONE)
401  .setDescription("");
402 
403  p.addPar ("BRTD",0.035,&VDMOS::Model::brtd)
404  .setUnit(U_NONE)
406  .setDescription("");
407 
408  p.addPar ("CRTD",0.1472,&VDMOS::Model::crtd)
409  .setUnit(U_NONE)
411  .setDescription("");
412 
413  p.addPar ("DRTD",0.0052,&VDMOS::Model::drtd)
414  .setUnit(U_NONE)
416  .setDescription("");
417 
418  p.addPar ("NRTD",0.115,&VDMOS::Model::nrtd)
419  .setUnit(U_NONE)
421  .setDescription("");
422 
423  p.addPar ("N2",1.0,&VDMOS::Model::n2)
424  .setUnit(U_NONE)
426  .setDescription("");
427 
428  p.addPar ("XQC",0.6,&VDMOS::Model::xqc)
429  .setUnit(U_NONE)
430  .setCategory(CAT_PROCESS)
431  .setDescription("Charge partitioning factor");
432 
433  p.addPar ("MCV",10.0,&VDMOS::Model::mcv)
434  .setUnit(U_NONE)
435  .setCategory(CAT_GEOMETRY)
436  .setDescription("Transition width parameter used by the charge partitioning scheme");
437 
438  p.addPar ("VFB",0.0,&VDMOS::Model::vfb)
439  .setUnit(U_VOLT)
440  .setCategory(CAT_VOLT)
441  .setDescription("Flat band voltage");
442 
443  p.addPar ("ALPHA",0.0,&VDMOS::Model::alpha)
444  .setUnit(U_NONE)
445  .setCategory(CAT_NONE)
446  .setDescription("Parameter accounting for the threshold dependence on the channel potential");
447 
448  p.addPar ("LS",35e-9,&VDMOS::Model::ls)
449  .setUnit(U_NONE)
451  .setDescription("");
452 
453  p.addPar ("RSUB",0.0,&VDMOS::Model::rsub)
454  .setUnit(U_NONE)
456  .setDescription("");
457 
458  p.addPar ("VP",0.0,&VDMOS::Model::vp)
459  .setGivenMember(&VDMOS::Model::vpGiven)
460  .setUnit(U_NONE)
462  .setDescription("");
463 
464  p.addPar ("AI",2e+9,&VDMOS::Model::ai)
465  .setUnit(U_NONE)
467  .setDescription("");
468 
469  p.addPar ("BI",8e+8,&VDMOS::Model::bi)
470  .setUnit(U_NONE)
472  .setDescription("");
473 
474  p.addPar ("DELMAX",0.9,&VDMOS::Model::delmax)
475  .setUnit(U_NONE)
477  .setDescription("");
478 
479  p.addPar ("MD",2.0,&VDMOS::Model::md)
480  .setUnit(U_NONE)
482  .setDescription("");
483 
484  p.addPar ("DRIFTPARAMA",0.08,&VDMOS::Model::driftParamA)
485  .setUnit(U_OHM)
486  .setCategory(CAT_RES)
487  .setDescription("Drift region resistance intercept parameter");
488 
489  p.addPar ("DRIFTPARAMB",0.013,&VDMOS::Model::driftParamB)
490  .setUnit(U_OHMPV)
491  .setCategory(CAT_RES)
492  .setDescription("Drift region resistance slope parameter");
493 
494  p.addPar ("RDSSHUNT",0.0,&VDMOS::Model::rdsshunt)
495  .setUnit(U_OHM)
496  .setCategory(CAT_RES)
497  .setDescription("Drain-source shunt resistance");
498 
499  p.addPar ("D1IS",1e-14,&VDMOS::Model::D1DIOsatCur)
500  .setUnit(U_AMP)
501  .setCategory(CAT_CURRENT)
502  .setDescription("Drain-Source diode saturation current");
503 
504  p.addPar ("D1RS",0.0,&VDMOS::Model::D1DIOresist)
505  .setUnit(U_OHM)
506  .setCategory(CAT_RES)
507  .setDescription("Drain-source diode ohmic resistance");
508 
510  .setUnit(U_NONE)
511  .setCategory(CAT_PROCESS)
512  .setDescription("Drain-source diode emission coefficient");
513 
514  p.addPar ("D1TT",0.0,&VDMOS::Model::D1DIOtransitTime)
515  .setUnit(U_SECOND)
516  .setCategory(CAT_PROCESS)
517  .setDescription("Drain-source diode transit time");
518 
519  p.addPar ("D1CJO",0.0,&VDMOS::Model::D1DIOjunctionCap)
520  .setUnit(U_FARAD)
521  .setCategory(CAT_CAP)
522  .setDescription("Drain-source diode junction capacitance");
523 
524  p.addPar ("D1VJ",1.0,&VDMOS::Model::D1DIOjunctionPot)
525  .setUnit(U_VOLT)
526  .setCategory(CAT_VOLT)
527  .setDescription("Drain-source diode junction potential");
528 
530  .setUnit(U_NONE)
531  .setCategory(CAT_PROCESS)
532  .setDescription("Drain-source diode grading coefficient");
533 
535  .setUnit(U_EV)
536  .setCategory(CAT_PROCESS)
537  .setDescription("Drain-source diode activation energy");
538 
540  .setUnit(U_NONE)
541  .setCategory(CAT_PROCESS)
542  .setDescription("Drain-source diode sat. current temperature exponent");
543 
545  .setUnit(U_NONE)
546  .setCategory(CAT_CAP)
547  .setDescription("Drain-source diode forward bias depletion capacitance");
548 
551  .setUnit(U_VOLT)
552  .setCategory(CAT_VOLT)
553  .setDescription("Drain-source diode reverse breakdown voltage");
554 
556  .setUnit(U_AMP)
557  .setCategory(CAT_CURRENT)
558  .setDescription("Drain-source diode current at breakdown voltage");
559 
560  p.addPar ("D1IKF",0.0,&VDMOS::Model::D1DIOikf)
561  .setUnit(U_AMP)
562  .setCategory(CAT_CURRENT)
563  .setDescription("Drain-source diode high injection knee currrent");
564 
565  p.addPar ("D1ISR",0.0,&VDMOS::Model::D1DIOisr)
566  .setUnit(U_AMP)
567  .setCategory(CAT_CURRENT)
568  .setDescription("Drain-source diode recombination saturation current");
569 
570  p.addPar ("D1NR",2.0,&VDMOS::Model::D1DIOnr)
571  .setUnit(U_NONE)
572  .setCategory(CAT_PROCESS)
573  .setDescription("Drain-source diode recombination emission coefficient");
574 
575  p.addPar ("D1KF",0.0,&VDMOS::Model::D1DIOfNcoef)
576  .setUnit(U_NONE)
577  .setCategory(CAT_FLICKER)
578  .setDescription("Drain-source diode flicker noise coefficient");
579 
580  p.addPar ("D1AF",1.0,&VDMOS::Model::D1DIOfNexp)
581  .setUnit(U_NONE)
582  .setCategory(CAT_FLICKER)
583  .setDescription("Drain-source diode flicker noise exponent");
584 
585  p.addPar ("D1TNOM",300.15,&VDMOS::Model::D1DIOnomTemp)
586  .setUnit(U_DEGC)
587  .setCategory(CAT_TEMP)
588  .setDescription("Drain-source diode nominal temperature");
589 
590  // Set up non-double precision variables:
591  p.addPar ("TPG",1,&VDMOS::Model::gateType)
592  .setUnit(U_NONE)
593  .setCategory(CAT_MATERIAL)
594  .setDescription("Gate material type (-1 = same as substrate, 0 = aluminum,1 = opposite of substrate)");
595 
596  p.addPar ("CV",1,&VDMOS::Model::cv)
597  .setUnit(U_NONE)
598  .setCategory(CAT_CONTROL)
599  .setDescription("Charge model storage selector");
600 
601  p.addPar ("CVE",1,&VDMOS::Model::cve)
602  .setUnit(U_NONE)
603  .setCategory(CAT_CONTROL)
604  .setDescription("Meyer-like capacitor model selector");
605 
606  p.addPar ("FPE",1,&VDMOS::Model::fpe)
607  .setUnit(U_NONE)
608  .setCategory(CAT_CONTROL)
609  .setDescription("Charge partitioning scheme selector");
610 
611  p.addPar ("ISUBMOD",0,&VDMOS::Model::isubmod)
612  .setUnit(U_NONE)
614  .setDescription("");
615 
617 }
618 
619 std::vector< std::vector<int> > Instance::jacStamp_DC_SC_GC;
620 std::vector< std::vector<int> > Instance::jacStamp_DC_GC;
621 std::vector< std::vector<int> > Instance::jacStamp_SC_GC;
622 std::vector< std::vector<int> > Instance::jacStamp_DC_SC;
623 std::vector< std::vector<int> > Instance::jacStamp_GC;
624 std::vector< std::vector<int> > Instance::jacStamp_SC;
625 std::vector< std::vector<int> > Instance::jacStamp_DC;
626 std::vector< std::vector<int> > Instance::jacStamp;
627 
628 std::vector<int> Instance::jacMap_DC_SC_GC;
629 std::vector<int> Instance::jacMap_DC_GC;
630 std::vector<int> Instance::jacMap_SC_GC;
631 std::vector<int> Instance::jacMap_DC_SC;
632 std::vector<int> Instance::jacMap_GC;
633 std::vector<int> Instance::jacMap_SC;
634 std::vector<int> Instance::jacMap_DC;
635 std::vector<int> Instance::jacMap;
636 
637 std::vector< std::vector<int> > Instance::jacMap2_DC_SC_GC;
638 std::vector< std::vector<int> > Instance::jacMap2_DC_GC;
639 std::vector< std::vector<int> > Instance::jacMap2_SC_GC;
640 std::vector< std::vector<int> > Instance::jacMap2_DC_SC;
641 std::vector< std::vector<int> > Instance::jacMap2_GC;
642 std::vector< std::vector<int> > Instance::jacMap2_SC;
643 std::vector< std::vector<int> > Instance::jacMap2_DC;
644 std::vector< std::vector<int> > Instance::jacMap2;
645 
646 // duplicating, but with RD1RS nonzero
647 std::vector< std::vector<int> > Instance::jacStamp_D1C_DC_SC_GC;
648 std::vector< std::vector<int> > Instance::jacStamp_D1C_DC_GC;
649 std::vector< std::vector<int> > Instance::jacStamp_D1C_SC_GC;
650 std::vector< std::vector<int> > Instance::jacStamp_D1C_DC_SC;
651 std::vector< std::vector<int> > Instance::jacStamp_D1C_GC;
652 std::vector< std::vector<int> > Instance::jacStamp_D1C_SC;
653 std::vector< std::vector<int> > Instance::jacStamp_D1C_DC;
654 std::vector< std::vector<int> > Instance::jacStamp_D1C;
655 
656 std::vector<int> Instance::jacMap_D1C_DC_SC_GC;
657 std::vector<int> Instance::jacMap_D1C_DC_GC;
658 std::vector<int> Instance::jacMap_D1C_SC_GC;
659 std::vector<int> Instance::jacMap_D1C_DC_SC;
660 std::vector<int> Instance::jacMap_D1C_GC;
661 std::vector<int> Instance::jacMap_D1C_SC;
662 std::vector<int> Instance::jacMap_D1C_DC;
663 std::vector<int> Instance::jacMap_D1C;
664 
665 std::vector< std::vector<int> > Instance::jacMap2_D1C_DC_SC_GC;
666 std::vector< std::vector<int> > Instance::jacMap2_D1C_DC_GC;
667 std::vector< std::vector<int> > Instance::jacMap2_D1C_SC_GC;
668 std::vector< std::vector<int> > Instance::jacMap2_D1C_DC_SC;
669 std::vector< std::vector<int> > Instance::jacMap2_D1C_GC;
670 std::vector< std::vector<int> > Instance::jacMap2_D1C_SC;
671 std::vector< std::vector<int> > Instance::jacMap2_D1C_DC;
672 std::vector< std::vector<int> > Instance::jacMap2_D1C;
673 
674 //--------------------- Class Model ----------------------------
675 
676 //-----------------------------------------------------------------------------
677 // Function : Model::Model
678 // Purpose :
679 // Special Notes :
680 // Scope : public
681 // Creator : Tom Russo, Component Information and Models
682 // Creation Date : 2/26/01
683 //-----------------------------------------------------------------------------
685  const Configuration & configuration,
686  const ModelBlock & MB,
687  const FactoryBlock & factory_block)
688  : DeviceModel(MB, configuration.getModelParameters(), factory_block),
689  dtype(CONSTNMOS),
690  gateType(0),
691 
692  l0(2e-6),
693  w0(2e-5),
694  tnom(getDeviceOptions().tnom),
695  latDiff(0.0),
696  jctSatCurDensity(0.0),
697  jctSatCur(0.0),
698  drainResistance(0.0),
699  gateResistance(0.0),
700  sourceResistance(0.0),
701  sheetResistance(0.0),
702  gateSourceOverlapCapFactor(0.0),
703  gateDrainOverlapCapFactor(0.0),
704  gateBulkOverlapCapFactor(0.0),
705  oxideCapFactor(0.0),
706  vt0(0.0),
707  capBD(0.0),
708  capBS(0.0),
709  timeScale(0.0),
710  bulkCapFactor(0.0),
711  sideWallCapFactor(0.0),
712  bulkJctPotential(0.0),
713  bulkJctBotGradingCoeff(0.0),
714  bulkJctSideGradingCoeff(0.0),
715  fwdCapDepCoeff(0.0),
716  phi(0.0),
717  gamma(0.0),
718  lambda(0.0),
719  substrateDoping(0.0),
720  surfaceStateDensity(0.0),
721  oxideThickness(0.0),
722  surfaceMobility(0.0),
723  surfaceMobility0(0.0),
724 
725  capBDGiven(0),
726  capBSGiven(0),
727  bulkCapFactorGiven(0),
728  sideWallCapFactorGiven(0),
729  vpGiven(0),
730 
731  maxDriftVel(0.0),
732  junctionDepth(0.0),
733  rdi(0.0),
734  rsi(0.0),
735 
736  delta(5.0),
737  eta(1.32),
738  m(4.0),
739  mc(3.0),
740  sigma0(0.048),
741  vsigmat(1.7),
742  vsigma(0.2),
743  theta(0.0),
744  gammas0(0.5),
745  gammal0(0.0),
746  lgammas(0.0),
747  wgammas(0.0),
748  lgammal_(0.0),
749  wgammal(0.0),
750  kacc(0.0),
751  gb(0.0),
752  knit(0.0),
753  nitd(0.0),
754  nits(0.0),
755  mm(0.5),
756  k(0.0),
757  deltaSqr(0.0),
758  kvt(0.0),
759  mdtemp(0.0),
760  kvs(0.0),
761  tvs(0.0),
762  mth(0.0),
763  artd(0.0),
764  brtd(0.035),
765  crtd(0.1472),
766  drtd(0.0052),
767  nrtd(0.115),
768  n2(1.0),
769 
770  xqc(0.6),
771  mcv(10.0),
772  vfb(0.0),
773  fpe(1),
774  alpha(1.05),
775 
776  cv(1),
777  cve(1),
778  ls(35e-9),
779  rsub(0.0),
780  vp(0.0),
781  ai(2e9),
782  bi(8e8),
783  delmax(0.9),
784  md(2.0),
785  isubmod(0),
786 
787  kaccd(0.0),
788  kaccs(0.0),
789  invm(0.0),
790  invmc(0.0),
791  invmd(0.0),
792 
793  fact1(0.0),
794  vtnom(0.0),
795  egfet1(0.0),
796  pbfact1(0.0),
797 
798  driftParamA(0.08),
799  driftParamB(0.013),
800  rdsshunt(0.0),
801 
802  D1DIOsatCur(0.0),
803  D1DIOresist(0.0),
804  D1DIOconductance(0.0),
805  D1DIOemissionCoeff(0.0),
806  D1DIOtransitTime(0.0),
807  D1DIOjunctionCap(0.0),
808  D1DIOjunctionPot(0.0),
809  D1DIOgradingCoeff(0.0),
810  D1DIOactivationEnergy(0.0),
811  D1DIOsaturationCurrentExp(0.0),
812  D1DIOdepletionCapCoeff(0.0),
813  D1DIObreakdownVoltage(0.0),
814  D1DIObreakdownCurrent(0.0),
815  D1DIOf2(0.0),
816  D1DIOf3(0.0),
817  D1DIOnomTemp(0.0),
818  D1DIOfNcoef(0.0),
819  D1DIOfNexp(0.0),
820  D1DIOikf(0.0),
821  D1DIOisr(0.0),
822  D1DIOnr(0.0),
823  D1DIObreakdownVoltageGiven(0)
824 
825 {
826 
827  if (getType() != "")
828  {
829  if (getType() == "NMOS") {
830  dtype = CONSTNMOS;
831  }
832  else if (getType() == "PMOS") {
833  dtype = CONSTPMOS;
834  }
835  else
836  {
837  UserError0(*this) << "Could not recognize the type for model " << getName();
838  }
839  }
840 
841 
842  // Set params to constant default values:
843  setDefaultParams ();
844 
845  // Set params according to .model line and constant defaults from metadata:
846  setModParams (MB.params);
847 
848  // Set any non-constant parameter defaults:
849  if (!given("TNOM"))
851  if (!given("L0"))
853  if (!given("W0"))
855  if (!given("ALPHA"))
856  alpha = (dtype == CONSTNMOS) ? 1.05 : 1.3;
857  if (capBD != 0)
858  capBDGiven = true;
859  if (capBS != 0)
860  capBSGiven = true;
861 
862  // Calculate any parameters specified as expressions:
863 
865 
866  // calculate dependent (ie computed) params and check for errors:
867  if (given("U0"))
868  {
869  if (given("UO"))
870  UserError0(*this) << "Both uo and u0 have been specified and, which is not allowed";
871  else
872  UserWarning0(*this) << "Surface mobility has been specified as u0 instead of uo, uo is the preferred syntax";
873 
875  }
876 
877  if(eta == 0)
878  {
879  UserError0(*this) << "ETA cannot be zero for level 18";
880  }
881 
882  if(driftParamA == 0 && driftParamB == 0)
883  {
884  UserError0(*this) << "Both driftParamA and driftParamB cannot be zero";
885  }
886 
887  if(cv != 1 && cv != 2)
888  {
889  UserError0(*this) << "Model error: use cv=1 (Meyer's model) or cv=2 (Meyer-like Model)";
890  }
891 
892  if(fpe != 1 && fpe != 2 && fpe != 3)
893  {
894  UserError0(*this) << "Model error: use fpe=1 (and xqc), fpe=2 (and xqc,mcv) or fpe=3";
895  }
896 
897  if(mcv <= 1)
898  {
899  UserError0(*this) << "Model error: charge conservation requires mcv > 1";
900  }
901 
902  if(xqc > 1 || xqc < 0.5)
903  {
904  UserError0(*this) << "Model error: charge conservation requires 0.5 <= xqc <= 1.0";
905  }
906 
907  // New version of VDMOSModel. We have to set lambda to zero
908  if(isubmod) lambda = 0.0;
909 
910  // limit diode #1 parameters
911  // limit grading coeff to max of 0.9
912  if(D1DIOgradingCoeff > 0.9)
913  {
914  UserWarning0(*this) << "grading coefficient too large, limited to 0.9";
915  D1DIOgradingCoeff = 0.9;
916  }
917 
918  // limit activation energy to min of 0.1
919  if(D1DIOactivationEnergy < 0.1)
920  {
921  UserWarning0(*this) << "activation energy too small, limited to 0.1";
922  D1DIOactivationEnergy = 0.1;
923  }
924 
925  // limit depletion cap coeff to max of 0.95
926  if(D1DIOdepletionCapCoeff > 0.95)
927  {
928  UserWarning0(*this) << "coefficient Fc too large, limited to 0.95";
929  D1DIOdepletionCapCoeff = 0.95;
930  }
931 
932  processParams ();
933 
934  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS))
935  {
936  Xyce::dout() << " " << std::endl;
937  Xyce::dout() << "l0 " << l0 << std::endl;
938  Xyce::dout() << "w0 " << w0 << std::endl;
939  Xyce::dout() << "vt0 = " << vt0 << std::endl;
940  Xyce::dout() << "gamma = " << gamma << std::endl;
941  Xyce::dout() << "lambda = " << lambda << std::endl;
942  Xyce::dout() << "phi = " << phi << std::endl;
943  Xyce::dout() << "drainResistance = " << drainResistance << std::endl;
944  Xyce::dout() << "gateResistance = " << gateResistance << std::endl;
945  Xyce::dout() << "sourceResistance = " << sourceResistance << std::endl;
946  Xyce::dout() << "capBD = " << capBD << std::endl;
947  Xyce::dout() << "capBS = " << capBS << std::endl;
948  Xyce::dout() << "timeScale = " << timeScale << std::endl;
949  Xyce::dout() << "jctSatCur = " << jctSatCur << std::endl;
950  Xyce::dout() << "bulkJctPotential = " << bulkJctPotential << std::endl;
951  Xyce::dout() << "gateSourceOverlapCapFactor = " << gateSourceOverlapCapFactor << std::endl;
952  Xyce::dout() << "gateDrainOverlapCapFactor = " << gateDrainOverlapCapFactor << std::endl;
953  Xyce::dout() << "gateBulkOverlapCapFactor = " << gateBulkOverlapCapFactor << std::endl;
954  Xyce::dout() << "sheetResistance = " << sheetResistance << std::endl;
955  Xyce::dout() << "bulkCapFactor = " << bulkCapFactor << std::endl;
956  Xyce::dout() << "bulkJctBotGradingCoeff = " << bulkJctBotGradingCoeff << std::endl;
957  Xyce::dout() << "sideWallCapFactor = " << sideWallCapFactor << std::endl;
958  Xyce::dout() << "bulkJctSideGradingCoeff = " << bulkJctSideGradingCoeff << std::endl;
959  Xyce::dout() << "jctSatCurDensity = " << jctSatCurDensity << std::endl;
960  Xyce::dout() << "oxideThickness = " << oxideThickness << std::endl;
961  Xyce::dout() << "oxideCapFactor = " << oxideCapFactor << std::endl;
962  Xyce::dout() << "latDiff = " << latDiff << std::endl;
963  Xyce::dout() << "surfaceMobility = " << surfaceMobility << std::endl;
964  Xyce::dout() << "fwdCapDepCoeff = " << fwdCapDepCoeff << std::endl;
965  Xyce::dout() << "substrateDoping = " << substrateDoping << std::endl;
966  Xyce::dout() << "gateType = " << gateType << std::endl;
967  Xyce::dout() << "surfaceStateDensity = " << surfaceStateDensity << std::endl;
968  Xyce::dout() << "tnom = " << tnom << std::endl;
969  Xyce::dout() << "driftParamA = " << driftParamA << std::endl;
970  Xyce::dout() << "driftParamB = " << driftParamB << std::endl;
971  Xyce::dout() << "rdsshunt = " << rdsshunt << std::endl;
972  }
973 
974 }
975 
976 //-----------------------------------------------------------------------------
977 // Function : Model::~Model
978 // Purpose : destructor
979 // Special Notes :
980 // Scope : public
981 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
982 // Creation Date : 3/16/00
983 //-----------------------------------------------------------------------------
985 {
986  std::vector<Instance*>::iterator iter;
987  std::vector<Instance*>::iterator first = instanceContainer.begin();
988  std::vector<Instance*>::iterator last = instanceContainer.end();
989 
990  for (iter=first; iter!=last; ++iter)
991  {
992  delete (*iter);
993  }
994 
995 }
996 
997 
998 //-----------------------------------------------------------------------------
999 // Function : Model::printOutInstances
1000 // Purpose : debugging tool.
1001 // Special Notes :
1002 // Scope : public
1003 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1004 // Creation Date : 4/03/00
1005 //-----------------------------------------------------------------------------
1006 std::ostream &Model::printOutInstances(std::ostream &os) const
1007 {
1008  std::vector<Instance*>::const_iterator iter;
1009  std::vector<Instance*>::const_iterator first = instanceContainer.begin();
1010  std::vector<Instance*>::const_iterator last = instanceContainer.end();
1011 
1012  int i;
1013  os << std::endl;
1014  os << " name model name Parameters" << std::endl;
1015  for (i=0, iter=first; iter!=last; ++iter, ++i)
1016  {
1017  os << " " << i << ": " << (*iter)->getName() << "\t";
1018  os << getName();
1019  os << std::endl;
1020  }
1021 
1022  os << std::endl;
1023 
1024  return os;
1025 }
1026 
1027 //-----------------------------------------------------------------------------
1028 // Function : Model::forEachInstance
1029 // Purpose :
1030 // Special Notes :
1031 // Scope : public
1032 // Creator : David Baur
1033 // Creation Date : 2/4/2014
1034 //-----------------------------------------------------------------------------
1035 /// Apply a device instance "op" to all instances associated with this
1036 /// model
1037 ///
1038 /// @param[in] op Operator to apply to all instances.
1039 ///
1040 ///
1041 void Model::forEachInstance(DeviceInstanceOp &op) const /* override */
1042 {
1043  for (std::vector<Instance *>::const_iterator it = instanceContainer.begin(); it != instanceContainer.end(); ++it)
1044  op(*it);
1045 }
1046 
1047 
1048 //-----------------------------------------------------------------------------
1049 // Function : Model::processParams
1050 // Purpose :
1051 // Special Notes :
1052 // Scope : public
1053 // Creator : pmc
1054 // Creation Date : 2/17/2004
1055 //-----------------------------------------------------------------------------
1057 {
1058  double wkfngs(0.0);
1059  double wkfng(0.0);
1060  double fermig(0.0);
1061  double fermis(0.0);
1062  double kt1(0.0);
1063  double arg1(0.0);
1064 
1067  kt1 = CONSTboltz * tnom;
1068  egfet1 = 1.16-(7.02e-4*tnom*tnom)/(tnom+1108);
1069  arg1 = -egfet1/(kt1+kt1)+1.1150877/(CONSTboltz*(CONSTREFTEMP+CONSTREFTEMP));
1070  pbfact1= -2*vtnom *(1.5*log(fact1)+CONSTQ*arg1);
1071 
1072 // assuming silicon - make definition for epsilon of silicon
1073 //#define EPSSIL (11.7 * 8.854214871e-12) = CONSTEPSSI (the Xyce value)
1074 
1075  oxideCapFactor = 3.9 * 8.854214871e-12/oxideThickness;
1076  if(given("NSUB"))
1077  {
1078  if(substrateDoping*1e6 >1.45e16)
1079  {
1080  if(!given("PHI"))
1081  {
1082  phi = 2*vtnom* log(substrateDoping*1e6/1.45e16);
1083  phi = std::max(.1,phi);
1084  }
1085  fermis = dtype*0.5*phi;
1086  wkfng = 3.2;
1087  if(gateType != 0)
1088  {
1089  fermig = dtype *gateType*0.5*egfet1;
1090  wkfng = 3.25 + 0.5 * egfet1 - fermig;
1091  }
1092  wkfngs = wkfng - (3.25 + 0.5 * egfet1 + fermis);
1093  gamma = sqrt(2 * 11.70 * 8.854214871e-12 * CONSTQ
1095  if(!given("GAMMAS0")) gammas0 = gamma;
1096 
1097  if(!given("VTO"))
1098  {
1099  if(!given("NSS")) surfaceStateDensity=0;
1100  if(!given("VFB"))
1101  {
1103  }
1104  vt0 = vfb + dtype * (gamma * sqrt(phi)+ phi);
1105  }
1106  else
1107  {
1108  if(!given("VFB")) vfb = vt0 - dtype * (gamma * sqrt(phi)+ phi);
1109  }
1110  }
1111  else
1112  {
1113  UserError0(*this) << "SubstrateDoping is 0";
1114 
1115  substrateDoping = 0;
1116  }
1117  }
1118 
1119  if(!given("CJ"))
1120  {
1122  /(2*bulkJctPotential));
1123  }
1124  if(!given("VFB"))
1125  {
1126  if(fabs(vfb) < 1e-12)
1127  {
1128  vfb = 1e-12;
1129  }
1130  else
1131  {
1132  //vfb = vfb;
1133  }
1134  }
1135 
1136  kaccd = kacc/(1+pow(k*nitd,mm));
1137  kaccs = kacc/(1+pow(k*nits,mm));
1138  deltaSqr = delta*delta;
1139 
1140  if(D1DIOresist == 0)
1141  D1DIOconductance = 0;
1142  else
1144 
1145  double D1xfc = log(1-D1DIOdepletionCapCoeff);
1146  D1DIOf2 = exp((1 + D1DIOgradingCoeff)*D1xfc);
1148 
1149  return true;
1150 }
1151 
1152 //----------------------------------------------------------------------------
1153 // Function : Model::processInstanceParams
1154 // Purpose :
1155 // Special Notes :
1156 // Scope : public
1157 // Creator : Dave Shirely, PSSI
1158 // Creation Date : 03/23/06
1159 //----------------------------------------------------------------------------
1161 {
1162 
1163  std::vector<Instance*>::iterator iter;
1164  std::vector<Instance*>::iterator first = instanceContainer.begin();
1165  std::vector<Instance*>::iterator last = instanceContainer.end();
1166 
1167  for (iter=first; iter!=last; ++iter)
1168  {
1169  (*iter)->processParams();
1170  }
1171 
1172  return true;
1173 }
1174 
1175 //------------------ Class Instance ----------------------------
1176 
1177 
1178 //-----------------------------------------------------------------------------
1179 // Function : Instance::Instance
1180 // Purpose : instance block constructor
1181 // Special Notes :
1182 // Scope : public
1183 // Creator : Tom Russo, Component Information and Models
1184 // Creation Date : 3/21/01
1185 //-----------------------------------------------------------------------------
1187  const Configuration & configuration,
1188  const InstanceBlock & IB,
1189  Model & Miter,
1190  const FactoryBlock & factory_block)
1191  : DeviceInstance(IB, configuration.getInstanceParameters(), factory_block),
1192  model_(Miter),
1193  dNode(0),
1194  gNode(0),
1195  sNode(0),
1196  bNode(0),
1197  dNodePrime(0),
1198  sNodePrime(0),
1199  dDriftNode(0),
1200  l(getDeviceOptions().defl),
1201  w(getDeviceOptions().defw),
1202  drainArea(getDeviceOptions().defad),
1203  sourceArea(getDeviceOptions().defas),
1204  drainSquares(1.0),
1205  sourceSquares(1.0),
1206  drainPerimeter(0.0),
1207  sourcePerimeter(0.0),
1208  sourceCond(0.0),
1209  gateCond(0.0),
1210  drainCond(0.0),
1211  draindriftCond(0.0),
1212  vt(0.0),
1213  gammas(0.0),
1214  gammal(0.0),
1215  gchi0(0.0),
1216  vtoo(0.0),
1217  vthLimit(0.0),
1218  numberParallel(1),
1219 
1220  temp(getDeviceOptions().temp.getImmutableValue<double>()),
1221  tSurfMob(0.0),
1222  tPhi(0.0),
1223  tVto(0.0),
1224  tSatCur(0.0),
1225  tSatCurDens(0.0),
1226  tCbd(0.0),
1227  tCbs(0.0),
1228  tCj(0.0),
1229  tCjsw(0.0),
1230  tBulkPot(0.0),
1231  tDepCap(0.0),
1232  tVbi(0.0),
1233 
1234  von(0.0),
1235  vdsat(0.0),
1236  vddsat(0.0),
1237  sourceVcrit(0.0),
1238  drainVcrit(0.0),
1239  draindriftVcrit(0.0),
1240  cd(0.0),
1241  cdd(0.0),
1242  gmbs(0.0),
1243  gm(0.0),
1244  gddd(0.0),
1245  dIdd_dVd(0.0),
1246  gds(0.0),
1247  gdds(0.0),
1248  gdsshunt(0.0),
1249  gbd(0.0),
1250  gbs(0.0),
1251  cbd(0.0),
1252  Cbd(0.0),
1253  Cbdsw(0.0),
1254  cbs(0.0),
1255  Cbs(0.0),
1256  Cbssw(0.0),
1257  f2d(0.0),
1258  f3d(0.0),
1259  f4d(0.0),
1260  f2s(0.0),
1261  f3s(0.0),
1262  f4s(0.0),
1263  mode(1),
1264  mode_low(0.0),
1265  mode_high(0.0),
1266  off(0),
1267  dNodePrimeSet(0),
1268  sNodePrimeSet(0),
1269  limitedFlag(false),
1270  //calculated quantities
1271  EffectiveLength(0),
1272  DrainSatCur(0),
1273  SourceSatCur(0),
1274  GateSourceOverlapCap(0),
1275  GateDrainOverlapCap(0),
1276  GateBulkOverlapCap(0),
1277  OxideCap(0),
1278 
1279  Isource(0.0),
1280  Igate(0.0),
1281  Idrain(0.0),
1282  Idraindrift(0.0),
1283  Irdsshunt(0.0),
1284  Ird1rs(0.0),
1285  mm1(0.0),
1286  dmm1vgs(0.0),
1287  dmm1vds(0.0),
1288  dmm1vbs(0.0),
1289  cdrain(0.0),
1290  cdraindrift(0.0),
1291 
1292  D1DIOoff1(0),
1293  D1DIOarea(1.0),
1294  D1DIOinitCond(0.0),
1295  D1DIOtemp(0.0),
1296  D1DIOtJctPot(0.0),
1297  D1DIOtJctCap(0.0),
1298  D1DIOtDepCap(0.0),
1299  D1DIOtSatCur(0.0),
1300  D1DIOtSatRCur(0.0),
1301  D1DIOtVcrit(0.0),
1302  D1DIOtF1(0.0),
1303  D1DIOtBrkdwnV(0.0),
1304  D1gspr(0.0),
1305  D1gd(0.0),
1306  D1cdeq(0.0),
1307  D1vt(0.0),
1308  D1vte(0.0),
1309  D1capd(0.0),
1310 
1311  D1DIOcapCharge(0.0),
1312  D1DIOcapCurrent(0.0),
1313 
1314  Vd (0.0),
1315  Vs (0.0),
1316  Vg (0.0),
1317  Vb (0.0),
1318  Vdp (0.0),
1319  Vgp (0.0),
1320  Vsp (0.0),
1321  Vdd (0.0),
1322  Vddp (0.0),
1323  Vddd (0.0),
1324  Vdddp (0.0),
1325  Vssp (0.0),
1326  Vbsp (0.0),
1327  Vbdp (0.0),
1328  Vggp (0.0),
1329  Vgpsp (0.0),
1330  Vgpdp (0.0),
1331  Vgpb (0.0),
1332  Vdpsp (0.0),
1333  Vbdd (0.0),
1334  D1vd (0.0),
1335  Capgs (0.0),
1336  Capgdd(0.0),
1337  Capgb (0.0),
1338  Gm (0.0),
1339  Gmbs (0.0),
1340  revsum(0.0),
1341  nrmsum(0.0),
1342 
1343  // solution,rhs vector:
1344  // local indices
1345  li_Drain (-1),
1346  li_DrainPrime (-1),
1347  li_Source (-1),
1348  li_SourcePrime (-1),
1349  li_Gate (-1),
1350  li_GatePrime (-1),
1351  li_Bulk (-1),
1352  li_DrainDrift (-1),
1353  li_D1Prime (-1),
1354 
1355  // jacobian matrix:
1356  // matrix and vector pointers:
1357  // drain row
1358  ADrainEquDrainNodeOffset (-1),
1359  ADrainEquSourceNodeOffset (-1),
1360  ADrainEquDrainDriftNodeOffset (-1),
1361  ADrainEquD1PrimeNodeOffset (-1),
1362  // gate row
1363  AGateEquGateNodeOffset (-1),
1364  AGateEquGatePrimeNodeOffset (-1),
1365  // source row
1366  ASourceEquDrainNodeOffset (-1),
1367  ASourceEquSourceNodeOffset (-1),
1368  ASourceEquSourcePrimeNodeOffset (-1),
1369  // bulk row
1370  ABulkEquBulkNodeOffset (-1),
1371  ABulkEquDrainPrimeNodeOffset (-1),
1372  ABulkEquGatePrimeNodeOffset (-1),
1373  ABulkEquSourcePrimeNodeOffset (-1),
1374  // drain' row
1375  ADrainPrimeEquBulkNodeOffset (-1),
1376  ADrainPrimeEquDrainPrimeNodeOffset (-1),
1377  ADrainPrimeEquGatePrimeNodeOffset (-1),
1378  ADrainPrimeEquSourcePrimeNodeOffset (-1),
1379  ADrainPrimeEquDrainDriftNodeOffset (-1),
1380  // gate' row
1381  AGatePrimeEquGateNodeOffset (-1),
1382  AGatePrimeEquBulkNodeOffset (-1),
1383  AGatePrimeEquDrainPrimeNodeOffset (-1),
1384  AGatePrimeEquGatePrimeNodeOffset (-1),
1385  AGatePrimeEquSourcePrimeNodeOffset (-1),
1386  // source' row
1387  ASourcePrimeEquSourceNodeOffset (-1),
1388  ASourcePrimeEquBulkNodeOffset (-1),
1389  ASourcePrimeEquDrainPrimeNodeOffset (-1),
1390  ASourcePrimeEquGatePrimeNodeOffset (-1),
1391  ASourcePrimeEquSourcePrimeNodeOffset (-1),
1392  // drain drift row
1393  ADrainDriftEquDrainNodeOffset (-1),
1394  ADrainDriftEquDrainPrimeNodeOffset (-1),
1395  ADrainDriftEquDrainDriftNodeOffset (-1),
1396  // D1Prime row
1397  AD1PrimeEquDrainNodeOffset (-1),
1398  AD1PrimeEquSourceNodeOffset (-1),
1399  AD1PrimeEquD1PrimeNodeOffset (-1),
1400 
1402  // F-matrix pointers:
1403  f_DrainEquDrainNodePtr(0),
1404  f_DrainEquSourceNodePtr(0),
1405  f_DrainEquDrainDriftNodePtr(0),
1406  f_DrainEquD1PrimeNodePtr(0),
1407 
1408  f_GateEquGateNodePtr(0),
1409  f_GateEquGatePrimeNodePtr(0),
1410 
1411  f_SourceEquDrainNodePtr(0),
1412  f_SourceEquSourceNodePtr(0),
1413  f_SourceEquSourcePrimeNodePtr(0),
1414  f_SourceEquD1PrimeNodePtr(0),
1415 
1416  f_BulkEquBulkNodePtr(0),
1417  f_BulkEquDrainPrimeNodePtr(0),
1418  f_BulkEquGatePrimeNodePtr(0),
1419  f_BulkEquSourcePrimeNodePtr(0),
1420 
1421  f_DrainPrimeEquBulkNodePtr(0),
1422  f_DrainPrimeEquDrainPrimeNodePtr(0),
1423  f_DrainPrimeEquGatePrimeNodePtr(0),
1424  f_DrainPrimeEquSourcePrimeNodePtr(0),
1425  f_DrainPrimeEquDrainDriftNodePtr(0),
1426 
1427  f_GatePrimeEquGateNodePtr(0),
1428  f_GatePrimeEquBulkNodePtr(0),
1429  f_GatePrimeEquDrainPrimeNodePtr(0),
1430  f_GatePrimeEquGatePrimeNodePtr(0),
1431  f_GatePrimeEquSourcePrimeNodePtr(0),
1432 
1433  f_SourcePrimeEquSourceNodePtr(0),
1434  f_SourcePrimeEquBulkNodePtr(0),
1435  f_SourcePrimeEquDrainPrimeNodePtr(0),
1436  f_SourcePrimeEquGatePrimeNodePtr(0),
1437  f_SourcePrimeEquSourcePrimeNodePtr(0),
1438 
1439  f_DrainDriftEquDrainNodePtr(0),
1440  f_DrainDriftEquDrainPrimeNodePtr(0),
1441  f_DrainDriftEquDrainDriftNodePtr(0),
1442 
1443  f_D1PrimeEquDrainNodePtr(0),
1444  f_D1PrimeEquSourceNodePtr(0),
1445  f_D1PrimeEquD1PrimeNodePtr(0),
1446 
1447  // Q-matrix pointers:
1448  q_DrainEquDrainNodePtr(0),
1449  q_DrainEquSourceNodePtr(0),
1450  q_DrainEquDrainDriftNodePtr(0),
1451  q_DrainEquD1PrimeNodePtr(0),
1452 
1453  q_GateEquGateNodePtr(0),
1454  q_GateEquGatePrimeNodePtr(0),
1455 
1456  q_SourceEquDrainNodePtr(0),
1457  q_SourceEquSourceNodePtr(0),
1458  q_SourceEquSourcePrimeNodePtr(0),
1459  q_SourceEquD1PrimeNodePtr(0),
1460 
1461  q_BulkEquBulkNodePtr(0),
1462  q_BulkEquDrainPrimeNodePtr(0),
1463  q_BulkEquGatePrimeNodePtr(0),
1464  q_BulkEquSourcePrimeNodePtr(0),
1465 
1466  q_DrainPrimeEquBulkNodePtr(0),
1467  q_DrainPrimeEquDrainPrimeNodePtr(0),
1468  q_DrainPrimeEquGatePrimeNodePtr(0),
1469  q_DrainPrimeEquSourcePrimeNodePtr(0),
1470  q_DrainPrimeEquDrainDriftNodePtr(0),
1471 
1472  q_GatePrimeEquGateNodePtr(0),
1473  q_GatePrimeEquBulkNodePtr(0),
1474  q_GatePrimeEquDrainPrimeNodePtr(0),
1475  q_GatePrimeEquGatePrimeNodePtr(0),
1476  q_GatePrimeEquSourcePrimeNodePtr(0),
1477 
1478  q_SourcePrimeEquSourceNodePtr(0),
1479  q_SourcePrimeEquBulkNodePtr(0),
1480  q_SourcePrimeEquDrainPrimeNodePtr(0),
1481  q_SourcePrimeEquGatePrimeNodePtr(0),
1482  q_SourcePrimeEquSourcePrimeNodePtr(0),
1483 
1484  q_DrainDriftEquDrainNodePtr(0),
1485  q_DrainDriftEquDrainPrimeNodePtr(0),
1486  q_DrainDriftEquDrainDriftNodePtr(0),
1487 
1488  q_D1PrimeEquDrainNodePtr(0),
1489  q_D1PrimeEquSourceNodePtr(0),
1490  q_D1PrimeEquD1PrimeNodePtr(0),
1491 #endif
1492 
1493  vbdd(0.0),
1494  vbs(0.0),
1495  vgpdd(0.0),
1496  vgps(0.0),
1497  vdds(0.0),
1498 
1499  vbdd_old(0.0),
1500  vbs_old(0.0),
1501  vgpdd_old(0.0),
1502  vgps_old(0.0),
1503  vdds_old(0.0),
1504  D1vd_old(0.0),
1505 
1506  vbdd_orig(0.0),
1507  vbs_orig(0.0),
1508  vgpdd_orig(0.0),
1509  vgps_orig(0.0),
1510  vdds_orig(0.0),
1511  D1vd_orig(0.0),
1512 
1513  capgs(0.0),
1514  qgs(0.0),
1515  cqgs(0.0),
1516  capgdd(0.0),
1517  qgdd(0.0),
1518  cqgdd(0.0),
1519  capgb(0.0),
1520  qgb(0.0),
1521  cqgb(0.0),
1522  capbd(0.0),
1523  qbd(0.0),
1524  cqbd(0.0),
1525  capbs(0.0),
1526  qbs(0.0),
1527  cqbs(0.0),
1528 
1529  // local indices
1530  li_state_vbdd(-1),
1531  li_state_vbs(-1),
1532  li_state_vgps(-1),
1533  li_state_vdds(-1),
1534  li_state_D1vd(-1),
1535  li_state_qgs(-1),
1536  li_state_qgdd(-1),
1537  li_state_qgb(-1),
1538  li_state_capgs(-1),
1539  li_state_capgdd(-1),
1540  li_state_capgb(-1),
1541  li_state_qbd(-1),
1542  li_state_qbs(-1),
1543  //li_state_cqgs(-1),
1544  //li_state_cqgdd(-1),
1545  //li_state_cqgb(-1),
1546  //li_state_cqbd(-1),
1547  //li_state_cqbs(-1),
1548 
1549  li_state_D1DIOcapCharge(-1),
1550  //li_state_D1DIOcapCurrent(-1),
1551  li_state_von(-1),
1552  li_store_dev_id(-1),
1553  li_store_dev_ig(-1),
1554  li_store_dev_is(-1),
1555  li_store_dev_ib(-1)
1556 {
1557  //numStateVars = 21;
1558  numStateVars = 15;
1559  numExtVars = 4;
1560  numLeadCurrentStoreVars = 4; // drain, gate, source & base lead currents
1561 
1562  devConMap.resize(4);
1563  devConMap[0] = 1;
1564  devConMap[1] = 2;
1565  devConMap[2] = 1;
1566  devConMap[3] = 3;
1567 
1568  if( jacStamp.empty() )
1569  {
1570  // stamp for RS!=0, RD!=0, RG!=0, RD1RS!=0
1571  jacStamp_D1C_DC_SC_GC.resize(9);
1572  jacStamp_D1C_DC_SC_GC[0].resize(4); // Drain row
1573  jacStamp_D1C_DC_SC_GC[0][0]=0; // d-d
1574  jacStamp_D1C_DC_SC_GC[0][1]=2; // d-s
1575  jacStamp_D1C_DC_SC_GC[0][2]=7; // d-dd
1576  jacStamp_D1C_DC_SC_GC[0][3]=8; // d-d1'
1577  jacStamp_D1C_DC_SC_GC[1].resize(2); // Gate row
1578  jacStamp_D1C_DC_SC_GC[1][0]=1; // g-g
1579  jacStamp_D1C_DC_SC_GC[1][1]=5; // g-g'
1580  jacStamp_D1C_DC_SC_GC[2].resize(4); // Source row
1581  jacStamp_D1C_DC_SC_GC[2][0]=0; // s-d
1582  jacStamp_D1C_DC_SC_GC[2][1]=2; // s-s
1583  jacStamp_D1C_DC_SC_GC[2][2]=6; // s-s'
1584  jacStamp_D1C_DC_SC_GC[2][3]=8; // s-d1'
1585  jacStamp_D1C_DC_SC_GC[3].resize(4); // Bulk row
1586  jacStamp_D1C_DC_SC_GC[3][0]=3; // b-b
1587  jacStamp_D1C_DC_SC_GC[3][1]=4; // b-d'
1588  jacStamp_D1C_DC_SC_GC[3][2]=5; // b-g'
1589  jacStamp_D1C_DC_SC_GC[3][3]=6; // b-s'
1590  jacStamp_D1C_DC_SC_GC[4].resize(5); // Drain' row
1591  jacStamp_D1C_DC_SC_GC[4][0]=3; // d'-b
1592  jacStamp_D1C_DC_SC_GC[4][1]=4; // d'-d'
1593  jacStamp_D1C_DC_SC_GC[4][2]=5; // d'-g'
1594  jacStamp_D1C_DC_SC_GC[4][3]=6; // d'-s'
1595  jacStamp_D1C_DC_SC_GC[4][4]=7; // d'-dd
1596  jacStamp_D1C_DC_SC_GC[5].resize(5); // Gate' row
1597  jacStamp_D1C_DC_SC_GC[5][0]=1; // g'-g
1598  jacStamp_D1C_DC_SC_GC[5][1]=3; // g'-b
1599  jacStamp_D1C_DC_SC_GC[5][2]=4; // g'-d'
1600  jacStamp_D1C_DC_SC_GC[5][3]=5; // g'-g'
1601  jacStamp_D1C_DC_SC_GC[5][4]=6; // g'-s'
1602  jacStamp_D1C_DC_SC_GC[6].resize(5); // Source' row
1603  jacStamp_D1C_DC_SC_GC[6][0]=2; // s'-s
1604  jacStamp_D1C_DC_SC_GC[6][1]=3; // s'-b
1605  jacStamp_D1C_DC_SC_GC[6][2]=4; // s'-d'
1606  jacStamp_D1C_DC_SC_GC[6][3]=5; // s'-g'
1607  jacStamp_D1C_DC_SC_GC[6][4]=6; // s'-s'
1608  jacStamp_D1C_DC_SC_GC[7].resize(3); // DrainDrift row
1609  jacStamp_D1C_DC_SC_GC[7][0]=0; // dd-d
1610  jacStamp_D1C_DC_SC_GC[7][1]=4; // dd-d'
1611  jacStamp_D1C_DC_SC_GC[7][2]=7; // dd-dd
1612  jacStamp_D1C_DC_SC_GC[8].resize(3); // D1'pos row
1613  jacStamp_D1C_DC_SC_GC[8][0]=0; // d1'-d
1614  jacStamp_D1C_DC_SC_GC[8][1]=2; // d1'-s
1615  jacStamp_D1C_DC_SC_GC[8][2]=8; // d1'-d1'
1616  jacMap_D1C_DC_SC_GC.clear();
1617 
1618  // First, generate the stamp for when RD1RS is zero --- it will be used
1619  // for all the other calculations for that case.
1620  // When RD1RS is zero, d1' is just S
1624  8, 2, 9);
1625 
1626  // Now do the cases where RD1RS is nonzero, but one of the others is
1630  6, 2, 9); // s' becomes same as s
1634  7, 4, 9); // dd becomes same as d'
1638  5, 1, 9); // g' becomes same as g
1639 
1640  // s' is already same as s here, so dd has become 6 instead of 7, and
1641  // needs to be mapped onto d' (which is still 4)
1644 
1645  // g' is same as g here, making s' 5 instead of 6, map it onto s now:
1648  // g' is same as g here, making dd 6 instead of 7, map it onto d' now:
1651 
1652  // g'=g and s'=s, so dd is 5 instead of 7, map it onto d'
1654  jacStamp_D1C, jacMap_D1C, jacMap2_D1C, 5, 4, 9);
1655 
1656  // Now do the ones where RD1RS is zero, starting from the case where
1657  // everything *but* RD1RS is nonzero. These follow exactly the same
1658  // pattern as above, but starting from _DC_SC_GC instead of _D1C_DC_SC_GC
1659  // This only works out this way because d1' is the last row, and nothing
1660  // gets moved up when it goes away.
1661  // I won't repeat the comments for each block.
1662  // only one row goes away:
1669 
1670  // s'==s already
1672  jacStamp_GC, jacMap_GC, jacMap2_GC, 6, 4, 9);
1673 
1674  // g'=g already
1676  jacStamp_DC, jacMap_DC, jacMap2_DC, 5, 2, 9);
1678  jacStamp_SC, jacMap_SC, jacMap2_SC, 6, 4, 9);
1679 
1680  // s'=s, g'=g already
1682  jacStamp, jacMap, jacMap2, 5, 4, 9);
1683  }
1684 
1685 
1686  // Set params to constant default values:
1687  setDefaultParams ();
1688 
1689  // Set params according to instance line and constant defaults from metadata:
1690  setParams (IB.params);
1691 
1692  // Set any non-constant parameter defaults:
1693  if (!given("TEMP"))
1694  temp = getDeviceOptions().temp.getImmutableValue<double>();
1695  if (!given("L"))
1696  l = model_.l0;
1697  if (!given("W"))
1698  w = model_.w0;
1699 
1700  // Calculate any parameters specified as expressions:
1702 
1703  // calculate dependent (ie computed) params and check for errors:
1704 
1705  // call updateTemp which may change model parameters if
1706  // temperature effects are present
1707  processParams ();
1708 
1709  if(model_.drainResistance != 0)
1710  {
1712  }
1713  else if (model_.given("RSH"))
1714  {
1715  if(model_.sheetResistance != 0)
1716  {
1717  drainCond =
1719  }
1720  else
1721  {
1722  drainCond = 0;
1723  }
1724  }
1725  else
1726  {
1727  drainCond = 0;
1728  }
1729  if(model_.sourceResistance != 0)
1730  {
1732  }
1733  else if (model_.given("RSH"))
1734  {
1735  if(model_.sheetResistance != 0)
1736  {
1737  sourceCond =
1739  }
1740  else
1741  {
1742  sourceCond = 0;
1743  }
1744  }
1745  else
1746  {
1747  sourceCond = 0;
1748  }
1749 
1750  if(model_.given("RG"))
1751  {
1752  if(model_.gateResistance != 0)
1753  {
1755  }
1756  else
1757  {
1758  gateCond = 0;
1759  }
1760  }
1761  else
1762  {
1763  gateCond = 0;
1764  }
1765 
1766  w *= numberParallel;
1767 
1769  if(gddd != 0) gddd = 1.0/gddd;
1770  draindriftCond = gddd;
1771 
1772  if(l - 2 * model_.latDiff <=0)
1773  {
1774  UserError0(*this) << "Effective channel length less than zero.";
1775  }
1776 
1777  // values for diode #1
1778  D1DIOarea = 1;
1779 
1785 
1786  numIntVars = 1 + (((sourceCond == 0.0)?0:1)+((drainCond == 0.0) ? 0:1)
1787  + ((gateCond == 0.0) ? 0:1))
1788  + ((model_.D1DIOconductance == 0.0) ? 0:1);
1789 }
1790 
1791 //-----------------------------------------------------------------------------
1792 // Function : Instance::~Instance
1793 // Purpose : destructor
1794 // Special Notes :
1795 // Scope : public
1796 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1797 // Creation Date : 3/16/00
1798 //-----------------------------------------------------------------------------
1800 {
1801 }
1802 
1803 //-----------------------------------------------------------------------------
1804 // Function : Instance::registerLIDs
1805 // Purpose :
1806 // Special Notes :
1807 // Scope : public
1808 // Creator : Robert Hoekstra, Computational Sciences
1809 // Creation Date : 6/21/02
1810 //-----------------------------------------------------------------------------
1811 void Instance::registerLIDs( const std::vector<int> & intLIDVecRef,
1812  const std::vector<int> & extLIDVecRef )
1813 {
1814  AssertLIDs(intLIDVecRef.size() == numIntVars);
1815  AssertLIDs(extLIDVecRef.size() == numExtVars);
1816 
1817  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS))
1818  {
1819  Xyce::dout() << section_divider << std::endl;
1820  Xyce::dout() << " In Instance::register LIDs\n\n";
1821  Xyce::dout() << " name = " << getName() << std::endl;
1822  Xyce::dout() << " number of internal variables: " << numIntVars << std::endl;
1823  Xyce::dout() << " number of external variables: " << numExtVars << std::endl;
1824  }
1825 
1826  // copy over the global ID lists.
1827  intLIDVec = intLIDVecRef;
1828  extLIDVec = extLIDVecRef;
1829 
1830  // now use these lists to obtain the indices into the
1831  // linear algebra entities. This assumes an order.
1832  // For the matrix indices, first do the rows.
1833 
1834  li_Drain = extLIDVec[0];
1835  li_Gate = extLIDVec[1];
1836  li_Source = extLIDVec[2];
1837  li_Bulk = extLIDVec[3];
1838 
1839  int intLoc = 0;
1840  li_DrainPrime = intLIDVec[intLoc++];
1841 
1842  if( gateCond )
1843  li_GatePrime = intLIDVec[intLoc++];
1844  else
1846 
1847  if( sourceCond )
1848  li_SourcePrime = intLIDVec[intLoc++];
1849  else
1851 
1852  if( drainCond )
1853  li_DrainDrift = intLIDVec[intLoc++];
1854  else
1856 
1857  if( model_.D1DIOconductance )
1858  li_D1Prime = intLIDVec[intLoc++];
1859  else
1861 
1862  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS))
1863  {
1864  Xyce::dout() << "\n variable local indices:\n";
1865  Xyce::dout() << " li_Drain = " << li_Drain << std::endl;
1866  Xyce::dout() << " li_Source = " << li_Source << std::endl;
1867  Xyce::dout() << " li_Gate = " << li_Gate << std::endl;
1868  Xyce::dout() << " li_Bulk = " << li_Bulk << std::endl;
1869  Xyce::dout() << " li_DrainPrime = " << li_DrainPrime << std::endl;
1870  Xyce::dout() << " li_GatePrime = " << li_GatePrime << std::endl;
1871  Xyce::dout() << " li_SourcePrime = " << li_SourcePrime << std::endl;
1872  Xyce::dout() << " li_DrainDrift = " << li_DrainDrift << std::endl;
1873  Xyce::dout() << " li_D1Prime = " << li_D1Prime << std::endl;
1874 
1875  Xyce::dout() << section_divider << std::endl;
1876  }
1877 
1878 }
1879 
1880 //-----------------------------------------------------------------------------
1881 // Function : Instance::loadNodeSymbols
1882 // Purpose :
1883 // Special Notes :
1884 // Scope : public
1885 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
1886 // Creation Date : 05/13/05
1887 //-----------------------------------------------------------------------------
1888 void Instance::loadNodeSymbols(Util::SymbolTable &symbol_table) const
1889 {
1890  addInternalNode(symbol_table, li_DrainPrime, getName(), "drainprime");
1891 
1892  if (li_GatePrime != li_Gate )
1893  addInternalNode(symbol_table, li_GatePrime, getName(), "gateprime");
1894 
1895  if (li_SourcePrime != li_Source )
1896  addInternalNode(symbol_table, li_SourcePrime, getName(), "sourceprime");
1897 
1898  if (li_DrainDrift != li_DrainPrime )
1899  addInternalNode(symbol_table, li_DrainDrift, getName(), "draindrift");
1900 
1901  if (li_D1Prime != li_Source )
1902  addInternalNode(symbol_table, li_D1Prime, getName(), "d1pos");
1903 
1904  if (loadLeadCurrent)
1905  {
1906  addStoreNode(symbol_table, li_store_dev_id, getName(), "DEV_ID");
1907  addStoreNode(symbol_table, li_store_dev_is, getName(), "DEV_IS");
1908  addStoreNode(symbol_table, li_store_dev_ig, getName(), "DEV_IG");
1909  addStoreNode(symbol_table, li_store_dev_ib, getName(), "DEV_IB");
1910  }
1911 }
1912 
1913 //-----------------------------------------------------------------------------
1914 // Function : Instance::registerStateLIDs
1915 // Purpose :
1916 // Special Notes :
1917 // Scope : public
1918 // Creator : Robert Hoekstra, Computational Sciences
1919 // Creation Date : 6/21/02
1920 //-----------------------------------------------------------------------------
1921 void Instance::registerStateLIDs( const std::vector<int> & staLIDVecRef )
1922 {
1923  AssertLIDs(staLIDVecRef.size() == numStateVars);
1924 
1925  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS))
1926  {
1927  Xyce::dout() << std::endl;
1928  Xyce::dout() << section_divider << std::endl;
1929  Xyce::dout() << " In Instance::registerStateLIDs\n\n";
1930  Xyce::dout() << " name = " << getName() << std::endl;
1931  Xyce::dout() << " Number of State LIDs: " << numStateVars << std::endl;
1932  }
1933 
1934  // Copy over the global ID lists:
1935  staLIDVec = staLIDVecRef;
1936 
1937  int lid=0;
1938  li_state_vbdd = staLIDVec[lid++];
1939  li_state_vbs = staLIDVec[lid++];
1940  li_state_vgps = staLIDVec[lid++];
1941  li_state_vdds = staLIDVec[lid++];
1942  li_state_D1vd = staLIDVec[lid++];
1943 
1944  li_state_qgs = staLIDVec[lid++];
1945  //li_state_cqgs = staLIDVec[lid++];
1946  li_state_qgdd = staLIDVec[lid++];
1947  //li_state_cqgdd = staLIDVec[lid++];
1948  li_state_qgb = staLIDVec[lid++];
1949  //li_state_cqgb = staLIDVec[lid++];
1950 
1951  li_state_capgs = staLIDVec[lid++];
1952  li_state_capgdd = staLIDVec[lid++];
1953  li_state_capgb = staLIDVec[lid++];
1954 
1955  li_state_qbd = staLIDVec[lid++];
1956  //li_state_cqbd = staLIDVec[lid++];
1957  li_state_qbs = staLIDVec[lid++];
1958  //li_state_cqbs = staLIDVec[lid++];
1959 
1961  //li_state_D1DIOcapCurrent = staLIDVec[lid++];
1962 
1963  li_state_von = staLIDVec[lid++];
1964 
1965  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS))
1966  {
1967  Xyce::dout() << " State local indices:" << std::endl;
1968  Xyce::dout() << std::endl;
1969 
1970  Xyce::dout() << " li_state_vbdd = " << li_state_vbdd << "\n";
1971  Xyce::dout() << " li_state_vbs = " << li_state_vbs << "\n";
1972  Xyce::dout() << " li_state_vgps = " << li_state_vgps << "\n";
1973  Xyce::dout() << " li_state_vdds = " << li_state_vdds << "\n";
1974  Xyce::dout() << " li_state_qgs = " << li_state_qgs << "\n";
1975  //Xyce::dout() << " li_state_cqgs = " << li_state_cqgs << "\n";
1976  Xyce::dout() << " li_state_capgs = " << li_state_capgs << "\n";
1977  Xyce::dout() << " li_state_capgdd = " << li_state_capgdd << "\n";
1978  Xyce::dout() << " li_state_capgb = " << li_state_capgb << "\n";
1979  Xyce::dout() << " li_state_qgdd = " << li_state_qgdd << "\n";
1980  //Xyce::dout() << " li_state_cqgdd = " << li_state_cqgdd << "\n";
1981  Xyce::dout() << " li_state_qgb = " << li_state_qgb << "\n";
1982  //Xyce::dout() << " li_state_cqgb = " << li_state_cqgb << "\n";
1983  Xyce::dout() << " li_state_qbs = " << li_state_qbs << "\n";
1984  //Xyce::dout() << " li_state_cqbs = " << li_state_cqbs << "\n";
1985  Xyce::dout() << " li_state_qbd = " << li_state_qbd << "\n";
1986  //Xyce::dout() << " li_state_cqbd = " << li_state_cqbd << "\n";
1987  Xyce::dout() << " li_state_D1DIOcapCharge = " << li_state_D1DIOcapCharge << "\n";
1988  //Xyce::dout() << " li_state_D1DIOcapCurrent = " << li_state_D1DIOcapCurrent << "\n";
1989  Xyce::dout() << " li_state_von = " << li_state_von << "\n";
1990 
1991  Xyce::dout() << std::endl;
1992  Xyce::dout() << section_divider << std::endl;
1993  }
1994 
1995 }
1996 
1997 
1998 //-----------------------------------------------------------------------------
1999 // Function : Instance::registerStoreLIDs
2000 // Purpose :
2001 // Special Notes :
2002 // Scope : public
2003 // Creator : Richard Schiek, Electrical Systems Modeling
2004 // Creation Date : 4/23/2013
2005 //-----------------------------------------------------------------------------
2006 void Instance::registerStoreLIDs( const std::vector<int> & stoLIDVecRef )
2007 {
2008  AssertLIDs(stoLIDVecRef.size() == getNumStoreVars());
2009 
2010  // Copy over the global ID lists:
2011  stoLIDVec = stoLIDVecRef;
2012  if( loadLeadCurrent )
2013  {
2014  li_store_dev_id = stoLIDVec[0];
2015  li_store_dev_ig = stoLIDVec[1];
2016  li_store_dev_is = stoLIDVec[2];
2017  li_store_dev_ib = stoLIDVec[3];
2018  }
2019 }
2020 
2021 //-----------------------------------------------------------------------------
2022 // Function : Instance::jacobianStamp
2023 // Purpose :
2024 // Special Notes :
2025 // Scope : public
2026 // Creator : Robert Hoekstra, Computational Sciences
2027 // Creation Date : 9/3/02
2028 //-----------------------------------------------------------------------------
2029 const std::vector< std::vector<int> > & Instance::jacobianStamp() const
2030 {
2031  if( drainCond != 0.0 && gateCond != 0.0 && sourceCond != 0.0)
2033  else if( drainCond != 0.0 && gateCond != 0.0 && sourceCond == 0.0 )
2035  else if( drainCond == 0.0 && gateCond != 0.0 && sourceCond != 0.0 )
2037  else if( drainCond != 0.0 && gateCond == 0.0 && sourceCond != 0.0 )
2039  else if( drainCond == 0.0 && gateCond != 0.0 && sourceCond == 0.0 )
2041  else if( drainCond != 0.0 && gateCond == 0.0 && sourceCond == 0.0 )
2043  else if( drainCond == 0.0 && gateCond == 0.0 && sourceCond != 0.0 )
2045  else if( drainCond == 0.0 && gateCond == 0.0 && sourceCond == 0.0 )
2047  else
2049 }
2050 
2051 //-----------------------------------------------------------------------------
2052 // Function : Instance::registerJacLIDs
2053 // Purpose :
2054 // Special Notes :
2055 // Scope : public
2056 // Creator : Robert Hoekstra, Computational Sciences
2057 // Creation Date : 9/3/02
2058 //-----------------------------------------------------------------------------
2059 void Instance::registerJacLIDs( const std::vector< std::vector<int> > & jacLIDVec )
2060 {
2061  DeviceInstance::registerJacLIDs( jacLIDVec );
2062  std::vector<int> map;
2063  std::vector< std::vector<int> > map2;
2064  double d1dcond=(model_.D1DIOconductance);
2065  if (gateCond != 0.0)
2066  {
2067  if (drainCond != 0.0)
2068  {
2069  if (sourceCond != 0.0)
2070  {
2071  if (d1dcond != 0)
2072  {
2073  map = jacMap_D1C_DC_SC_GC;
2074  map2 = jacMap2_D1C_DC_SC_GC;
2075  }
2076  else
2077  {
2078  map = jacMap_DC_SC_GC;
2079  map2 = jacMap2_DC_SC_GC;
2080  }
2081  }
2082  else
2083  {
2084  if (d1dcond != 0)
2085  {
2086  map = jacMap_D1C_DC_GC;
2087  map2 = jacMap2_D1C_DC_GC;
2088  }
2089  else
2090  {
2091  map = jacMap_DC_GC;
2092  map2 = jacMap2_DC_GC;
2093  }
2094  }
2095  }
2096  else
2097  {
2098  if (sourceCond != 0.0)
2099  {
2100  if (d1dcond != 0)
2101  {
2102  map = jacMap_D1C_SC_GC;
2103  map2 = jacMap2_D1C_SC_GC;
2104  }
2105  else
2106  {
2107  map = jacMap_SC_GC;
2108  map2 = jacMap2_SC_GC;
2109  }
2110  }
2111  else
2112  {
2113  if (d1dcond != 0)
2114  {
2115  map = jacMap_D1C_GC;
2116  map2 = jacMap2_D1C_GC;
2117  }
2118  else
2119  {
2120  map = jacMap_GC;
2121  map2 = jacMap2_GC;
2122  }
2123  }
2124  }
2125  }
2126  else
2127  {
2128  if (drainCond != 0.0)
2129  {
2130  if (sourceCond != 0.0)
2131  {
2132  if (d1dcond != 0)
2133  {
2134  map = jacMap_D1C_DC_SC;
2135  map2 = jacMap2_D1C_DC_SC;
2136  }
2137  else
2138  {
2139  map = jacMap_DC_SC;
2140  map2 = jacMap2_DC_SC;
2141  }
2142  }
2143  else
2144  {
2145  if (d1dcond != 0)
2146  {
2147  map = jacMap_D1C_DC;
2148  map2 = jacMap2_D1C_DC;
2149  }
2150  else
2151  {
2152  map = jacMap_DC;
2153  map2 = jacMap2_DC;
2154  }
2155  }
2156  }
2157  else
2158  {
2159  if (sourceCond != 0.0)
2160  {
2161  if (d1dcond != 0)
2162  {
2163  map = jacMap_D1C_SC;
2164  map2 = jacMap2_D1C_SC;
2165  }
2166  else
2167  {
2168  map = jacMap_SC;
2169  map2 = jacMap2_SC;
2170  }
2171  }
2172  else
2173  {
2174  if (d1dcond != 0)
2175  {
2176  map = jacMap_D1C;
2177  map2 = jacMap2_D1C;
2178  }
2179  else
2180  {
2181  map = jacMap;
2182  map2 = jacMap2;
2183  }
2184  }
2185  }
2186  }
2187 
2188  ADrainEquDrainNodeOffset = jacLIDVec[map[0]][map2[0][0]];
2189  ADrainEquSourceNodeOffset = jacLIDVec[map[0]][map2[0][1]];
2190  ADrainEquDrainDriftNodeOffset = jacLIDVec[map[0]][map2[0][2]];
2191  ADrainEquD1PrimeNodeOffset = jacLIDVec[map[0]][map2[0][3]];
2192 
2193  AGateEquGateNodeOffset = jacLIDVec[map[1]][map2[1][0]];
2194  AGateEquGatePrimeNodeOffset = jacLIDVec[map[1]][map2[1][1]];
2195 
2196  ASourceEquDrainNodeOffset = jacLIDVec[map[2]][map2[2][0]];
2197  ASourceEquSourceNodeOffset = jacLIDVec[map[2]][map2[2][1]];
2198  ASourceEquSourcePrimeNodeOffset = jacLIDVec[map[2]][map2[2][2]];
2199  ASourceEquD1PrimeNodeOffset = jacLIDVec[map[2]][map2[2][3]];
2200 
2201  ABulkEquBulkNodeOffset = jacLIDVec[map[3]][map2[3][0]];
2202  ABulkEquDrainPrimeNodeOffset = jacLIDVec[map[3]][map2[3][1]];
2203  ABulkEquGatePrimeNodeOffset = jacLIDVec[map[3]][map2[3][2]];
2204  ABulkEquSourcePrimeNodeOffset = jacLIDVec[map[3]][map2[3][3]];
2205 
2206  ADrainPrimeEquBulkNodeOffset = jacLIDVec[map[4]][map2[4][0]];
2207  ADrainPrimeEquDrainPrimeNodeOffset = jacLIDVec[map[4]][map2[4][1]];
2208  ADrainPrimeEquGatePrimeNodeOffset = jacLIDVec[map[4]][map2[4][2]];
2209  ADrainPrimeEquSourcePrimeNodeOffset = jacLIDVec[map[4]][map2[4][3]];
2210  ADrainPrimeEquDrainDriftNodeOffset = jacLIDVec[map[4]][map2[4][4]];
2211 
2212  AGatePrimeEquGateNodeOffset = jacLIDVec[map[5]][map2[5][0]];
2213  AGatePrimeEquBulkNodeOffset = jacLIDVec[map[5]][map2[5][1]];
2214  AGatePrimeEquDrainPrimeNodeOffset = jacLIDVec[map[5]][map2[5][2]];
2215  AGatePrimeEquGatePrimeNodeOffset = jacLIDVec[map[5]][map2[5][3]];
2216  AGatePrimeEquSourcePrimeNodeOffset = jacLIDVec[map[5]][map2[5][4]];
2217 
2218  ASourcePrimeEquSourceNodeOffset = jacLIDVec[map[6]][map2[6][0]];
2219  ASourcePrimeEquBulkNodeOffset = jacLIDVec[map[6]][map2[6][1]];
2220  ASourcePrimeEquDrainPrimeNodeOffset = jacLIDVec[map[6]][map2[6][2]];
2221  ASourcePrimeEquGatePrimeNodeOffset = jacLIDVec[map[6]][map2[6][3]];
2222  ASourcePrimeEquSourcePrimeNodeOffset = jacLIDVec[map[6]][map2[6][4]];
2223 
2224  ADrainDriftEquDrainNodeOffset = jacLIDVec[map[7]][map2[7][0]];
2225  ADrainDriftEquDrainPrimeNodeOffset = jacLIDVec[map[7]][map2[7][1]];
2226  ADrainDriftEquDrainDriftNodeOffset = jacLIDVec[map[7]][map2[7][2]];
2227 
2228  AD1PrimeEquDrainNodeOffset = jacLIDVec[map[8]][map2[8][0]];
2229  AD1PrimeEquSourceNodeOffset = jacLIDVec[map[8]][map2[8][1]];
2230  AD1PrimeEquD1PrimeNodeOffset = jacLIDVec[map[8]][map2[8][2]];
2231 }
2232 
2233 //-----------------------------------------------------------------------------
2234 // Function : Instance::setupPointers
2235 // Purpose :
2236 // Special Notes :
2237 // Scope : public
2238 // Creator : Eric Keiter, SNL
2239 // Creation Date : 12/06/08
2240 //-----------------------------------------------------------------------------
2242 {
2243 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
2244  Linear::Matrix & dFdx = *(extData.dFdxMatrixPtr);
2245  Linear::Matrix & dQdx = *(extData.dQdxMatrixPtr);
2246 
2251 
2254 
2259 
2264 
2270 
2276 
2282 
2286 
2290 
2291 
2292 
2297 
2300 
2305 
2310 
2316 
2322 
2328 
2332 
2336 
2337 
2338 
2339 
2340 #endif
2341 }
2342 
2343 //-----------------------------------------------------------------------------
2344 // Function : Instance::loadDAEQVector
2345 //
2346 // Purpose : Loads the Q-vector contributions for a single
2347 // voltage source instance.
2348 //
2349 // Special Notes : The "Q" vector is part of a standard DAE formalism in
2350 // which the system of equations is represented as:
2351 //
2352 // f(x) = dQ(x)/dt + F(x) + B(t) = 0
2353 //
2354 // The "Q" vector contains charges and fluxes, mostly.
2355 // The voltage source will not make any contributions to Q,
2356 // so this function does nothing.
2357 //
2358 // Scope : public
2359 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
2360 // Creation Date : 06/02/05
2361 //-----------------------------------------------------------------------------
2363 {
2364  double * qVec = extData.daeQVectorRawPtr;
2365  double * dQdxdVp = extData.dQdxdVpVectorRawPtr;
2366 
2367  double coef_Jdxp(0.0);
2368 
2369  double ceqbs(0.0), ceqbd(0.0), ceqgb(0.0), ceqgs(0.0), ceqgdd(0.0); // 3f5 vars
2370  double Dtype =model_.dtype;
2371 
2372  ceqbs = Dtype*qbs;
2373  ceqbd = Dtype*qbd;
2374  ceqgb = Dtype*qgb;
2375  ceqgs = Dtype*qgs;
2376  ceqgdd = Dtype*qgdd;
2377 
2378  qVec[li_Bulk] += (ceqbs + ceqbd - ceqgb);
2379  qVec[li_DrainPrime] += -(ceqbd + ceqgdd);
2380  qVec[li_GatePrime] += (ceqgs+ceqgdd+ceqgb);
2381  qVec[li_SourcePrime] += (-(ceqbs + ceqgs));
2382  qVec[li_D1Prime] += D1DIOcapCharge;
2383  qVec[li_Drain] += -D1DIOcapCharge;
2384  // voltlim section:
2385  if (!origFlag)
2386  {
2387  // bulk
2388  coef_Jdxp = Dtype*(-Capgb*(vgps-vgps_orig-vbs+vbs_orig)
2389  + (+Capgb)*(vbdd-vbdd_orig)
2390  + (+capbs)*(vbs-vbs_orig));
2391  dQdxdVp[li_Bulk] += coef_Jdxp;
2392 
2393  // drain-prime
2394  coef_Jdxp = Dtype*(-Capgdd*(vgpdd-vgpdd_orig)-capbd*(vbdd-vbdd_orig));
2395  dQdxdVp[li_DrainPrime] += coef_Jdxp;
2396 
2397  // gate-prime
2398  coef_Jdxp = Dtype*(Capgdd*(vgpdd-vgpdd_orig)+Capgs*(vgps-vgps_orig)+
2400  dQdxdVp[li_GatePrime] += coef_Jdxp;
2401 
2402  // source-prime
2403  coef_Jdxp = Dtype*(-Capgs*(vgps-vgps_orig)-capbs*(vbs-vbs_orig));
2404  dQdxdVp[li_SourcePrime] += coef_Jdxp;
2405 
2406  coef_Jdxp = -D1capd*(D1vd-D1vd_orig);
2407  dQdxdVp[li_D1Prime] -= coef_Jdxp;
2408  dQdxdVp[li_Drain] += coef_Jdxp;
2409  }
2410 
2411  if( loadLeadCurrent )
2412  {
2413  double * storeLeadQ = extData.storeLeadCurrQCompRawPtr;
2414 
2415  storeLeadQ[li_store_dev_id] = -D1DIOcapCharge;
2416  storeLeadQ[li_store_dev_is] = 0;
2417  storeLeadQ[li_store_dev_ig] = 0;
2418  storeLeadQ[li_store_dev_ib] = (ceqbs + ceqbd - ceqgb);
2419  // case where optional nodes become external nodes:
2420  if( !gateCond )
2421  {
2422  // G' is G
2423  storeLeadQ[li_store_dev_ig] += (ceqgs+ceqgdd+ceqgb);
2424  }
2425 
2426  if( !sourceCond )
2427  {
2428  // S' is S
2429  storeLeadQ[li_store_dev_is] += (-(ceqbs + ceqgs));
2430  }
2431 
2432  if( !drainCond )
2433  {
2434  // Ddrift is D'
2435  }
2436 
2437  if( !model_.D1DIOconductance )
2438  {
2439  // D1' is S
2440  storeLeadQ[li_store_dev_is] += D1DIOcapCharge;
2441  }
2442  }
2443 
2444  return true;
2445 }
2446 
2447 //-----------------------------------------------------------------------------
2448 // Function : Instance::loadDAEFVector
2449 //
2450 // Purpose : Loads the F-vector contributions for a single
2451 // VDMOS instance.
2452 //
2453 // Special Notes :
2454 //
2455 // Scope : public
2456 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
2457 // Creation Date : 06/02/05
2458 //-----------------------------------------------------------------------------
2460 {
2461  double * fVec = extData.daeFVectorRawPtr;
2462  double * dFdxdVp = extData.dFdxdVpVectorRawPtr;
2463 
2464  double coef_Jdxp(0.0), gd_Jdxp(0.0);
2465  double gmin1 = getDeviceOptions().gmin;
2466  double Dtype =model_.dtype;
2467 
2468  double ceqbs = Dtype*(cbs+cqbs);
2469  double ceqbd = Dtype*(cbd+cqbd);
2470 
2471  double D1current = Dtype*D1cdeq; // don't add in the capacitor stuff
2472  fVec[li_Drain] += (Idraindrift + Irdsshunt - D1current);
2473 
2474  if (Igate != 0.0)
2475  {
2476  fVec[li_Gate] += Igate;
2477  fVec[li_GatePrime] += -Igate;
2478  }
2479  fVec[li_Source] += (Isource - Irdsshunt + Ird1rs);
2480  fVec[li_Bulk] += (ceqbs + ceqbd);
2481  fVec[li_DrainPrime] += (-Idrain-(ceqbd - cdreq));
2482  fVec[li_SourcePrime] += (-Isource-(ceqbs + cdreq));
2483  fVec[li_DrainDrift] += (Idrain - Idraindrift);
2484  fVec[li_D1Prime] += (D1current - Ird1rs);
2485 
2486  // limiter terms:
2487  if (!origFlag)
2488  {
2489  // bulk
2490  coef_Jdxp = Dtype*( + ((gbd-gmin1))*(vbdd-vbdd_orig)
2491  + ((gbs-gmin1))*(vbs-vbs_orig));
2492  dFdxdVp[li_Bulk] += coef_Jdxp;
2493 
2494  // drain-prime
2495  coef_Jdxp = Dtype*(-((gbd-gmin1))*
2497  +Gm*((mode>0)?(vgps-vgps_orig):(vgpdd-vgpdd_orig))
2498  +Gmbs*((mode>0)?(vbs-vbs_orig):(vbdd-vbdd_orig)));
2499  dFdxdVp[li_DrainPrime] += coef_Jdxp;
2500 
2501  // source prime
2502  coef_Jdxp = Dtype*(-((gbs-gmin1))*(vbs-vbs_orig)
2503  -gdds*(vdds-vdds_orig)
2504  -Gm*((mode>0)?(vgps-vgps_orig):(vgpdd-vgpdd_orig))
2505  -Gmbs*((mode>0)?(vbs-vbs_orig):(vbdd-vbdd_orig)));
2506  dFdxdVp[li_SourcePrime] += coef_Jdxp;
2507 
2508  gd_Jdxp = -D1gd * (D1vd-D1vd_orig);
2509  dFdxdVp[li_Drain] += gd_Jdxp;
2510  dFdxdVp[li_D1Prime] -= gd_Jdxp;
2511  }
2512 
2513  if( loadLeadCurrent )
2514  {
2515  double * storeLeadF = extData.nextStoVectorRawPtr;
2516 
2517  storeLeadF[li_store_dev_id] = (Idraindrift + Irdsshunt - D1current);
2518  storeLeadF[li_store_dev_is] = (Isource - Irdsshunt + Ird1rs);
2519  storeLeadF[li_store_dev_ig] = 0;
2520  storeLeadF[li_store_dev_ib] = (ceqbs + ceqbd);
2521  // case where optional nodes become external nodes:
2522  if (Igate != 0.0)
2523  {
2524  storeLeadF[li_store_dev_ig] += Igate;
2525  }
2526  if( !gateCond )
2527  {
2528  // G' is G
2529  storeLeadF[li_store_dev_ig] += -Igate;
2530  }
2531 
2532  if( !sourceCond )
2533  {
2534  // S' is S
2535  storeLeadF[li_store_dev_is] += (-Isource-(ceqbs + cdreq));
2536  }
2537 
2538  if( !drainCond )
2539  {
2540  // Ddrift is D'
2541  }
2542 
2543  if( !model_.D1DIOconductance )
2544  {
2545  // D1' is S
2546  storeLeadF[li_store_dev_is] += (D1current - Ird1rs);
2547  }
2548  }
2549 
2550  return true;
2551 }
2552 
2553 //-----------------------------------------------------------------------------
2554 // Function : Instance::loadDAEdQdx
2555 //
2556 // Purpose : Loads the Q-vector contributions for a single
2557 // VDMOS instance.
2558 //
2559 // Special Notes : The "Q" vector is part of a standard DAE formalism in
2560 // which the system of equations is represented as:
2561 //
2562 // f(x) = dQ(x)/dt + F(x) - B(t) = 0
2563 //
2564 // Scope : public
2565 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
2566 // Creation Date : 06/02/05
2567 //-----------------------------------------------------------------------------
2569 {
2570  Linear::Matrix & dQdx = *(extData.dQdxMatrixPtr);
2571 
2572  if (getSolverState().dcopFlag) return true;
2573 
2578 
2582 
2587 
2591 
2594 
2597 
2598  return true;
2599 }
2600 
2601 //-----------------------------------------------------------------------------
2602 // Function : Instance::loadDAEdFdx ()
2603 //
2604 // Purpose : Loads the F-vector contributions for a single
2605 // instance.
2606 //
2607 // Special Notes : The F-vector is an algebraic constaint.
2608 //
2609 // Scope : public
2610 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
2611 // Creation Date : 06/02/05
2612 //-----------------------------------------------------------------------------
2614 {
2615  Linear::Matrix & dFdx = *(extData.dFdxMatrixPtr);
2616 
2621 
2622  if (gateCond != 0)
2623  {
2626  }
2627 
2631 
2632  if (sourceCond != 0.0)
2633  {
2635  }
2636 
2639 
2641 
2646  if (drainCond != 0.0)
2647  {
2649  }
2650 
2651  if (gateCond != 0)
2652  {
2655  }
2656 
2657  if (sourceCond != 0.0)
2658  {
2660  }
2665 
2667 
2668  if (drainCond != 0.0)
2669  {
2671  }
2673 
2677 
2678  return true;
2679 }
2680 
2681 //-----------------------------------------------------------------------------
2682 // Function : Instance::updateIntermediateVars
2683 // Purpose :
2684 // Special Notes :
2685 // Scope : public
2686 // Creator : pmc
2687 // Creation Date : 2/17/2004
2688 //-----------------------------------------------------------------------------
2690 {
2691  bool bsuccess = true;
2692 
2693  double * solVec = extData.nextSolVectorRawPtr;
2694 
2695  // 3f5 likes to use the same variable names in local variables and in
2696  // structures. Messes with us! Define some local versions with capitals
2697  // instead
2698 #define ISUBMOD model_.isubmod
2699 
2700  double Von(0.0);
2701  double Vddsat(0.0);
2702  //
2703  double evbs(0.0);
2704  double evbdd(0.0);
2705  double sarg(0.0);
2706  double sargsw(0.0);
2707  double arg(0.0);
2708  int Check(1);
2709 
2710  double capgs_old(0.0);
2711  double capgdd_old(0.0);
2712  double capgb_old(0.0);
2713 
2714  // diode #1
2715  double D1temp(0.0);
2716  double D1power(0.0);
2717  double D1arg(0.0);
2718  double D1cd(0.0);
2719  double D1cdr(0.0);
2720  double D1csat(0.0);
2721  double D1csatr(0.0);
2722  double D1czero(0.0);
2723  double D1czof2(0.0);
2724  double D1evd(0.0);
2725  double D1evr(0.0);
2726  double D1isr(0.0);
2727  double D1evrev(0.0);
2728  double D1gdr(0.0);
2729  double D1sarg(0.0);
2730  double D1vdtemp(0.0);
2731  double D1vtr(0.0);
2732  int D1Check(0);
2733 
2734  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
2735  {
2736  Xyce::dout() << subsection_divider << std::endl;
2737  Xyce::dout() <<" Instance::updateIntermediateVars.\n"<<std::endl;
2738  Xyce::dout() <<" name = " << getName() << std::endl;
2739  Xyce::dout() <<" Model name = " << model_.getName() << std::endl;
2740  Xyce::dout() <<" dtype is " << model_.dtype << std::endl;
2741  Xyce::dout() << std::endl;
2742  Xyce::dout().width(25); Xyce::dout().precision(17); Xyce::dout().setf(std::ios::scientific);
2743  }
2744 
2745  if( (tSatCurDens == 0) || (drainArea == 0) || (sourceArea == 0))
2746  {
2747  DrainSatCur = tSatCur;
2749  }
2750  else
2751  {
2754  }
2755 
2756  // we need our solution variables for any of this stuff
2757 
2758  Vd = 0.0;
2759  Vs = 0.0;
2760  Vg = 0.0;
2761  Vb = 0.0;
2762  Vdp = 0.0;
2763  Vgp = 0.0;
2764  Vsp = 0.0;
2765  Vdd = 0.0;
2766 
2767  Vd = solVec[li_Drain];
2768  Vg = solVec[li_Gate];
2769  Vs = solVec[li_Source];
2770  Vb = solVec[li_Bulk];
2771  Vsp = solVec[li_SourcePrime];
2772  Vgp = solVec[li_GatePrime];
2773  Vdp = solVec[li_DrainPrime];
2774  Vdd = solVec[li_DrainDrift];
2775 
2776  // node for diode #1
2777  Vd1p = solVec[li_D1Prime];
2778 
2781  D1vtr = model_.D1DIOnr * D1vt;
2782 
2783  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
2784  {
2785  Xyce::dout() << " " << std::endl;
2786  Xyce::dout() << " Solution vector: " << std::endl;
2787  Xyce::dout() << " Vg = " << Vg << std::endl;
2788  Xyce::dout() << " Vd = " << Vd << std::endl;
2789  Xyce::dout() << " Vs = " << Vs << std::endl;
2790  Xyce::dout() << " Vb = " << Vb << std::endl;
2791  Xyce::dout() << " Vdp = " << Vdp << std::endl;
2792  Xyce::dout() << " Vgp = " << Vgp << std::endl;
2793  Xyce::dout() << " Vsp = " << Vsp << std::endl;
2794  Xyce::dout() << " Vdd = " << Vdd << std::endl;
2795  Xyce::dout() << " Vd1p= " << Vd1p << std::endl;
2796  }
2797 
2798  // voltage drops
2799  Vddp = Vd - Vdp;
2800  Vddd = Vd - Vdd;
2801  Vdddp = Vdd - Vdp;
2802  Vssp = Vs - Vsp;
2803  Vbsp = Vb - Vsp;
2804  Vbdp = Vb - Vdp;
2805 
2806  Vggp = Vg - Vgp;
2807  Vgpsp = Vgp - Vsp;
2808  Vgpdp = Vgp - Vdp;
2809  Vgpb = Vgp - Vb;
2810  Vdpsp = Vdp - Vsp;
2811 
2812  // setup corrections for dtype.
2813  vbs = model_.dtype * Vbsp;
2814  vgps = model_.dtype * Vgpsp;
2815  vdds = model_.dtype * Vdpsp;
2816  D1vd = model_.dtype * (Vd1p - Vd);
2817 
2818  vbdd = vbs-vdds;
2819  vgpdd = vgps-vdds;
2820 
2821  // set up the orig voltages.
2822  origFlag = 1;
2823  limitedFlag = false;
2824  vgps_orig = vgps;
2825  vdds_orig = vdds;
2826  vbs_orig = vbs;
2827  vbdd_orig = vbdd;
2828  vgpdd_orig = vgpdd;
2829  D1vd_orig = D1vd;
2830 
2831  if (getSolverState().newtonIter == 0)
2832  {
2833 
2835  {
2837  {
2838  Linear::Vector * flagSolVectorPtr = extData.flagSolVectorPtr;
2839  if ((*flagSolVectorPtr)[li_Drain] == 0 || (*flagSolVectorPtr)[li_Gate] == 0 ||
2840  (*flagSolVectorPtr)[li_Source] == 0 || (*flagSolVectorPtr)[li_SourcePrime] ||
2841  (*flagSolVectorPtr)[li_DrainDrift] == 0 || (*flagSolVectorPtr)[li_GatePrime] ||
2842  (*flagSolVectorPtr)[li_DrainPrime] || (*flagSolVectorPtr)[li_Bulk] )
2843  {
2844  vbs = -1;
2845  vgps = model_.dtype*tVto;
2846  vdds = 0;
2847  vbdd = vbs-vdds;
2848  vgpdd = vgps-vdds;
2849  }
2850  }
2851  else
2852  {
2853  vbs = -1;
2854  vgps = model_.dtype*tVto;
2855  vdds = 0;
2856  vbdd = vbs-vdds;
2857  vgpdd = vgps-vdds;
2858  }
2859  }
2860 
2861  ////////////////////////////////////////////////////////////////////////
2862  // Note: the "old" variables should be the values for the previous
2863  // Newton iteration. That previous newton iteration could
2864  // have happened in a previous time step...
2866  {
2872  Von = model_.dtype *
2874  }
2875  else
2876  { // there's no history
2877  vbs_old = vbs;
2878  vbdd_old = vbdd;
2879  vgps_old = vgps;
2880  vdds_old = vdds;
2881  D1vd_old = D1vd;
2882  }
2883  }
2884  else
2885  {
2891  Von = model_.dtype *
2894  }
2895 
2896  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
2897  {
2898  Xyce::dout() << "After Von first set, ";
2899  Xyce::dout() << " von = " << von << " Von = " << Von << std::endl;
2900  }
2901 
2902  ////////////////////////////////////////////
2903  // SPICE-type Voltage Limiting (PINNING)
2904  ////////////////////////////////////////////
2905  if (getDeviceOptions().voltageLimiterFlag)
2906  {
2907  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
2908  {
2909  Xyce::dout() << " checking whether to limit voltages "<< std::endl;
2910  Xyce::dout() << " Von = " << Von << std::endl;
2911  Xyce::dout() << " before limiting: " << std::endl;
2912  Xyce::dout() << " vgpdd = " << vgpdd << " vgpdd_old = " << vgpdd_old << std::endl;
2913  Xyce::dout() << " vgps = " << vgps << " vgps_old = " << vgps_old << std::endl;
2914  Xyce::dout() << " vdds = " << vdds << " vdds_old = " << vdds_old << std::endl;
2915  Xyce::dout() << " vbs = " << vbs << " vbs_old = " << vbs_old << std::endl;
2916  Xyce::dout() << " vbdd = " << vbdd << " vbdd_old = " << vbdd_old << std::endl;
2917  }
2918 
2919  if (vdds_old >= 0)
2920  {
2921  vgps = devSupport.fetlim( vgps, vgps_old, Von);
2922  vdds = vgps - vgpdd;
2923  vdds = devSupport.limvds( vdds, vdds_old);
2924  vgpdd = vgps - vdds;
2925  }
2926  else
2927  {
2929  vdds = vgps - vgpdd;
2930  vdds = -devSupport.limvds( -vdds, -vdds_old );
2931  vgps = vgpdd + vdds;
2932  }
2933 
2934  if (vdds >= 0.0)
2935  {
2936  vbs = devSupport.pnjlim( vbs, vbs_old, vt, sourceVcrit, &Check);
2937  vbdd = vbs - vdds;
2938  }
2939  else
2940  {
2941  vbdd = devSupport.pnjlim( vbdd, vbdd_old, vt, drainVcrit, &Check);
2942  vbs = vbdd + vdds;
2943  }
2944 
2945  if (Check == 1) limitedFlag=true;
2946  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
2947  {
2948  Xyce::dout() << " After limiting: " << std::endl;
2949  Xyce::dout() << " vgpdd = " << vgpdd << std::endl;
2950  Xyce::dout() << " vgps = " << vgps << std::endl;
2951  Xyce::dout() << " vdds = " << vdds << std::endl;
2952  Xyce::dout() << " vbs = " << vbs << std::endl;
2953  Xyce::dout() << " vbdd = " << vbdd << std::endl;
2954  }
2955 
2956  // limit diode junction voltage
2957  bool tmpGiven = false;
2958  tmpGiven = (model_.D1DIObreakdownVoltageGiven);
2959  if (tmpGiven && (D1vd < std::min(0.0, 10.0*D1vte-D1DIOtBrkdwnV) ))
2960  {
2961  D1vdtemp = -(D1vd + D1DIOtBrkdwnV);
2962  D1vdtemp = devSupport.pnjlim(D1vdtemp,
2964  D1DIOtVcrit, &D1Check);
2965  D1vd = -(D1vdtemp + D1DIOtBrkdwnV);
2966  }
2967  else
2968  {
2970  D1DIOtVcrit, &D1Check);
2971  }
2972 
2973  if (D1Check == 1) limitedFlag=true;
2974  } // voltage limiter flag
2975 
2976  ////
2977  // now all the preliminaries are over - we can start doing the
2978  // real work
2979  ////
2980  vbdd = vbs - vdds;
2981  vgpdd = vgps - vdds;
2982  Vgpb = vgps - vbs;
2983 
2984  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
2985  {
2986  Xyce::dout() << " vbs = " << vbs << std::endl;
2987  Xyce::dout() << " vgps = " << vgps << std::endl;
2988  Xyce::dout() << " vdds = " << vdds << std::endl;
2989  Xyce::dout() << " vbdd = " << vbdd << std::endl;
2990  Xyce::dout() << " vgpdd= " << vgpdd << std::endl;
2991 
2992  Xyce::dout() << " Vddp = " << Vddp << std::endl;
2993  Xyce::dout() << " Vddd = " << Vddd << std::endl;
2994  Xyce::dout() << " Vdddp= " << Vdddp << std::endl;
2995  Xyce::dout() << " Vssp = " << Vssp << std::endl;
2996  Xyce::dout() << " Vbsp = " << Vbsp << std::endl;
2997  Xyce::dout() << " Vbdp = " << Vbdp << std::endl;
2998  Xyce::dout() << " Vggp = " << Vggp << std::endl;
2999  Xyce::dout() << " Vgpsp = " << Vgpsp << std::endl;
3000  Xyce::dout() << " Vgpdp = " << Vgpdp << std::endl;
3001  Xyce::dout() << " Vgpb = " << Vgpb << std::endl;
3002  Xyce::dout() << " Vdpsp= " << Vdpsp << std::endl;
3003 
3004  }
3005 
3006  // Now set the origFlag
3007  if (vgps_orig != vgps || vdds_orig != vdds ||
3008  vbs_orig != vbs || vbdd_orig != vbdd || vgpdd_orig != vgpdd) origFlag = 0;
3009 
3010  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3011  {
3012  if (origFlag == 0)
3013  {
3014  Xyce::dout() << " Something modified the voltages. " << std::endl;
3015  Xyce::dout() << " Voltage before after diff " << std::endl;
3016  Xyce::dout() << " vgps " << vgps_orig << " " << vgps << " " << vgps-vgps_orig << std::endl;
3017  Xyce::dout() << " vdds " << vdds_orig << " " << vdds << " " << vdds-vdds_orig << std::endl;
3018  Xyce::dout() << " vbs " << vbs_orig << " " << vbs << " " << vbs-vbs_orig << std::endl;
3019  Xyce::dout() << " vbdd " << vbdd_orig << " " << vbdd << " " << vbdd-vbdd_orig << std::endl;
3020  Xyce::dout() << " vgpdd " << vgpdd_orig << " " << vgpdd << " " << vgpdd-vgpdd_orig << std::endl;
3021  }
3022  }
3023 
3024 
3025  ////
3026  // bulk-source and bulk-drain diodes
3027  // here we just evaluate the ideal diode current and the
3028  // corresponding derivative (conductance).
3029  ////
3030 
3031  double pi(0.0);
3032  double end10(0.0);
3033  double end1(0.0);
3034  double end20(0.0);
3035  double end2(0.0);
3036  double end3(0.0);
3037  double j1(0.0);
3038  double fr1(0.0);
3039  double fr2(0.0);
3040  double d1j1(0.0);
3041  double fr3(0.0);
3042  double d2j1(0.0);
3043  double dj1(0.0);
3044 
3045  if(model_.artd != 0.0 && sourceArea != 0.0)
3046  {
3048  end1=1+exp(std::min(CONSTMAX_EXP_ARG,end10/vt));
3050  end2=1+exp(std::min(CONSTMAX_EXP_ARG,end20/vt));
3051  end3=(model_.crtd-model_.nrtd*vbs)/model_.drtd;
3052  pi = 3.1415927;
3053  j1=model_.artd*log(end1/end2)*(pi/2+atan(end3));
3054  j1=j1*sourceArea;
3055 
3056  fr1=(end1-1)*model_.nrtd/(end1*vt);
3057  fr2=(end2-1)*model_.nrtd/(end2*vt);
3058  d1j1=(fr1-fr2)*(pi/2+atan(end3));
3059  fr3=(-model_.nrtd/model_.drtd)/(1+end3*end3);
3060  d2j1=log(end1/end2)*fr3;
3061  dj1=model_.artd*(d1j1+d2j1);
3062  dj1=dj1*sourceArea;
3063  }
3064  else
3065  {
3066  j1 = 0.0;
3067  dj1 = 0.0;
3068  }
3069 
3070  if(vbs <= 0)
3071  {
3074  cbs = gbs*vbs+j1;
3075 // PMC for D1 diode testing
3076 // gbs = 0.0;
3077 // cbs = 0.0;
3078 
3079  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3080  {
3081  Xyce::dout() << "*******Setting cbs for vbs<=0 ******" << std::endl;
3082  Xyce::dout() << " vbs = " << vbs << std::endl;
3083  Xyce::dout() << " vt = " << vt << std::endl;
3084  Xyce::dout() << " n2 = " << model_.n2 << std::endl;
3085  Xyce::dout() << " SSC = " << SourceSatCur << std::endl;
3086  Xyce::dout() << " gbs = " << gbs << std::endl;
3087  Xyce::dout() << " cbs = " << cbs << std::endl;
3088  Xyce::dout() << " j1 = " << j1 << std::endl;
3089  }
3090  }
3091  else
3092  {
3093  evbs = exp(std::min(CONSTMAX_EXP_ARG,model_.n2*vbs/vt));
3095  + getDeviceOptions().gmin + dj1);
3096  cbs = numberParallel*(SourceSatCur*(evbs-1) + j1);
3097 // PMC for D1 diode testing
3098 // gbs = 0.0;
3099 // cbs = 0.0;
3100 
3101  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3102  {
3103  Xyce::dout() << "*******Setting cbs for vbs>0 ******" << std::endl;
3104  Xyce::dout() << " vbs = " << vbs << std::endl;
3105  Xyce::dout() << " vt = " << vt << std::endl;
3106  Xyce::dout() << " n2 = " << model_.n2 << std::endl;
3107  Xyce::dout() << " vbs/vt = " << vbs/vt << std::endl;
3108  Xyce::dout() << " Maxarg = " << CONSTMAX_EXP_ARG << std::endl;
3109  Xyce::dout() << " arg = " << std::min(CONSTMAX_EXP_ARG,
3110  model_.n2*vbs/vt) << std::endl;
3111  Xyce::dout() << " evbs = " << evbs << std::endl;
3112  Xyce::dout() << " gbs = " << gbs << std::endl;
3113  Xyce::dout() << " cbs = " << cbs << std::endl;
3114  }
3115  }
3116  if(vbdd <= 0)
3117  {
3120  cbd = gbd *vbdd;
3121 // PMC for D1 diode testing
3122 // gbd = 0.0;
3123 // cbd = 0.0;
3124 
3125  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3126  {
3127  Xyce::dout() << "*******Setting cbd for vbdd<=0 ******" << std::endl;
3128  Xyce::dout() << " vbdd = " << vbdd << std::endl;
3129  Xyce::dout() << " SSC = " << SourceSatCur << std::endl;
3130  Xyce::dout() << " vt = " << vt << std::endl;
3131  Xyce::dout() << " n2 = " << model_.n2 << std::endl;
3132  Xyce::dout() << " gbd = " << gbd << std::endl;
3133  Xyce::dout() << " cbd = " << cbd << std::endl;
3134  }
3135 
3136  }
3137  else
3138  {
3139  evbdd = exp(std::min(CONSTMAX_EXP_ARG,model_.n2*vbdd/vt));
3141  + getDeviceOptions().gmin);
3142  cbd = numberParallel*DrainSatCur*(evbdd-1);
3143 // PMC for D1 diode testing
3144 // gbd = 0.0;
3145 // cbd = 0.0;
3146 
3147  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3148  {
3149  Xyce::dout() << "*******Setting cbd for vbdd>0 ******" << std::endl;
3150  Xyce::dout() << " vbdd = " << vbdd << std::endl;
3151  Xyce::dout() << " vt = " << vt << std::endl;
3152  Xyce::dout() << " vbdd/vt = " << vbdd/vt << std::endl;
3153  Xyce::dout() << " Maxarg = " << CONSTMAX_EXP_ARG << std::endl;
3154  Xyce::dout() << " arg = " << std::min(CONSTMAX_EXP_ARG,
3155  model_.n2*vbdd/vt) << std::endl;
3156  Xyce::dout() << " evbdd = " << evbdd << std::endl;
3157  Xyce::dout() << " gbd = " << gbd << std::endl;
3158  Xyce::dout() << " cbd = " << cbd << std::endl;
3159  }
3160  }
3161 
3162  if (vdds >= 0)
3163  mode = 1;
3164  else
3165  mode = -1;
3166 
3167  double dvonvbs(0.0);
3168  double lvbs = mode == 1 ? vbs : vbdd;
3169  UCCMcvon(lvbs, von, dvonvbs);
3170  Von = model_.dtype * von; // von contains dtype, Von does not
3171 
3172  if(!ISUBMOD)
3173  {
3174  // Old version of model with corrected Jacobian terms
3175  UCCMmosa1(vdds>0?vgps:vgpdd,mode*vdds,dvonvbs, cdraindrift, vddsat);
3176  }
3177  else
3178  {
3179  // Version of model with substrate currents
3180  UCCMmosa2(vdds>0?vgps:vgpdd,mode*vdds,dvonvbs, cdraindrift, vddsat);
3181  }
3182  Vddsat = model_.dtype * vddsat;
3183 
3184  ////
3185  // COMPUTE EQUIVALENT DRAIN CURRENT SOURCE
3186  ////
3187 
3188  // Substrate current calculations
3189  ISUB = cdraindrift*mm1;
3190  GMSUB = gm*mm1 + cdraindrift*dmm1vgs;
3191  GDDSSUB = gdds*mm1 + cdraindrift*dmm1vds;
3192  GBSSUB = gmbs*mm1 + cdraindrift*dmm1vbs;
3193 
3194 // PMC ISUBMOD not fully implemented
3195 // - probably we don't have the right model parameters
3196 // cdrain = cdraindrift + ISUB;
3197 
3198  cdrain = cdraindrift;
3199  cd = mode*cdrain - cbd;
3200 
3201  // diode #1 values
3202  D1csat = D1DIOtSatCur * D1DIOarea;
3203  D1csatr = D1DIOtSatRCur * D1DIOarea;
3205  D1vtr = model_.D1DIOnr * D1vt;
3206  D1evd = D1arg = D1evrev = 0.0;
3207 
3208  // compute diode dc current and derivatives
3209  if(D1csatr != 0)
3210  {
3211  D1evr = exp(D1vd/D1vtr);
3212  D1temp = 1 - D1vd/D1DIOtJctPot;
3213  D1arg = D1temp*D1temp;
3214  D1power = pow(D1arg + 0.001,model_.D1DIOgradingCoeff/2);
3215  D1isr = D1csatr*D1power;
3216  D1cdr = D1isr*(D1evr - 1);
3217  D1gdr = -D1cdr*model_.D1DIOgradingCoeff*D1temp/
3218  (D1DIOtJctPot*(D1arg + 0.001)) + D1isr*D1evr/D1vtr;
3219  } else
3220  {
3221  D1cdr = 0;
3222  D1gdr = 0;
3223  }
3224 
3225  // add in diode bulk resistance
3226  if(D1vd >= -3*D1vte)
3227  {
3228  D1evd = exp(D1vd/D1vte);
3229  D1cd = D1csat*(D1evd-1)+D1cdr + getDeviceOptions().gmin*D1vd;
3230  D1gd = D1csat*D1evd/D1vte+D1gdr + getDeviceOptions().gmin;
3231  if(model_.D1DIOikf > 0) {
3232  D1arg = sqrt(model_.D1DIOikf/(model_.D1DIOikf + D1cd));
3233  D1gd = D1arg*D1gd*(1 - 0.5*D1cd/(model_.D1DIOikf + D1cd));
3234  D1cd = D1arg*D1cd;
3235  }
3236  }
3237  else if( (!(D1DIOtBrkdwnV)) || (D1vd >= -D1DIOtBrkdwnV) )
3238  {
3239  D1arg = 3*D1vte/(D1vd*CONSTe);
3240  D1arg = D1arg * D1arg * D1arg;
3241  D1cd = -D1csat*(1 + D1arg) + getDeviceOptions().gmin*D1vd + D1cdr;
3242  D1gd = D1csat*3*D1arg/D1vd + getDeviceOptions().gmin + D1gdr;
3243  }
3244  else
3245  {
3246  D1evrev= exp(-(D1DIOtBrkdwnV + D1vd)/D1vte);
3247  D1cd = -D1csat*D1evrev + getDeviceOptions().gmin*D1vd + D1cdr;
3248  D1gd = D1csat*D1evrev/D1vte + getDeviceOptions().gmin + D1gdr;
3249  }
3250  D1cdeq = D1cd;
3251 
3252  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3253  {
3254  Xyce::dout() << " " << std::endl;
3255  Xyce::dout() << " dtype = " << model_.dtype << std::endl;
3256  Xyce::dout() << " D1vd = " << D1vd << std::endl;
3257  Xyce::dout() << " D1vtr = " << D1vtr << std::endl;
3258  Xyce::dout() << " D1csat = " << D1csat << std::endl;
3259  Xyce::dout() << " D1csatr = " << D1csatr << std::endl;
3260  Xyce::dout() << " D1gspr = " << D1gspr << std::endl;
3261  Xyce::dout() << " D1cdr = " << D1cdr << std::endl;
3262  Xyce::dout() << " D1gdr = " << D1gdr << std::endl;
3263  Xyce::dout() << " D1vte = " << D1vte << std::endl;
3264  Xyce::dout() << " D1evd = " << D1evd << std::endl;
3265  Xyce::dout() << " D1arg = " << D1arg << std::endl;
3266  Xyce::dout() << " D1evrev = " << D1evrev << std::endl;
3267  Xyce::dout() << " D1DIOtBrkdwnV = " << D1DIOtBrkdwnV << std::endl << std::endl;
3268  Xyce::dout() << " D1cd = " << D1cd << std::endl;
3269  Xyce::dout() << " D1gd = " << D1gd << std::endl;
3270  Xyce::dout() << " D1cdeq = " << D1cdeq << std::endl;
3271  }
3272 
3273  // in 3f5 this is all in a block conditioned on CKTmode, but since
3274  // it's valid for MODETRAN and MODETRANOP we'll just always do it
3275 
3276  ////
3277  // * now we do the hard part of the bulk-drain and bulk-source
3278  // * diode - we evaluate the non-linear capacitance and
3279  // * charge
3280  // *
3281  // * the basic equations are not hard, but the implementation
3282  // * is somewhat long in an attempt to avoid log/exponential
3283  // * evaluations
3284  ////
3285  ////
3286  // * charge storage elements
3287  // *
3288  // *.. bulk-drain and bulk-source depletion capacitances
3289  ////
3290  // I took out all the CAPBYPASS stuff, and the
3291  // unnecessary curly braces that wind up there if you do
3292 
3293  // can't bypass the diode capacitance calculations
3294  if(Cbs != 0 || Cbssw != 0 )
3295  {
3296  if (vbs < tDepCap)
3297  {
3298  arg=1-vbs/tBulkPot;
3299  ////
3300  // * the following block looks somewhat long and messy,
3301  // * but since most users use the default grading
3302  // * coefficients of .5, and sqrt is MUCH faster than an
3303  // * exp(log()) we use this special case code to buy time.
3304  // * (as much as 10% of total job time!)
3305  ////
3307  {
3308  if(model_.bulkJctBotGradingCoeff == .5)
3309  {
3310  sarg = sargsw = 1/sqrt(arg);
3311  }
3312  else
3313  {
3314  sarg = sargsw = exp(-model_.bulkJctBotGradingCoeff*log(arg));
3315  }
3316  }
3317  else
3318  {
3319  if(model_.bulkJctBotGradingCoeff == .5)
3320  {
3321  sarg = 1/sqrt(arg);
3322  }
3323  else
3324  {
3325  sarg = exp(-model_.bulkJctBotGradingCoeff*log(arg));
3326  }
3328  {
3329  sargsw = 1/sqrt(arg);
3330  }
3331  else
3332  {
3333  sargsw =exp(-model_.bulkJctSideGradingCoeff* log(arg));
3334  }
3335  }
3336  qbs = tBulkPot*(Cbs*(1-arg*sarg)/(1-model_.bulkJctBotGradingCoeff)
3337  +Cbssw*(1-arg*sargsw)/(1-model_.bulkJctSideGradingCoeff));
3338  capbs=Cbs*sarg + Cbssw*sargsw;
3339  }
3340  else
3341  {
3342  qbs = f4s + vbs*(f2s+vbs*(f3s/2));
3343  capbs=f2s+f3s*vbs;
3344  }
3345  }
3346  else
3347  {
3348  qbs = 0;
3349  capbs=0;
3350  }
3351 
3352  //// can't bypass the diode capacitance calculations
3353  if(Cbd != 0 || Cbdsw != 0 )
3354  {
3355 
3356  if (vbdd < tDepCap)
3357  {
3358  arg=1-vbdd/tBulkPot;
3359  ////
3360  // * the following block looks somewhat long and messy,
3361  // * but since most users use the default grading
3362  // * coefficients of .5, and sqrt is MUCH faster than an
3363  // * exp(log()) we use this special case code to buy time.
3364  // * (as much as 10% of total job time!)
3365  ////
3366  if(model_.bulkJctBotGradingCoeff == .5 &&
3368  {
3369  sarg = sargsw = 1/sqrt(arg);
3370  }
3371  else
3372  {
3373  if(model_.bulkJctBotGradingCoeff == .5)
3374  {
3375  sarg = 1/sqrt(arg);
3376  }
3377  else
3378  {
3379  sarg = exp(-model_.bulkJctBotGradingCoeff*log(arg));
3380  }
3382  {
3383  sargsw = 1/sqrt(arg);
3384  }
3385  else
3386  {
3387  sargsw =exp(-model_.bulkJctSideGradingCoeff*log(arg));
3388  }
3389  }
3390  qbd =
3391  tBulkPot*(Cbd*(1-arg*sarg)
3393  +Cbdsw*(1-arg*sargsw)
3395  capbd=Cbd*sarg + Cbdsw*sargsw;
3396  }
3397  else
3398  {
3399  qbd = f4d + vbdd * (f2d + vbdd * f3d/2);
3400  capbd= f2d + vbdd * f3d;
3401  }
3402  }
3403  else
3404  {
3405  qbd = 0;
3406  capbd = 0;
3407  }
3408 
3409  // charge storage elements
3410  D1czero = D1DIOtJctCap*D1DIOarea;
3411  if (D1vd < D1DIOtDepCap)
3412  {
3413  D1arg = 1-D1vd/model_.D1DIOjunctionPot;
3414  D1sarg = exp(-model_.D1DIOgradingCoeff*log(D1arg));
3416  model_.D1DIOjunctionPot*D1czero*
3417  (1-D1arg*D1sarg)/(1-model_.D1DIOgradingCoeff);
3418  D1capd = model_.D1DIOtransitTime*D1gd + D1czero*D1sarg;
3419  }
3420  else
3421  {
3422  D1czof2 = D1czero/model_.D1DIOf2;
3423  D1DIOcapCharge = model_.D1DIOtransitTime*D1cd+D1czero*
3424  D1DIOtF1+D1czof2*(model_.D1DIOf3*
3429  D1capd = model_.D1DIOtransitTime*D1gd+D1czof2*
3432  }
3433 
3434  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3435  {
3436  Xyce::dout() << " " << std::endl;
3437  Xyce::dout() << " Going into qmeyer..." << std::endl;
3438  Xyce::dout() << " Mode is " << mode << std::endl;
3439  Xyce::dout() << " Args are vgps = " << vgps << " vgpdd = " << vgpdd << std::endl;
3440  Xyce::dout() << " Vgpb = " << Vgpb << " Von = " << Von << " Vddsat = " << Vddsat << std::endl;
3441  Xyce::dout() << " tPhi = " << tPhi << " OxideCap = " << OxideCap << std::endl;
3442  }
3443  if (mode > 0)
3444  {
3445  if(model_.cve == 1)
3446  {
3447  UCCMqmeyer (vgps,vgpdd,Vgpb,Von,Vddsat,
3449  }
3450  else
3451  {
3452  // Ward-like model
3454  }
3455  }
3456  else
3457  {
3458  if(model_.cve == 1)
3459  {
3460  UCCMqmeyer (vgpdd,vgps,Vgpb,Von,Vddsat,
3462  }
3463  else
3464  {
3465  // Ward-like model
3467  }
3468  }
3469 
3470  if((getSolverState().dcopFlag))
3471  {
3472  Capgs = 2 * capgs + GateSourceOverlapCap ;
3474  Capgb = 2 * capgb + GateBulkOverlapCap ;
3475  }
3476  else
3477  {
3478  capgs_old = (*extData.currStaVectorPtr)[li_state_capgs];
3479  capgdd_old = (*extData.currStaVectorPtr)[li_state_capgdd];
3480  capgb_old = (*extData.currStaVectorPtr)[li_state_capgb];
3481 
3482  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3483  {
3484  Xyce::dout() << "Doing meyer back averaging..."<< std::endl;
3485  Xyce::dout() << " capgs = " << capgs << " capgs_old = " << capgs_old << std::endl;
3486  Xyce::dout() << " capgdd = " << capgdd << " capgdd_old = " << capgdd_old << std::endl;
3487  Xyce::dout() << " capgb = " << capgb << " capgb_old = " << capgb_old << std::endl;
3488  }
3489  Capgs = ( capgs+capgs_old + GateSourceOverlapCap );
3490  Capgdd = ( capgdd+capgdd_old + GateDrainOverlapCap );
3491  Capgb = ( capgb+capgb_old + GateBulkOverlapCap );
3492  }
3493  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3494  {
3495  Xyce::dout() << "Capgs = " << Capgs << std::endl;
3496  Xyce::dout() << "Capgdd = " << Capgdd << std::endl;
3497  Xyce::dout() << "Capgb = " << Capgb << std::endl;
3498  Xyce::dout() << "capgs = " << capgs << std::endl;
3499  Xyce::dout() << "capgdd = " << capgdd << std::endl;
3500  Xyce::dout() << "capgb = " << capgb << std::endl;
3501  }
3502  Capgs *= (Capgs < 0.0)?-1.0:1.0;
3503  Capgdd *= (Capgdd < 0.0)?-1.0:1.0;
3504  Capgb *= (Capgb < 0.0)?-1.0:1.0;
3505 
3506  // Voltage-dependent drain-drift conductance
3507  double absV = fabs(Vddd);
3509  if(gddd != 0) gddd = 1.0/gddd;
3510  draindriftCond = gddd;
3511 
3512  // need a more precise derivative than gddd for the jacobian.
3513  // w.r.t. Vd and Vdd.
3514  double dVddd_dVd (1.0);
3515  double d_absVddd_dVd (0.0);
3516  if (Vddd > 0.0)
3517  {
3518  d_absVddd_dVd = 1.0;
3519  }
3520  else if (Vddd < 0.0)
3521  {
3522  d_absVddd_dVd = -1.0;
3523  }
3524 
3525  dIdd_dVd = dVddd_dVd * gddd -
3526  Vddd * (gddd*gddd)*(model_.driftParamB)*d_absVddd_dVd;
3527 
3528  // now do the shundt resistor
3529  if(model_.rdsshunt == 0.0)
3530  {
3531  gdsshunt = 0.0;
3532  }
3533  else
3534  {
3536  }
3537 
3538  Idrain = drainCond * Vdddp;
3539  Igate = gateCond * Vggp;
3540  Isource = sourceCond * Vssp;
3542  Irdsshunt = gdsshunt * (Vd - Vs);
3543  // Rd1rs current:
3544  Ird1rs = D1gspr*(Vs - Vd1p);
3545 
3546  if (mode >= 0) // Normal mode
3547  {
3548  Gm = gm; // (xnrm-xrev)*gm in 3f5
3549  Gmbs = gmbs; // (xnrm-xrev)*gmbs in 3f5
3550  nrmsum = Gm+Gmbs; // xnrm*(gm+gmbs)
3551  revsum = 0; // xrev*(gm+gmbs)
3553  }
3554  else
3555  {
3556  Gm = -gm;
3557  Gmbs = -gmbs;
3558  nrmsum = 0;
3559  revsum = -(Gm+Gmbs); // because Gm and Gmbs already have - in them!
3560  cdreq = -(model_.dtype)*cdrain;
3561  }
3562 
3563  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3564  {
3565  Xyce::dout() << " Done with Instance::updateIntermediateVars " << std::endl;
3566  Xyce::dout() << " mode = " << mode << std::endl;
3567  Xyce::dout() << " Idrain = " << Idrain << std::endl;
3568  Xyce::dout() << " Igate = " << Igate << std::endl;
3569  Xyce::dout() << " Isource = " << Isource << std::endl;
3570  Xyce::dout() << " Idraindrift = " << Idraindrift << std::endl;
3571  Xyce::dout() << " Ird1rs = " << Ird1rs << std::endl;
3572  Xyce::dout() << " dIdd_dVd = " << dIdd_dVd << std::endl;
3573  Xyce::dout() << " gddd = " << gddd << std::endl;
3574  Xyce::dout() << " Irdsshunt = " << Irdsshunt << std::endl;
3575  Xyce::dout() << " D1DIOcapCurrent = " << D1DIOcapCurrent << std::endl;
3576  Xyce::dout() << " cbd = " << cbd << std::endl;
3577  Xyce::dout() << " cbs = " << cbs << std::endl;
3578  Xyce::dout() << " qbd = " << qbd << std::endl;
3579  Xyce::dout() << " qbs = " << qbs << std::endl;
3580  Xyce::dout() << " cdrain = " << cdrain << std::endl;
3581  Xyce::dout() << " cdraindrift = " << cdraindrift << std::endl;
3582  Xyce::dout() << " cdreq = " << cdreq << std::endl;
3583  Xyce::dout() << " gdds = " << gdds << std::endl;
3584  Xyce::dout() << " gdsshunt = " << gdsshunt << std::endl;
3585  Xyce::dout() << " gm = " << gm << std::endl;
3586  Xyce::dout() << " gmbs = " << gmbs << std::endl;
3587  Xyce::dout() << " Gm = " << Gm << std::endl;
3588  Xyce::dout() << " Gmbs = " << Gmbs << std::endl;
3589  }
3590 
3591  /// CURRENTS to load into RHS:
3592 
3593  // current out of drain is
3594  // Idraindrift + Irdsshunt - Dtype*(D1cdeq + D1DIOcapCurrent)
3595 
3596  // current out of gate:
3597  // dtype*( (deriv of qgs) + (deriv of qgdd) + (deriv of qgb))
3598 
3599  // the current *out of* the source is
3600  // Isource - Irdsshunt + Ird1rs
3601 
3602  // current out of bulk is
3603  // dtype*(deriv of qbd) + dtype*cbd + dtype*cbs + dtype*(deriv of qbs)
3604  // - dtype*(deriv of qgb)
3605 
3606  // current out of drain' is
3607  // -Idrain - dtype*(deriv of qgd) - (deriv of qbd) - dtype*cbd +
3608  // mode*dtype*cdrain
3609 
3610  // the current out of the source' is
3611  // -Isource - dtype*(deriv of qgs) - dtype*cbs - (deriv of qbs) -
3612  // mode*dtype*cdrain -Irdsshunt
3613 
3614  // the current out of the draindrift node is
3615  // Idrain - IdraindriftDtype*(D1cdeq + D1DIOcapCurrent)
3616 
3617  // the current out of the d1' pos node is
3618  // Dtype*(D1cdeq + D1DIOcapCurrent) - Ird1rs
3619 
3620  return bsuccess;
3621 }
3622 
3623 //-----------------------------------------------------------------------------
3624 // Function : Instance::updateTemperature
3625 // Purpose :
3626 // Special Notes :
3627 // Scope : public
3628 // Creator : pmc
3629 // Creation Date : 2/17/2004
3630 //-----------------------------------------------------------------------------
3631 bool Instance::updateTemperature ( const double & temp_tmp)
3632 {
3633  bool bsuccess = true;
3634  // mos3temp vars
3635  double czbd(0.0); // zero voltage bulk-drain capacitance
3636  double czbdsw(0.0); // zero voltage bulk-drain sidewall capacitance
3637  double czbs(0.0); // zero voltage bulk-source capacitance
3638  double czbssw(0.0); // zero voltage bulk-source sidewall capacitance
3639  double arg(0.0); // 1 - fc
3640  double sarg(0.0); // (1-fc) ^^ (-mj)
3641  double sargsw(0.0); // (1-fc) ^^ (-mjsw)
3642  double ratio,ratio4(0.0);
3643  double fact2(0.0);
3644  double kt(0.0);
3645  double egfet(0.0);
3646  double pbfact(0.0);
3647  double capfact(0.0);
3648  double phio(0.0);
3649  double pbo(0.0);
3650  double gmanew,gmaold(0.0);
3651  // end of mos3temp stuff
3652 
3653 //variables related to the diode 1
3654  double D1xfc(0.0);
3655  double D1vte_loc(0.0);
3656  double D1cbv(0.0);
3657  double D1xbv(0.0);
3658  double D1xcbv(0.0);
3659  double D1tol(0.0);
3660  double D1vt_loc(0.0);
3661  register int D1iter(0.0);
3662  double D1egfet1(0.0),D1arg1(0.0),D1fact1(0.0),D1pbfact1(0.0),D1pbo(0.0),D1gmaold(0.0);
3663  double D1fact2(0.0),D1pbfact(0.0),D1arg(0.0),D1egfet(0.0),D1gmanew(0.0);
3664  double reltol(0.0);
3665 
3666  double tnom(0.0);
3667  double VMAX(0.0);
3668 
3669 // oxide dielectric permitivity
3670 //#define EPSSIO2 3.453130e-11 == CONSTEPSOX (Xyce value)
3671 
3672  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3673  {
3674  Xyce::dout() << subsection_divider << std::endl;
3675  Xyce::dout() << " Instance::Begin of updateTemperature. \n";
3676  Xyce::dout() <<" name = " << getName() << std::endl;
3677  Xyce::dout() << std::endl;
3678  }
3679 
3680  // first set the instance temperature to the new temperature:
3681  if (temp_tmp != -999.0) temp = temp_tmp;
3683  {
3684  // make sure interpolation doesn't take any resistance negative
3689  if(model_.D1DIOresist < 0) model_.D1DIOresist = 0;
3690 
3691  // some params may have changed during interpolation
3693  }
3694 
3695  VMAX = model_.maxDriftVel;
3696  tnom = model_.tnom;
3697  ratio = temp/tnom;
3698  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3699  {
3700  Xyce::dout() << " Temperature = "<< temp << std::endl;
3701  Xyce::dout() << " tnom = " << tnom << std::endl;
3702  Xyce::dout() << " ratio = " << ratio << std::endl;
3703  }
3704 
3705  vt = temp * CONSTKoverQ;
3706  ratio = temp/tnom;
3707  fact2 = temp/CONSTREFTEMP;
3708  kt = temp * CONSTboltz;
3709  egfet = 1.16-(7.02e-4*temp*temp)/(temp+1108);
3710  arg = -egfet/(kt+kt)+1.1150877/(CONSTboltz*(CONSTREFTEMP+CONSTREFTEMP));
3711  pbfact = -2*vt *(1.5*log(fact2)+CONSTQ*arg);
3712 
3713  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3714  {
3715  Xyce::dout() << " fact1 = " << model_.fact1 << std::endl;
3716  Xyce::dout() << " vtnom = " << model_.vtnom << std::endl;
3717  Xyce::dout() << " egfet1 = " << model_.egfet1 << std::endl;
3718  Xyce::dout() << " pbfact1= " << model_.pbfact1 << std::endl;
3719  Xyce::dout() << " vt = " << vt << std::endl;
3720  Xyce::dout() << " ratio = " << ratio << std::endl;
3721  Xyce::dout() << " fact2 = " << fact2 << std::endl;
3722  Xyce::dout() << " kt = " << kt << std::endl;
3723  Xyce::dout() << " egfet = " << egfet << std::endl;
3724  Xyce::dout() << " arg = " << arg << std::endl;
3725  Xyce::dout() << " pbfact = " << pbfact << std::endl;
3726  }
3727 
3728  ratio4 = ratio * sqrt(ratio);
3729  tSurfMob = 1.e-4*model_.surfaceMobility/ratio4;
3730  phio= (model_.phi-model_.pbfact1)/model_.fact1;
3731  tPhi = fact2 * phio + pbfact;
3732  tVbi = model_.vt0 - model_.dtype *
3733  (model_.gamma* sqrt(model_.phi))+.5*(model_.egfet1-egfet)
3734  + model_.dtype*.5* (tPhi-model_.phi);
3735  tVto = tVbi + model_.dtype * model_.gamma * sqrt(tPhi);
3736  if(model_.mdtemp != 0) tVto = model_.vt0 -
3737  model_.kvt*(temp-tnom);
3739  exp(-egfet/vt+model_.egfet1/model_.vtnom);
3741  exp(-egfet/vt+model_.egfet1/model_.vtnom);
3743  gmaold = (model_.bulkJctPotential-pbo)/pbo;
3744  capfact = 1/(1+model_.bulkJctBotGradingCoeff*
3745  (4e-4*(model_.tnom-CONSTREFTEMP)-gmaold));
3746  tCbd = model_.capBD * capfact;
3747  tCbs = model_.capBS * capfact;
3748  tCj = model_.bulkCapFactor * capfact;
3749  capfact = 1/(1+model_.bulkJctSideGradingCoeff*
3750  (4e-4*(model_.tnom-CONSTREFTEMP)-gmaold));
3751  tCjsw = model_.sideWallCapFactor * capfact;
3752  tBulkPot = fact2 * pbo+pbfact;
3753  gmanew = (tBulkPot-pbo)/pbo;
3754  capfact = (1+model_.bulkJctBotGradingCoeff *
3755  (4e-4*(temp-CONSTREFTEMP)-gmanew));
3756  tCbd *= capfact;
3757  tCbs *= capfact;
3758  tCj *= capfact;
3759  capfact = (1+model_.bulkJctSideGradingCoeff *
3760  (4e-4*(temp-CONSTREFTEMP)-gmanew));
3761  tCjsw *= capfact;
3763 
3764  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
3765  {
3766  Xyce::dout() << " ratio4 = " << ratio4 << std::endl;
3767  Xyce::dout() << " tSurfMob = " << tSurfMob << std::endl;
3768  Xyce::dout() << " phio = " << phio << std::endl;
3769  Xyce::dout() << " tPhi = " << tPhi << std::endl;
3770  Xyce::dout() << " tVbi = " << tVbi << std::endl;
3771  Xyce::dout() << " tVto = " << tVto << std::endl;
3772  Xyce::dout() << " tSatCur = " << tSatCur << std::endl;
3773  Xyce::dout() << " tSatCurDens = " << tSatCurDens << std::endl;
3774  Xyce::dout() << " pbo = " << pbo << std::endl;
3775  Xyce::dout() << " gmaold = " << gmaold << std::endl;
3776  Xyce::dout() << " tBulkPot = " << tBulkPot << std::endl;
3777  Xyce::dout() << " gmanew = " << gmanew << std::endl;
3778  Xyce::dout() << " capfact = " << capfact << std::endl;
3779  Xyce::dout() << " tCbd = " << tCbd << std::endl;
3780  Xyce::dout() << " tCbs = " << tCbs << std::endl;
3781  Xyce::dout() << " tCj = " << tCj << std::endl;
3782  Xyce::dout() << " capfact = " << capfact << std::endl;
3783  Xyce::dout() << " tCjsw = " << tCjsw << std::endl;
3784  Xyce::dout() << " tDepCap = " << tDepCap << std::endl;
3785  }
3786 
3787  if( (model_.jctSatCurDensity == 0) || (drainArea == 0) ||
3788  (sourceArea == 0) )
3789  {
3791  vt*log(vt/(CONSTroot2*model_.jctSatCur));
3792  }
3793  else
3794  {
3795  drainVcrit = vt * log( vt / (CONSTroot2 *
3797  sourceVcrit = vt * log( vt / (CONSTroot2 *
3799  }
3800  if(model_.capBDGiven)
3801  {
3802  czbd = tCbd;
3803  }
3804  else
3805  {
3807  {
3808  czbd=tCj*drainArea;
3809  }
3810  else
3811  {
3812  czbd=0;
3813  }
3814  }
3816  {
3817  czbdsw= tCjsw * drainPerimeter;
3818  }
3819  else
3820  {
3821  czbdsw=0;
3822  }
3823  arg = 1-model_.fwdCapDepCoeff;
3824  sarg = exp( (-model_.bulkJctBotGradingCoeff) * log(arg) );
3825  sargsw = exp( (-model_.bulkJctSideGradingCoeff) * log(arg) );
3826  Cbd = czbd;
3827  Cbdsw = czbdsw;
3828  f2d = czbd*(1-model_.fwdCapDepCoeff*
3829  (1+model_.bulkJctBotGradingCoeff))* sarg/arg
3830  + czbdsw*(1-model_.fwdCapDepCoeff*
3831  (1+model_.bulkJctSideGradingCoeff))*sargsw/arg;
3832  f3d = czbd * model_.bulkJctBotGradingCoeff * sarg/arg/tBulkPot
3833  + czbdsw * model_.bulkJctSideGradingCoeff *
3834  sargsw/arg/tBulkPot;
3835  f4d = czbd*tBulkPot*(1-arg*sarg)/(1-model_.bulkJctBotGradingCoeff) +
3836  czbdsw*tBulkPot*(1-arg*sargsw)/
3838  f3d/2*(tDepCap*tDepCap) - tDepCap * f2d;
3839  if(model_.capBSGiven)
3840  {
3841  czbs=tCbs;
3842  }
3843  else
3844  {
3846  {
3847  czbs=tCj*sourceArea;
3848  }
3849  else
3850  {
3851  czbs=0;
3852  }
3853  }
3855  {
3856  czbssw = tCjsw * sourcePerimeter;
3857  }
3858  else
3859  {
3860  czbssw=0;
3861  }
3862  arg = 1-model_.fwdCapDepCoeff;
3863  sarg = exp( (-model_.bulkJctBotGradingCoeff) * log(arg) );
3864  sargsw = exp( (-model_.bulkJctSideGradingCoeff) * log(arg) );
3865  Cbs = czbs;
3866  Cbssw = czbssw;
3867  f2s = czbs*(1-model_.fwdCapDepCoeff*
3868  (1+model_.bulkJctBotGradingCoeff))*sarg/arg +
3869  czbssw*(1-model_.fwdCapDepCoeff*
3870  (1+model_.bulkJctSideGradingCoeff))*sargsw/arg;
3871  f3s = czbs * model_.bulkJctBotGradingCoeff * sarg/arg/tBulkPot
3872  + czbssw * model_.bulkJctSideGradingCoeff *
3873  sargsw/arg /tBulkPot;
3874  f4s = czbs*tBulkPot*(1-arg*sarg)/(1-model_.bulkJctBotGradingCoeff)
3875  + czbssw*tBulkPot*(1-arg*sargsw)/
3877  (tDepCap*tDepCap) - tDepCap * f2s;
3878 
3880 
3881  if(!model_.vpGiven)
3882  {
3883  vp = VMAX*(l-2*model_.latDiff)/tSurfMob;
3884  if(model_.dtype == CONSTNMOS) vp *= 2.0;
3885  } else
3886  vp = model_.vp;
3887 
3888  gchi0 = CONSTQ*w/(l-2*model_.latDiff);
3889  gammas =
3891  model_.wgammas*(1-model_.w0/w);
3892  gammal =
3894  model_.wgammal*(1-model_.w0/w);
3895  if(gammal != 0.0) {
3896  vthLimit = gammas/(2*gammal);
3898  } else
3899  vthLimit = Util::MachineDependentParams::DoubleMax();
3900  vtoo = tVto*model_.dtype+ gammal*tPhi- gammas*sqrt(tPhi);
3901 
3902  // quantities related to diode #1
3903 
3904  D1DIOtemp = temp;
3905  model_.D1DIOnomTemp = tnom;
3906  D1vt_loc = CONSTKoverQ * D1DIOtemp;
3907  D1xfc = log(1-model_.D1DIOdepletionCapCoeff);
3908 
3909 // this part gets really ugly - I don't know how to explain these equations
3910 
3911  D1fact2 = D1DIOtemp/CONSTREFTEMP;
3912  D1egfet = 1.16-(7.02e-4*D1DIOtemp*D1DIOtemp)/(D1DIOtemp+1108);
3913  D1arg = -D1egfet/(2*CONSTboltz*D1DIOtemp) + 1.1150877/
3914  (CONSTboltz*(CONSTREFTEMP+CONSTREFTEMP));
3915  D1pbfact = -2*D1vt_loc*(1.5*log(D1fact2)+CONSTQ*D1arg);
3916  D1egfet1 = 1.16 - (7.02e-4*model_.D1DIOnomTemp*
3918  D1arg1 = -D1egfet1/(CONSTboltz*2*model_.D1DIOnomTemp) +
3919  1.1150877/(2*CONSTboltz*CONSTREFTEMP);
3920  D1fact1 = model_.D1DIOnomTemp/CONSTREFTEMP;
3921  D1pbfact1 = -2*CONSTKoverQ*model_.D1DIOnomTemp*
3922  (1.5*log(D1fact1)+CONSTQ*D1arg1);
3923  D1pbo = (model_.D1DIOjunctionPot-D1pbfact1)/D1fact1;
3924  D1gmaold = (model_.D1DIOjunctionPot -D1pbo)/D1pbo;
3927  (400e-6*(model_.D1DIOnomTemp-CONSTREFTEMP)-D1gmaold) );
3928  D1DIOtJctPot = D1pbfact+D1fact2*D1pbo;
3929  D1gmanew = (D1DIOtJctPot-D1pbo)/D1pbo;
3931  (400e-6*(D1DIOtemp-CONSTREFTEMP)-D1gmanew);
3932  D1DIOtSatCur = model_.D1DIOsatCur*exp( ((D1DIOtemp/
3933  model_.D1DIOnomTemp)-1)*
3935  (model_.D1DIOemissionCoeff*D1vt_loc)+
3938  log(D1DIOtemp/model_.D1DIOnomTemp) );
3940 
3941  // the defintion of f1, just recompute after temperature adjusting
3942  // all the variables used in it
3943  D1DIOtF1=D1DIOtJctPot*(1-exp((1-model_.D1DIOgradingCoeff)*D1xfc))/
3945 
3946  // same for Depletion Capacitance
3948 
3949  // and Vcrit
3950  D1vte_loc=model_.D1DIOemissionCoeff*D1vt_loc;
3951  D1DIOtVcrit=D1vte_loc*log(D1vte_loc/(CONSTroot2*D1DIOtSatCur));
3952 
3953  // and now to copute the breakdown voltage, again using
3954  // temperature adjusted basic parameters
3956  {
3958  if (D1cbv < D1DIOtSatCur*model_.D1DIObreakdownVoltage/D1vt_loc)
3959  {
3960  D1cbv=D1DIOtSatCur*model_.D1DIObreakdownVoltage/D1vt_loc;
3961  Xyce::dout() << " breakdown current increased to " << D1cbv <<
3962  "to resolve incompatability " <<
3963  "with specified saturation current" << std::endl;
3965  } else
3966  {
3967  reltol = 1e-3;
3968  D1tol=reltol*D1cbv;
3970  D1vt_loc*log(1+D1cbv/D1DIOtSatCur);
3971  for(D1iter=0; D1iter<25; ++D1iter)
3972  {
3973  D1xbv=model_.D1DIObreakdownVoltage-D1vt_loc*log(D1cbv/
3974  D1DIOtSatCur+1-D1xbv/D1vt_loc);
3975  D1xcbv=D1DIOtSatCur*(exp((model_.D1DIObreakdownVoltage
3976  -D1xbv)/D1vt_loc)-1+D1xbv/D1vt_loc);
3977  if (fabs(D1xcbv-D1cbv) <= D1tol) goto matched;
3978  }
3979  Xyce::dout() << " unable to match forward and reverse diode regions: D1bv = "
3980  << D1xbv << " D1ibv = " << D1xcbv << std::endl;
3981  }
3982  matched:
3983  D1DIOtBrkdwnV = D1xbv;
3984  }
3985 
3986  return bsuccess;
3987 }
3988 
3989 //-----------------------------------------------------------------------------
3990 // Function : Instance::updatePrimaryState
3991 // Purpose :
3992 // Special Notes :
3993 // Scope : public
3994 // Creator : Tom Russo, Component Information and Models
3995 // Creation Date : 02/28/01
3996 //-----------------------------------------------------------------------------
3998 {
3999  bool bsuccess = true;
4000  double vgs1(0.0), vgdd1(0.0), vbs1(0.0),vgb1(0.0), vdds1(0.0);
4001  double * staVector = extData.nextStaVectorRawPtr;
4002  double * currStaVector = extData.currStaVectorRawPtr;
4003 
4004  bool tmpBool = updateIntermediateVars ();
4005  bsuccess = bsuccess && tmpBool;
4006 
4007  // voltage drops:
4008  staVector[li_state_vbdd] = vbdd;
4009  staVector[li_state_vbs] = vbs;
4010  staVector[li_state_vgps] = vgps;
4011  staVector[li_state_vdds] = vdds;
4012  staVector[li_state_D1vd] = D1vd;
4013 
4014  // now the meyer capacitances
4015  // we didn't calculate these charges in update IntermediateVars
4016  // but we did calculate the voltage drops and capacitances.
4017  // first store the capacitances themselves:
4018  staVector[li_state_capgs] = capgs;
4019  staVector[li_state_capgdd] = capgdd;
4020  staVector[li_state_capgb] = capgb;
4021 
4022  // now the charges
4023  // BE CAREFUL! We can only do Q=CV for DCOP! Otherwise it's
4024  // supposed to be *INTEGRATED*:
4025  // Q = int(t0,t1)C(V)*dV --- and we approximate that by
4026  // Q(t1)-Q(t0) = CBar*(V(t1)-V(t0)) where CBar is the average.
4027  // Now with Meyer back averaging, Capxx is the average between the last
4028  // time step and this one. So we gotta do the right thing for non-DCOP
4029  // when backaverage is on.
4030 
4031  if((getSolverState().dcopFlag))
4032  {
4033  qgs = Capgs *vgps;
4034  qgdd = Capgdd*vgpdd;
4035  qgb = Capgb *Vgpb;
4036  }
4037  else
4038  {
4039  // get the ones from last time step
4040  qgs = currStaVector[li_state_qgs];
4041  qgdd = currStaVector[li_state_qgdd];
4042  qgb = currStaVector[li_state_qgb];
4043  // get the voltage drops, too
4044  vgs1 = currStaVector[li_state_vgps];
4045  vbs1 = currStaVector[li_state_vbs];
4046  vdds1 = currStaVector[li_state_vdds];
4047 
4048  vgb1 = vgs1-vbs1;
4049  vgdd1 = vgs1-vdds1;
4050 
4051  // NOW we can calculate the charge update
4052  qgs += Capgs*(vgps-vgs1);
4053  qgdd += Capgdd*(vgpdd-vgdd1);
4054  qgb += Capgb*((vgps-vbs)-vgb1);
4055  }
4056 
4057  staVector[li_state_qgs] = qgs;
4058  staVector[li_state_qgdd] = qgdd;
4059  staVector[li_state_qgb] = qgb;
4060 
4061  // and the diode parasitic capacitors
4062  // these charges were set in updateIntermediateVars
4063  staVector[li_state_qbd] = qbd;
4064  staVector[li_state_qbs] = qbs;
4065 
4066 // diode #1 quatntities
4067 
4069 
4070  // In the case of a charge, need to handle this:
4071  // Ensure dQ/dt = 0 for initial step after dcOP by setting the history equal
4072  // to the current value
4073  if( !getSolverState().dcopFlag && getSolverState().initTranFlag && !getSolverState().newtonIter )
4074  {
4075  currStaVector[li_state_D1DIOcapCharge] = D1DIOcapCharge;
4076  }
4077 
4078  staVector[ li_state_von ] = von;
4079 
4080  return bsuccess;
4081 }
4082 
4083 //-----------------------------------------------------------------------------
4084 // Function : Instance::updateSecondaryState
4085 // Purpose :
4086 // Special Notes :
4087 // Scope : public
4088 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
4089 // Creation Date : 01/09/01
4090 //-----------------------------------------------------------------------------
4092 {
4093  double * staDerivVector = (extData.nextStaDerivVectorRawPtr);
4094 
4095  cqgs = staDerivVector[li_state_qgs];
4096  cqgdd = staDerivVector[li_state_qgdd];
4097  cqgb = staDerivVector[li_state_qgb];
4098  cqbd = staDerivVector[li_state_qbd];
4099  cqbs = staDerivVector[li_state_qbs];
4100 
4101  D1DIOcapCurrent = staDerivVector[li_state_D1DIOcapCharge];
4102 
4103  double Dtype =model_.dtype;
4104  double D1current = Dtype*(D1cdeq + D1DIOcapCurrent);
4105 
4106  // Even though these do not look exactly like the RHS loads, they do account
4107  // for the currents into the leads even with the lead resistors in place
4108  // --- since the lead resistances are the only thing between the primed node
4109  // and the lead, the currents are still correct AT CONVERGED SOLUTIONS
4110 
4111  /*
4112  lead current load has been moved to F & Q loads
4113  LeadCurrentD = -model_.dtype*(-mode*cdrain+cbd+cqgdd+cqbd) - D1current-Irdsshunt;
4114  LeadCurrentG = model_.dtype*(cqgs+cqgb+cqgdd);
4115  LeadCurrentS = model_.dtype*(-mode*cdrain+cbs-cqgs-cqbs) + D1current+Irdsshunt;
4116  LeadCurrentB = model_.dtype*(cbd-cbs-cqgb+cqbd+cqbs);
4117  */
4118 
4119  return true;
4120 }
4121 
4122 
4123 //-----------------------------------------------------------------------------
4124 // Function : Instance::processParams
4125 // Purpose :
4126 // Special Notes :
4127 // Scope : public
4128 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
4129 // Creation Date : 6/03/02
4130 //-----------------------------------------------------------------------------
4132 {
4133 
4134  if (DEBUG_DEVICE && isActive(Diag::DEVICE_PARAMETERS) && getSolverState().debugTimeFlag)
4135  {
4136  Xyce::dout() << " L = " << l << std::endl;
4137  Xyce::dout() << " W = " << w<< std::endl;
4138  Xyce::dout() << " drainArea = " << drainArea<< std::endl;
4139  Xyce::dout() << " sourceArea = " << sourceArea<< std::endl;
4140  Xyce::dout() << " drainSquares = " << drainSquares<< std::endl;
4141  Xyce::dout() << " sourceSquares = " << sourceSquares<< std::endl;
4142  Xyce::dout() << " drainPerimeter = " << drainPerimeter<< std::endl;
4143  Xyce::dout() << " sourcePerimeter = " << sourcePerimeter<< std::endl;
4144  Xyce::dout() << " drainCond = " << drainCond<< std::endl;
4145  Xyce::dout() << " sourceCond = " << sourceCond << std::endl;
4146  Xyce::dout() << " draindriftCond = " << draindriftCond<< std::endl;
4147  Xyce::dout() << " temp = " << temp<< std::endl;
4148  }
4149 
4150  // now set the temperature related stuff
4152 
4153  return true;
4154 }
4155 
4156 // Additional Declarations
4157 
4158 
4159 //-----------------------------------------------------------------------------
4160 // Function : Instance::UCCMqmeyer
4161 // Purpose : Compute the MOS overlap capacitances as functions of the
4162 // device terminal voltages
4163 //
4164 // Special Notes :
4165 //
4166 // Scope : public
4167 // Creator : pmc
4168 // Creation Date : 02/17/2004
4169 //-----------------------------------------------------------------------------
4171  double xvgs, // initial voltage gate-source (mode > 0)
4172  double xvgdd, // initial voltage gate-drain (mode > 0)
4173  double vgb, // initial voltage gate-bulk
4174  double von_local, // threshold voltage
4175  double vddsat_local, // saturation drain voltage
4176  double & capgs_local, // non-constant portion of g-s overlap capacitance
4177  double & capgdd_local,// non-constant portion of g-d overlap capacitance
4178  double & capgb_local, // non-constant portion of g-b overlap capacitance
4179  double phi,
4180  double cox // oxide capactiance
4181  )
4182 {
4183  double vdds_local(0.0);
4184  double vddif(0.0);
4185  double vddif1(0.0);
4186  double vddif2(0.0);
4187  double vgst(0.0);
4188  double etavt(0.0);
4189  double cgc(0.0);
4190  double x(0.0);
4191 
4192 #define W w
4193 #define L EffectiveLength
4194 #define ETA model_.eta
4195 #define TOX model_.oxideThickness
4196 
4197  vgst = xvgs-von_local;
4198  if (vgst <= -phi)
4199  {
4200  capgb_local = cox/2;
4201  capgs_local = 0;
4202  capgdd_local = 0;
4203  }
4204  else if (vgst <= -phi/2)
4205  {
4206  capgb_local = -vgst*cox/(2*phi);
4207  capgs_local = 0;
4208  capgdd_local = 0;
4209  }
4210  else if (vgst <= 0)
4211  {
4212  capgb_local = -vgst*cox/(2*phi);
4213  capgs_local = vgst*cox/(1.5*phi)+cox/3;
4214  capgdd_local = 0;
4215  }
4216  else
4217  {
4218  vdds_local = xvgs-xvgdd;
4219  capgb_local = 0;
4220  if (vddsat_local <= vdds_local)
4221  {
4222  capgs_local = cox/3;
4223  capgdd_local = 0;
4224  }
4225  else
4226  {
4227  vddif = 2.0*vddsat_local-vdds_local;
4228  vddif1 = vddsat_local-vdds_local/*-1.0e-12*/;
4229  vddif2 = vddif*vddif;
4230  capgdd_local = cox*(1.0-vddsat_local*vddsat_local/vddif2)/3;
4231  capgs_local = cox*(1.0-vddif1*vddif1/vddif2)/3;
4232  }
4233  }
4234 
4235  // at this point we have unmodified Xyce values
4236  // the following is the UCCM adjustment
4237 
4238  if(model_.cv == 2 && vddsat_local != 0)
4239  {
4240  vdds_local = fabs(xvgs-xvgdd);
4241  vdds_local = vdds_local/pow(1+pow(vdds_local/vddsat_local,model_.mc),
4242  1.0/model_.mc);
4243  vddif = 2.0*vddsat_local-vdds_local;
4244  vddif1 = vddsat_local-vdds_local;
4245  vddif2 = vddif*vddif;
4246  etavt = ETA*vt;
4247  x = vgst/etavt;
4248  cgc = L*W/(TOX/CONSTEPSSI+etavt/(CONSTQ*n0)*exp(-x));
4249  capgdd_local = cgc*(1.0-vddsat_local*vddsat_local/vddif2)/3;
4250  capgs_local = cgc*(1.0-vddif1*vddif1/vddif2)/3;
4251  }
4252 
4253  return true;
4254 }
4255 
4256 //-----------------------------------------------------------------------------
4257 // Function : Instance::UCCMMeyercap
4258 // Purpose :
4259 //
4260 //
4261 // Special Notes :
4262 //
4263 // Scope : public
4264 // Creator : pmc
4265 // Creation Date : 02/17/2004
4266 //-----------------------------------------------------------------------------
4268  double xvgs, // initial voltage gate-source (mode > 0)
4269  double xvgdd, // initial voltage gate-drain (mode > 0)
4270  double vgb, // initial voltage gate-bulk
4271  double & cgs,
4272  double & cgd,
4273  double & cgb
4274  )
4275 {
4276  double vbs_local(0.0);
4277  double vdds_local(0.0);
4278  double zetanb(0.0);
4279  double xi(0.0);
4280  double xisqrt(0.0);
4281  double nsd(0.0);
4282  double nss(0.0);
4283  double tnss(0.0);
4284  double tnsd(0.0);
4285  double DeltaVT(0.0);
4286  double DnsdVgs(0.0);
4287  double DnsdVgd(0.0);
4288  double DnsdVgb(0.0);
4289  double DnssVgs(0.0);
4290  double DnssVgd(0.0);
4291  double DnssVgb(0.0);
4292  double DndepVgs(0.0);
4293  double DndepVgd(0.0);
4294  double DndepVgb(0.0);
4295  double UI(0.0);
4296  double VI(0.0);
4297  double DxiVgs(0.0);
4298  double DxiVgd(0.0);
4299  double DxiVgb(0.0);
4300  double DUIVgs(0.0);
4301  double DUIVgd(0.0);
4302  double DUIVgb(0.0);
4303  double DVIVgs(0.0);
4304  double DVIVgd(0.0);
4305  double DVIVgb(0.0); // With use of the chain rule, this routine may
4306  double DqiVgs(0.0); // be used to calculate the Ward capacitors.
4307  double DqiVgd(0.0); //
4308  double DqiVgb(0.0); // Cs, gs=dqs/dVgs= DqiVgs*fp+qiVgs*dfp
4309  double DqbVgs(0.0); // fp: partitioning factor
4310  double DqbVgd(0.0);
4311  double DqbVgb(0.0);
4312  double DqgVgs(0.0);
4313  double DqgVgd(0.0);
4314  double DqgVgb(0.0);
4315  double SigmaVgs(0.0); // This Sigma-parameter acounts for the DIBL
4316  double SigmaVgd(0.0); // effect on the capasitors. The physical
4317  double etavt(0.0);
4318  double VFB(0.0);
4319  double mqWL(0.0);
4320  double gamma(0.0);
4321  double A(0.0);
4322  double ALPHA(0.0);
4323 
4324  vbs_local = xvgs-vgb;
4325  vdds_local = xvgs-xvgdd;
4327  VFB =model_.vfb;
4328  ALPHA =model_.alpha;
4329 
4330  if(xvgs > VFB+vbs_local)
4331  {
4332  gamma= model_.gammas0;
4333  mqWL=(-CONSTQ)*EffectiveLength*(w);
4334  etavt=model_.eta*vt;
4335  nss=2.0*n0*log(1+0.5*exp((xvgs-von)/etavt));
4336  nsd=2.0*n0*log(1+0.5*exp((xvgs-von-ALPHA*vdds_local)/etavt));
4337 
4338  if (nss<1e-36) {nss=1.0e-36;}
4339  if (nsd<1e-36) {nsd=1.0e-36;}
4340  zetanb = (1-ALPHA)/ALPHA;
4341  xi = gamma*gamma/(A*A)+4/(A*A)*(vgb-VFB)-4.0/A*nss;
4342  xisqrt=sqrt(xi);
4343  tnsd = etavt/nsd+A;
4344  tnss = etavt/nss+A;
4345  if (vbs_local <= 0) { DeltaVT= 0.5*gamma/sqrt(tPhi-vbs_local); }
4346  if ((vbs_local > 0)&&(vbs_local <= 2*tPhi)) { DeltaVT=0.5*gamma/sqrt(tPhi); }
4347  if ((vbs_local > 0)&&(vbs_local > 2*tPhi)) { DeltaVT=0; }
4348 
4349  DnsdVgs=(1+DeltaVT-ALPHA+SigmaVgs)/tnsd;
4350  DnsdVgd=(-SigmaVgd+ALPHA)/tnsd;
4351  DnsdVgb=-DeltaVT/tnsd;
4352 
4353  DnssVgs= (1+DeltaVT+SigmaVgs)/tnss;
4354  DnssVgd=-SigmaVgd/tnss;
4355  DnssVgb=-DeltaVT/tnss;
4356 
4357  UI= etavt/2.0*(nsd*nsd-nss*nss) + A/3.0*(nsd*nsd*nsd-nss*nss*nss);
4358  VI= etavt*(nsd-nss) + A/2.0*(nsd*nsd-nss*nss);
4359  DUIVgs = etavt*(nsd*DnsdVgs-nss*DnssVgs) +
4360  A*(nsd*nsd*DnsdVgs-nss*nss*DnssVgs);
4361  DUIVgd= etavt*(nsd*DnsdVgd-nss*DnssVgd) +
4362  A*(nsd*nsd*DnsdVgd-nss*nss*DnssVgd);
4363  DUIVgb= etavt*(nsd*DnsdVgb-nss*DnssVgb) +
4364  A*(nsd*nsd*DnsdVgb-nss*nss*DnssVgb);
4365 
4366  DVIVgs= etavt*(DnsdVgs-DnssVgs) + A*(nsd*DnsdVgs-nss*DnssVgs);
4367  DVIVgd= etavt*(DnsdVgd-DnssVgd) + A*(nsd*DnsdVgd-nss*DnssVgd);
4368  DVIVgb= etavt*(DnsdVgb-DnssVgb) + A*(nsd*DnsdVgb-nss*DnssVgb);
4369  if (VI != 0)
4370  {
4371  DqiVgs= mqWL*(DUIVgs*VI-DVIVgs*UI)/(VI*VI);
4372  DqiVgd= mqWL*(DUIVgd*VI-DVIVgd*UI)/(VI*VI);
4373  DqiVgb= mqWL*(DUIVgb*VI-DVIVgb*UI)/(VI*VI);
4374  }
4375  else
4376  {
4377  DqiVgs=mqWL*DnssVgs;
4378  DqiVgd=mqWL*DnssVgd;
4379  DqiVgb=mqWL*DnssVgb;
4380  }
4381 
4382  DxiVgs= -4.0/A*DnssVgs;
4383  DxiVgd= -4.0/A*DnssVgd;
4384  DxiVgb= 4.0/(A*A)-4/A*DnssVgb;
4385 
4386  if (xisqrt !=0)
4387  {DndepVgs= gamma/4.0*1.0/xisqrt*DxiVgs;
4388  DndepVgd= gamma/4.0*1.0/xisqrt*DxiVgd;
4389  DndepVgb= gamma/4.0*1.0/xisqrt*DxiVgb;
4390  }
4391  else
4392  {
4393  DndepVgs=0;
4394  DndepVgd=0;
4395  DndepVgb=0;
4396  }
4397 
4398  DqbVgs = mqWL*(DndepVgs-zetanb*DnssVgs)+zetanb*DqiVgs;
4399  DqbVgd = mqWL*(DndepVgd-zetanb*DnssVgd)+zetanb*DqiVgd;
4400  DqbVgb = mqWL*(DndepVgb-zetanb*DnssVgb)+zetanb*DqiVgb;
4401 
4402  DqgVgs= -DqiVgs-DqbVgs;
4403  DqgVgd= -DqiVgd-DqbVgd;
4404  DqgVgb= -DqiVgb-DqbVgb;
4405  }
4406  else
4407  {
4408  // Below flat band
4409  DqgVgs=0.0;
4410  DqgVgd=0.0;
4412  }
4413 
4414  if (DqgVgs<0) {DqgVgs=0;} // Because of truncation error cgb may turn
4415  if (DqgVgd<0) {DqgVgd=0;} // negative in extreme inversion, while cgd
4416  if (DqgVgb<0) {DqgVgb=0;} // and cgs may turn negative in accumulation
4417  cgs= 0.5*DqgVgs;
4418  cgd= 0.5*DqgVgd;
4419  cgb= 0.5*DqgVgb;
4420 
4421  return true;
4422 }
4423 
4424 //----------------------------------------------------------------------------
4425 // Function : Instance::UCCMCharges
4426 // Purpose :
4427 //
4428 //
4429 // Special Notes :
4430 //
4431 // Scope : public
4432 // Creator : pmc
4433 // Creation Date : 02/17/2004
4434 //----------------------------------------------------------------------------
4435 bool Instance::UCCMCharges( double vgs, double vgdd,
4436  double vgb, double & qD, double & qS, double & qB)
4437 {
4438  double VT(0.0); // Threshold voltage
4439  double nss(0.0); // Inversion charge density at source
4440  double nsd(0.0); // Inversion charge density at drain
4441  double etavth(0.0);
4442  double mqWL(0.0);
4443  double nsdsqr(0.0);
4444  double nsssqr(0.0);
4445  double arg1(0.0);
4446  double arg2(0.0);
4447  double qn(0.0); // Inversion charge
4448  double VFB(0.0); // Flat-band voltage
4449  double xisqrt(0.0);
4450  double ndeps(0.0); // Depletion charge density
4451  double zetanb(0.0);
4452  double fp1(0.0); // Partitioning factor
4453  double gamma(0.0);
4454  double ALPHA(0.0);
4455  double A(0.0);
4456  double cox(0.0); // oxide capasity per unit area
4457  double vdds_local(0.0);
4458 
4459  cox = model_.oxideCapFactor;
4460  gamma = model_.gammas0;
4461  A = CONSTQ/cox;
4462  ETA = model_.eta;
4463  ALPHA = model_.alpha;
4464  VFB = model_.vfb;
4465 
4466  vdds_local=vgs-vgdd;
4467 
4468  if (vgb > VFB)
4469  {
4470  VT = von;
4471  etavth = ETA*temp*CONSTKoverQ;
4472 
4473  nsd=2*n0*log(1+0.5*exp((vgs-VT-ALPHA*vdds_local)/etavth));
4474  nss=2*n0*log(1+0.5*exp((vgs-VT)/etavth));
4475 
4476  nsdsqr=nsd*nsd;
4477  nsssqr=nss*nss;
4478  mqWL=-CONSTQ*EffectiveLength*(w);
4479 
4480  arg1= 0.5*etavth*(nsdsqr-nsssqr)+1/3.0*A*(nsdsqr*nsd-nsssqr*nss);
4481  arg2= etavth*(nsd-nss) +0.5*A*(nsdsqr-nsssqr);
4482 
4483  if (arg2==0) //||(fabs(nss-nsd)<1e-12))
4484  { qn=mqWL*nss;}
4485  else {qn = mqWL*arg1/arg2;}
4486  xisqrt= sqrt(gamma*gamma/(A*A)+4.0/(A*A)*(vgb-VFB)-4.0/A*nss);
4487  ndeps = -gamma*gamma/(2.0*A)+ gamma/2.0*xisqrt;
4488  zetanb =(1-ALPHA)/ALPHA;
4489  mqWL=-CONSTQ*EffectiveLength*(w);
4490  qB =mqWL*(ndeps-zetanb*nss)+zetanb*qn;
4491 
4492  switch(model_.fpe)
4493  {
4494  case 1:
4495  fp1=model_.xqc;
4496  break;
4497 
4498  case 2:
4499  {
4500  double vsat(0.0);
4501  double beta(0.0);
4502  double a0(0.0);
4503  double vdsabs(0.0);
4504  double m(0.0);
4505  double fp1denum(0.0);
4506 
4507  vdsabs=fabs(vdds_local);
4508  vsat=vddsat;
4509  if (vgs-VT > 1e-36)
4510  {
4511  if((vsat > 1e-36)&&(vdds_local>1e-36))
4512  { beta=vdsabs/vsat;}else{beta=1e-36;}
4513  m = model_.mcv;
4514  a0 = (model_.xqc - 0.5)/vsat;
4515  fp1denum = exp(1/m*log(1+exp(m*log(beta))));
4516  if (fp1denum > 1e-36)
4517  {fp1= 0.5+ a0*vdsabs/fp1denum;}
4518  else {fp1=0.5;}
4519 
4520  } else {fp1=0.5;}
4521  }
4522  break;
4523  case 3:
4524  {
4525  double nd0(0.0);
4526  double ns0(0.0);
4527  double arg1_loc(0.0);
4528  double arg2_loc(0.0);
4529  double arg3_loc(0.0);
4530 
4531  nd0=nss/(A*etavth);
4532  ns0=nss/(A*etavth);
4533  arg1_loc= nd0*nd0*nd0*(1/3.0+nd0*(3.0/8.0+1/10.0*nd0))
4534  - (1.0/2.0+1.0/3.0*nd0)*(1.0+1/2.0*ns0)*nd0*nd0*ns0
4535  + ns0*ns0*ns0*(1.0/6.0+ns0*(5.0/24.0+1.0/15.0*ns0));
4536  arg2_loc = 1.0/2.0*(nd0*nd0-ns0*ns0)+1.0/3.0*(nd0*nd0*nd0-ns0*ns0*ns0);
4537  arg3_loc = nd0-ns0+1.0/2.0*(nd0*nd0-ns0*ns0);
4538  if ((arg2_loc==0)||(arg3_loc==0)) {fp1=0.5;}
4539  else { fp1 = 1-arg1_loc/(arg2_loc*arg3_loc);}
4540  }
4541  break;
4542  default:
4543  {
4544  UserWarning(*this) << "Partitioning model does not exist";
4545  return true;
4546  }
4547  }
4548  }
4549  else
4550  { // Below flatband
4551  qB = -EffectiveLength*(w)*cox*(vgb-VFB);
4552  qn = 0;
4553  }
4554 
4555  qS = -fp1*qn;
4556  qD = -(1-fp1)*qn;
4557  qB = -qB;
4558 
4559  return true;
4560 }
4561 
4562 //----------------------------------------------------------------------------
4563 // Function : Instance::UCCMcvon
4564 // Purpose : Compute the threshold voltage
4565 //
4566 // Special Notes : ERK. The "local" named variables (von_local, dvonvbs_local)
4567 // are named that way to avoid conflicts with the instance
4568 // variables of the same name.
4569 //
4570 // Scope : public
4571 // Creator : pmc
4572 // Creation Date : 02/17/2004
4573 //----------------------------------------------------------------------------
4574 bool Instance::UCCMcvon(double vbs_local,
4575  double & von_local, double & dvonvbs_local)
4576 {
4577  double PhiMinVbs = tPhi - vbs_local;
4578  double sarg(0.0);
4579  double dsrgdb(0.0);
4580 
4581 // vtoo, vthLimit calculated in updateTemp
4582 
4583  if(vtoo-vbs_local > vthLimit) {
4584  von_local = vtoo + gammal*vthLimit;
4585  dvonvbs_local = 0.0;
4586  return true;
4587  }
4588  if(PhiMinVbs > 0.0)
4589  {
4590  sarg = sqrt(PhiMinVbs);
4591  dsrgdb = -0.5/sarg;
4592  } else
4593  {
4594  sarg = 0;
4595  dsrgdb = 0;
4596  }
4597  von_local = vtoo + gammas*sarg - gammal*PhiMinVbs;
4598  dvonvbs_local = gammas*dsrgdb + gammal;
4599 
4600  return true;
4601 }
4602 
4603 //-----------------------------------------------------------------------------
4604 // Function : Instance::UCCMmosa1
4605 // Purpose : Compute currents and conductances in drain region
4606 //
4607 //
4608 // Special Notes :
4609 //
4610 // Scope : public
4611 // Creator : pmc
4612 // Creation Date : 02/17/2004
4613 //-----------------------------------------------------------------------------
4614 
4615 bool Instance::UCCMmosa1(double xvgs, double xvdds,
4616  double dvonvbs, double & cdraindrift_loc, double & vsate)
4617 {
4618  double etavt(0.0);
4619  double vgt(0.0);
4620  double vgt0(0.0);
4621  double sigma(0.0);
4622  double vgte(0.0);
4623  double isat(0.0);
4624  double mu(0.0);
4625  double ns(0.0);
4626  double a(0.0);
4627  double b(0.0);
4628  double d(0.0);
4629  double g(0.0);
4630  double h(0.0);
4631  double q(0.0);
4632  double r(0.0);
4633  double t(0.0);
4634  double u(0.0);
4635  double x(0.0),y(0.0),z(0.0);
4636 
4637  double gch(0.0);
4638  double gchi(0.0);
4639  double rt(0.0);
4640  double vl(0.0);
4641  double vl2(0.0);
4642  double ichoo(0.0);
4643 
4644  double dichoodvds(0.0);
4645  double dichoodvgs(0.0);
4646  double dichoodvbs(0.0);
4647  double dichooisat(0.0);
4648  double dichoodgch(0.0);
4649  double delgchgchi(0.0);
4650  double disatdvds(0.0);
4651  double disatdvgs(0.0);
4652  double disatdvbs(0.0);
4653  double dgchidvds(0.0);
4654  double dgchidvgs(0.0);
4655  double dgchidvbs(0.0);
4656  double disatvgte(0.0);
4657  double disatgchi(0.0);
4658  double dvgtedvgt(0.0);
4659  double dvgtdvds(0.0);
4660  double dvgtdvgs(0.0);
4661  double dvgtdvbs(0.0);
4662  double dnsdvgt(0.0);
4663  double dnsdvds(0.0);
4664  double dnsdvgs(0.0);
4665  double dnsdvbs(0.0);
4666  double dmudvgte(0.0);
4667  double dmudvds(0.0);
4668  double dmudvgs(0.0);
4669  double dmudvbs(0.0);
4670 
4671  static int output=0;
4672 
4673 #define ETA model_.eta
4674 #define RS model_.sourceResistance
4675 #define RD model_.drainResistance
4676 #define VSIGMAT model_.vsigmat
4677 #define VSIGMA model_.vsigma
4678 #define SIGMA0 model_.sigma0
4679 #define THETA model_.theta
4680 #define LAMBDA model_.lambda
4681 #define DELTA model_.delta
4682 #define VMAX model_.maxDriftVel
4683 #define TOX model_.oxideThickness
4684 
4685 #define EXP_MAX 150.0
4686 
4687  if(output == 1) {
4688  Xyce::dout() << " " << std::endl;
4689  Xyce::dout() << "ETA = " << ETA << std::endl;
4690  Xyce::dout() << "RS = " << RS << std::endl;
4691  Xyce::dout() << "RD = " << RD << std::endl;
4692  Xyce::dout() << "VSIGMAT = " << VSIGMAT << std::endl;
4693  Xyce::dout() << "VSIGMA = " << VSIGMA << std::endl;
4694  Xyce::dout() << "SIGMA0 = " << SIGMA0 << std::endl;
4695  Xyce::dout() << "THETA = " << THETA << std::endl;
4696  Xyce::dout() << "LAMBDA = " << LAMBDA << std::endl;
4697  Xyce::dout() << "DELTA = " << DELTA << std::endl;
4698  Xyce::dout() << "VMAX = " << VMAX << std::endl;
4699  Xyce::dout() << "TOX = " << TOX << std::endl;
4700  output=0;
4701  }
4702 
4703  etavt = ETA*vt;
4704  rt = RS+RD;
4705  vgt0 = xvgs - von;
4706  a = exp((vgt0-VSIGMAT)/VSIGMA);
4707  sigma = SIGMA0/(1+a);
4708  vgt = vgt0+sigma*xvdds;
4709  b = 0.5*vgt/vt-1;
4710  q = sqrt(model_.deltaSqr+b*b);
4711  vgte = vt*(1+b+1+q);
4712  u = 1+THETA*(vgte+2*von)/TOX;
4713  mu = tSurfMob/u;
4714 
4715  x = vgt/etavt;
4716  if (x > 50.0)
4717  {
4718  ns = n0*2.0*x;
4719  }
4720  else if (x < -30)
4721  {
4722  ns = n0*exp(x);
4723  }
4724  else
4725  {
4726  ns = 2.0*n0*log(1+0.5*exp(x));
4727  }
4728 
4729  if(ns < 1.0e-38) {
4730  cdraindrift_loc = 0.0;
4731  vsate = 0.0;
4732  gm = 0.0;
4733  gdds = 0.0;
4734  gmbs = 0.0;
4735  mm1 = 0.0;
4736  dmm1vgs = 0.0;
4737  dmm1vds = 0.0;
4738  dmm1vbs = 0.0;
4739  return true;
4740  }
4741 
4742  gchi = gchi0*mu*ns;
4743  gch = gchi/(1+gchi*rt);
4744  t = VMAX*L;
4745 // vl = t/mu;
4746  vl = t/tSurfMob;
4747  vl2 = vl*vl;
4748  d = sqrt(1+2*gchi*RS + vgte*vgte/vl2);
4749  r = gchi*vgte;
4750  h = 1+gchi*RS+d;
4751  isat = r/h;
4752  vsate = isat/gch;
4753  y = xvdds/(vsate);
4754 
4755  if(fabs(y) > EXP_MAX) z = 0;
4756  else z = 1/(cosh(y)*cosh(y));
4757 
4758  // drain current
4759  ichoo = isat*(1+LAMBDA*xvdds)*tanh(y);
4760 
4761  dvgtedvgt = 0.5*(1+b/q);
4762  dnsdvgt = n0/(etavt*(exp(-x)+0.5));
4763  dmudvgte = -tSurfMob*THETA/(TOX*u*u);
4764 
4765  dvgtdvds = sigma;
4766  dvgtdvgs = 1 - SIGMA0*a*xvdds/(VSIGMA*(1+a)*(1+a));
4767  dvgtdvbs = -dvonvbs*dvgtdvgs;
4768 
4769  dnsdvds = dnsdvgt*dvgtdvds;
4770  dnsdvgs = dnsdvgt*dvgtdvgs;
4771  dnsdvbs = dnsdvgt*dvgtdvbs;
4772  dmudvds = (dmudvgte*dvgtedvgt)*dvgtdvds;;
4773  dmudvgs = (dmudvgte*dvgtedvgt)*dvgtdvgs;;
4774  dmudvbs = (dmudvgte*dvgtedvgt)*dvgtdvbs + 2*dmudvgte*dvonvbs;
4775 
4776  dgchidvds = gchi0*(mu*dnsdvds + ns*dmudvds);
4777  dgchidvgs = gchi0*(mu*dnsdvgs + ns*dmudvgs);
4778  dgchidvbs = gchi0*(mu*dnsdvbs + ns*dmudvbs);
4779 
4780  disatvgte = gchi/h - gchi*vgte*vgte/vl2/(d*h*h);
4781  disatgchi = vgte/h - gchi*vgte*RS*(1+1/d)/(h*h);
4782 
4783  disatdvds = (disatvgte*dvgtedvgt)*dvgtdvds + disatgchi*dgchidvds;
4784  disatdvgs = (disatvgte*dvgtedvgt)*dvgtdvgs + disatgchi*dgchidvgs;
4785  disatdvbs = (disatvgte*dvgtedvgt)*dvgtdvbs + disatgchi*dgchidvbs;
4786 
4787  dichoodgch = (1+LAMBDA*xvdds)*xvdds*z;
4788  dichooisat = (1+LAMBDA*xvdds)*(tanh(y) - gch*xvdds*z/isat);
4789  g = 1+gchi*rt;
4790  delgchgchi = 1/(g*g);
4791 
4792  dichoodvds = dichooisat*disatdvds + (dichoodgch*delgchgchi)*dgchidvds
4793  + isat*LAMBDA*tanh(y) + (1+LAMBDA*xvdds)*gch*z;
4794  dichoodvgs = dichooisat*disatdvgs + (dichoodgch*delgchgchi)*dgchidvgs;
4795  dichoodvbs = dichooisat*disatdvbs + (dichoodgch*delgchgchi)*dgchidvbs;
4796 
4797 
4798  cdraindrift_loc = ichoo;
4799  gm = dichoodvgs;
4800  gdds = dichoodvds;
4801  gmbs = dichoodvbs;
4802 
4803 // the mm1 factors are for the ISUB calc which we don't use
4804  mm1 = 0.0;
4805  dmm1vgs = 0.0;
4806  dmm1vds = 0.0;
4807  dmm1vbs = 0.0;
4808 
4809  return true;
4810 }
4811 
4812 
4813 //-----------------------------------------------------------------------------
4814 // Function : Instance::UCCMmosa2
4815 // Purpose : Compute currents and conductances in
4816 // : drain and substrate regions
4817 //
4818 //
4819 // Special Notes :
4820 //
4821 // Scope : public
4822 // Creator : pmc
4823 // Creation Date : 02/17/2004
4824 //-----------------------------------------------------------------------------
4825 
4826 bool Instance::UCCMmosa2(double xvgs, double xvdds,
4827  double dvonvbs, double & cdraindrift_loc, double & vsate)
4828 {
4829  double etavt(0.0);
4830  double vgt(0.0);
4831  double vgt0(0.0);
4832  double sigma(0.0);
4833  double vgte(0.0);
4834  double isat(0.0);
4835  double mu(0.0);
4836  double ns(0.0);
4837  double a(0.0);
4838  double b(0.0);
4839  double c(0.0);
4840  double d(0.0);
4841  double e(0.0);
4842  double f(0.0);
4843  double g(0.0);
4844  double h(0.0);
4845  double p(0.0);
4846  double q(0.0);
4847  double r(0.0);
4848  double s(0.0);
4849  double t(0.0);
4850  double u(0.0);
4851  double x,y,z(0.0);
4852  double gch(0.0);
4853  double gchi(0.0);
4854  double rt(0.0);
4855  double vl(0.0);
4856  double vl2(0.0);
4857  double ichoo(0.0);
4858  double gmoo(0.0);
4859  double gdsoo(0.0);
4860  double gbsoo(0.0);
4861  double icho(0.0);
4862  double gmo(0.0);
4863  double gdso(0.0);
4864  double gbso(0.0);
4865  double vsat(0.0);
4866  double vdso(0.0);
4867  double vdse(0.0);
4868  double cox(0.0);
4869  double temp1(0.0);
4870  double qs(0.0);
4871  double dqsvgt(0.0);
4872  double dqsvbs(0.0);
4873  double delidgch(0.0);
4874  double delgchgchi(0.0);
4875  double delgchins(0.0);
4876  double dvgtevgt(0.0);
4877  double delidvsate(0.0);
4878  double delvsateisat(0.0);
4879  double delisatvgte(0.0);
4880  double delisatgchi(0.0);
4881  double delisatvl(0.0);
4882  double dgchivgt(0.0);
4883  double delvsategch(0.0);
4884  double delidvds(0.0);
4885  double dsigmavgs(0.0);
4886  double delnsvgt(0.0);
4887  double dvsatevgt(0.0);
4888  double dvgtvgs(0.0);
4889  double dmuvon(0.0);
4890  double dmuvgt(0.0);
4891  double dvgtvon(0.0);
4892  double delvsatqs(0.0);
4893  double delvsatmu(0.0);
4894  double dvsatvgt(0.0);
4895  double dvsatvbs(0.0);
4896  double delvdsevdso(0.0);
4897  double delvdsevsat(0.0);
4898  double dvdsevgs(0.0);
4899  double dvdsevds(0.0);
4900  double dvdsevbs(0.0);
4901 
4902  double del(0.0);
4903  double delmm1vdso(0.0);
4904  double ddelvgs(0.0);
4905  double ddelvds(0.0);
4906  double ddelvbs(0.0);
4907 
4908  // These "loc" variables are named such to avoid name conflicts with
4909  // instance variables.
4910  double mm1_loc(0.0);
4911  double dmm1vgs_loc(0.0);
4912  double dmm1vds_loc(0.0);
4913  double dmm1vbs_loc(0.0);
4914 
4915  static int output(0);
4916 
4917 #define ETA model_.eta
4918 #define RS model_.sourceResistance
4919 #define RD model_.drainResistance
4920 #define VSIGMAT model_.vsigmat
4921 #define VSIGMA model_.vsigma
4922 #define SIGMA0 model_.sigma0
4923 #define THETA model_.theta
4924 #define LAMBDA model_.lambda
4925 #define DELTA model_.delta
4926 #define ALPHA model_.alpha
4927 #define VMAX model_.maxDriftVel
4928 #define TOX model_.oxideThickness
4929 #define W w
4930 #define L EffectiveLength
4931 #define LS model_.ls
4932 #define M model_.m
4933 #define INVM 1.0/model_.m
4934 #define MC model_.mc
4935 #define INVMC 1.0/model_.mc
4936 #define RSUB model_.rsub
4937 #define AI model_.ai
4938 #define BI model_.bi
4939 #define MD model_.md
4940 #define INVMD 1.0/model_.md
4941 #define DELMAX model_.delmax
4942 
4943  if(output == 1) {
4944  Xyce::dout() << " " << std::endl;
4945  Xyce::dout() << "ETA = " << ETA << std::endl;
4946  Xyce::dout() << "RS = " << RS << std::endl;
4947  Xyce::dout() << "RD = " << RD << std::endl;
4948  Xyce::dout() << "VSIGMAT = " << VSIGMAT << std::endl;
4949  Xyce::dout() << "VSIGMA = " << VSIGMA << std::endl;
4950  Xyce::dout() << "SIGMA0 = " << SIGMA0 << std::endl;
4951  Xyce::dout() << "THETA = " << THETA << std::endl;
4952  Xyce::dout() << "LAMBDA = " << LAMBDA << std::endl;
4953  Xyce::dout() << "DELTA = " << DELTA << std::endl;
4954  Xyce::dout() << "ALPHA = " << ALPHA << std::endl;
4955  Xyce::dout() << "VMAX = " << VMAX << std::endl;
4956  Xyce::dout() << "TOX = " << TOX << std::endl;
4957  Xyce::dout() << "W = " << W << std::endl;
4958  Xyce::dout() << "L = " << L << std::endl;
4959  Xyce::dout() << "LS = " << LS << std::endl;
4960  Xyce::dout() << "M = " << M << std::endl;
4961  Xyce::dout() << "INVM = " << INVM << std::endl;
4962  Xyce::dout() << "MC = " << MC << std::endl;
4963  Xyce::dout() << "INVMC = " << INVMC << std::endl;
4964  Xyce::dout() << "RSUB = " << RSUB << std::endl;
4965  Xyce::dout() << "AI = " << AI << std::endl;
4966  Xyce::dout() << "BI = " << BI << std::endl;
4967  Xyce::dout() << "MD = " << MD << std::endl;
4968  Xyce::dout() << "INVMD = " << INVMD << std::endl;
4969  Xyce::dout() << "DELMAX = " << DELMAX << std::endl;
4970  Xyce::dout() << "VP = " << vp << std::endl;
4971  output=0;
4972  }
4973 
4974 
4975  etavt = ETA*vt;
4976  rt = RS+RD;
4977  vgt0 = xvgs - von;
4978  a = exp((vgt0-VSIGMAT)/VSIGMA);
4979  sigma = SIGMA0/(1+a);
4980  vgt = vgt0+sigma*xvdds;
4981  b = 0.5*vgt/vt-1;
4982  q = sqrt(model_.deltaSqr+b*b);
4983  vgte = vt*(1+b+1+q);
4984  u = 1+THETA*(vgte+2*von)/TOX;
4985  mu = tSurfMob/u;
4986 
4987  x = vgt/etavt;
4988  if (x > 50.0)
4989  {
4990  ns = n0*2.0*x;
4991  }
4992  else if (x < -30)
4993  {
4994  ns = n0*exp(x);
4995  }
4996  else
4997  {
4998  ns = 2.0*n0*log(1+0.5*exp(x));
4999  }
5000  //ns = 2.0*n0*log(1+0.5*c);
5001 
5002  if(ns < 1.0e-38)
5003  {
5004  cdraindrift_loc = 0.0;
5005  vsate = 0.0;
5006  gm = 0.0;
5007  gdds = 0.0;
5008  gmbs = 0.0;
5009 
5010  mm1 = mm1_loc = 0.0;
5011  dmm1vgs = dmm1vgs_loc = 0.0;
5012  dmm1vds = dmm1vds_loc = 0.0;
5013  dmm1vbs = dmm1vbs_loc = 0.0;
5014  return true;
5015  }
5016 
5017  gchi = gchi0*mu*ns;
5018  gch = gchi/(1+gchi*rt);
5019  t = VMAX*L;
5020  vl = t/mu;
5021  vl2 = vl*vl;
5022  d = sqrt(1+2*gchi*RS + vgte*vgte/vl2);
5023  r = gchi*vgte;
5024  isat = r/(1+gchi*RS+d);
5025  vsate = isat/gch;
5026  y = xvdds/(vsate);
5027  z = 0.5*(cosh(2*y)+1);
5028 
5029 // Old model
5030 // e = pow(xvdds/(*vsate),M);
5031 // f = pow(1+e,INVM);
5032 // delidgch = xvdds*(1+LAMBDA*xvdds)/f;
5033 
5034  ichoo = isat*(1+LAMBDA*xvdds)*tanh(y);
5035  delidgch = ichoo/gch;
5036  g = 1+gchi*rt;
5037  delgchgchi = 1/(g*g);
5038  delgchins = gchi0*mu;
5039  x = vgt/etavt;
5040  delnsvgt = n0/(etavt*(exp(-x) + 0.5));
5041  dvgtevgt = 0.5*(1+b/q);
5042 
5043 // delidvds = gch*(1+2*LAMBDA*xvdds)/f-
5044 // ichoo*pow(vdds/(*vsate),M-1)/((*vsate)*(1+e));
5045 // delidvsate = ichoo*e/((*vsate)*(1+e));
5046  delidvds = ichoo*LAMBDA/(1+LAMBDA*xvdds) + gch*(1+LAMBDA*xvdds)/z;
5047  delidvsate = gch*(1+LAMBDA*xvdds)*xvdds/(vsate)/z;
5048 
5049  delvsateisat = 1/gch;
5050  delvsategch = -(vsate)/gch;
5051  h = 1+gchi*RS+d;
5052  delisatgchi = (vgte*h - gchi*vgte*RS*(1+1/d))/(h*h);
5053  delisatvgte = (gchi*h - gchi*vgte*vgte/(vl2*d))/(h*h);
5054  delisatvl = isat*isat*vgte/(gchi*d*vl2*vl);
5055  s = THETA/(tSurfMob*TOX);
5056  dgchivgt = delgchins*delnsvgt-gchi*mu*s*dvgtevgt;
5057  dsigmavgs = -SIGMA0*a/(VSIGMA*((1+a)*(1+a)));
5058  p = delgchgchi*dgchivgt;
5059  dvsatevgt = delvsateisat*(delisatgchi*dgchivgt+
5060  (delisatvgte+delisatvl*t*s)*dvgtevgt)+ delvsategch*p;
5061  g = delidgch*p+delidvsate*dvsatevgt;
5062  dvgtvgs = (1+xvdds*dsigmavgs);
5063  gmoo = g*dvgtvgs;
5064  gdsoo = delidvds+g*sigma;
5065 
5066  dmuvon = -(2-dvgtevgt*dvgtvgs)*mu*mu*s;
5067  dmuvgt = -mu*mu*THETA*dvgtevgt/(tSurfMob*TOX);
5068  dvgtvon = -dvgtvgs;
5069 
5070  if(THETA == 0.0)
5071  {
5072  gbsoo = -(gmoo)*dvonvbs;
5073  }
5074  else
5075  {
5076  double dgchivon = -delgchins*delnsvgt*dvgtvgs+gchi*dmuvon/mu;
5077  p = delgchgchi*dgchivon;
5078  double dvsatevon = delvsateisat*(delisatgchi*dgchivon+
5079  delisatvgte*dvgtevgt*dvgtvon+delisatvl*(-vl2/t)*dmuvon)+delvsategch*p;
5080  gbsoo = (delidgch*p+delidvsate*dvsatevon)*dvonvbs;
5081  }
5082 
5083 // mosa2: New Version
5084 
5085  cox = CONSTEPSOX/TOX;
5086  temp1 = cox*RS/CONSTQ;
5087  qs = CONSTQ*(ns-ichoo*temp1);
5088  dqsvgt = CONSTQ*(delnsvgt-g*temp1);
5089  dqsvbs = CONSTQ*(delnsvgt*dvgtvon*dvonvbs-gbsoo*temp1);
5090 
5091  // A more precise vsat model for deltal and isub calculations
5092  if(model_.dtype == CONSTNMOS)
5093  {
5094  a = 2*qs*VMAX*L;
5095  b = qs*mu+2*cox*ALPHA*VMAX*L;
5096  temp1 = b*b;
5097  vsat = a/b;
5098  delvsatqs = 2*VMAX*L*(b-qs*mu)/temp1;
5099  delvsatmu = -a*qs/temp1;
5100  dvsatvgt = delvsatqs*dqsvgt+delvsatmu*dmuvgt;
5101  dvsatvbs = delvsatqs*dqsvbs+delvsatmu*dmuvon*dvonvbs;
5102  }
5103  else
5104  {
5105  a = VMAX*L/mu;
5106  b = 2/(cox*ALPHA*VMAX*L);
5107  c = sqrt(1+qs*mu*b);
5108  vsat = a*(c-1);
5109  delvsatqs = a*(0.5*b*mu/c);
5110  delvsatmu = -vsat/mu+a*(0.5*b*qs/c);
5111  dvsatvgt = delvsatqs*dqsvgt+delvsatmu*dmuvgt;
5112  dvsatvbs = delvsatqs*dqsvbs+delvsatmu*dmuvon*dvonvbs;
5113  }
5114 
5115  // Effective intrinsic drain-source voltage calculation
5116  vdso = xvdds-ichoo*rt;
5117  a = pow(vdso/vsat,MC);
5118  b = pow(1+a,INVMC);
5119  vdse = vdso/b;
5120  delvdsevdso = (1-a/(1+a))/b;
5121  delvdsevsat = vdse*a/(vsat*(1+a));
5122 
5123  // DeltaL calculation
5124  double deltal(0.0);
5125  double deldeltalvdse(0.0);
5126  double deldeltalvdso(0.0);
5127  double deldeltalmu(0.0);
5128  double ddeltalvgs(0.0);
5129  double ddeltalvds(0.0);
5130  double ddeltalvbs(0.0);
5131 
5132  a = 1+(vdso-vdse)/vp;
5133  b = W*mu*cox*vdse*RS/L;
5134  c = 1+vdse/vp+b;
5135  d = LS*log10(a);
5136  deltal = d/c;
5137  e = 1/(1-deltal/L);
5138  icho = ichoo*e;
5139  f = log(10.0);
5140  deldeltalvdse = (-LS/(f*a*vp)*c-d*(1/vp+W*mu*cox*RS/L))/(c*c);
5141  deldeltalvdso = LS/(c*a*f*vp);
5142  deldeltalmu = -d*W*cox*vdse*RS/(L*c*c);
5143  temp1 = -gmoo*rt;
5144  dvdsevgs = delvdsevsat*dvsatvgt*dvgtvgs+delvdsevdso*temp1;
5145  ddeltalvgs = deldeltalvdse*dvdsevgs + deldeltalmu*dmuvgt*dvgtvgs
5146  + deldeltalvdso*temp1;
5147  gmo = e*(gmoo+ichoo*ddeltalvgs*e/L);
5148  temp1 = 1-gdsoo*rt;
5149  dvdsevds = delvdsevsat*dvsatvgt*sigma+delvdsevdso*temp1;
5150  ddeltalvds = deldeltalvdse*dvdsevds + deldeltalmu*dmuvgt*sigma
5151  + deldeltalvdso*temp1;
5152  gdso = e*(gdsoo+ichoo*ddeltalvds*e/L);
5153  temp = -gbsoo*rt;
5154  dvdsevbs = delvdsevsat*dvsatvbs+delvdsevdso*temp1;
5155  ddeltalvbs = deldeltalvdse*dvdsevbs + deldeltalmu*dmuvon*dvonvbs
5156  + deldeltalvdso*temp1;;
5157  gbso = e*(gbsoo+ichoo*ddeltalvbs*e/L);
5158 
5159  a = AI/BI;
5160  b = vdso-vdse+vt;
5161  c = exp(-LS*BI/b);
5162  mm1_loc = a*b*c;
5163  if(RSUB != 0.0)
5164  {
5165  d = 0.5*W*(ALPHA-1)*cox*RSUB/L;
5166  del = d*mm1_loc*mu*vdse;
5167  f = pow(del/DELMAX,MD);
5168  g = pow(1+f,INVMD);
5169  h = (1-f/(1+f))/g;
5170  del = del/g;
5171  e = 1/(1-del);//1+del;
5172  }
5173  else
5174  {
5175  d = 0;
5176  h = 1;
5177  e = 1;
5178  }
5179 
5180  cdraindrift_loc = icho*e;
5181  delmm1vdso = a*c+mm1_loc*LS*BI/(b*b);
5182  dmm1vgs_loc = delmm1vdso*(-gmoo*rt-dvdsevgs);
5183  ddelvgs = h*d*(mu*vdse*dmm1vgs_loc + mm1_loc*vdse*dmuvgt*dvgtvgs
5184  + mu*mm1_loc*dvdsevgs);
5185  gm = e*(gmo+icho*ddelvgs*e); //gmo*e+icho*ddelvgs;
5186  dmm1vds_loc = delmm1vdso*(1-gdsoo*rt-dvdsevds);
5187  ddelvds = h*d*(mu*vdse*dmm1vds_loc + mm1_loc*vdse*dmuvgt*sigma
5188  + mu*mm1_loc*dvdsevds);
5189  gdds = e*(gdso+icho*ddelvds*e); //gdso*e+icho*ddelvds;
5190  dmm1vbs_loc = delmm1vdso*(-gbsoo*rt-dvdsevbs);
5191  ddelvbs = h*d*(mu*vdse*dmm1vbs_loc + mm1_loc*vdse*dmuvon*dvonvbs
5192  + mu*mm1_loc*dvdsevbs);
5193  gmbs = e*(gbso+icho*ddelvbs*e); //gbso*e+icho*ddelvbs;
5194 
5195 
5196  // saving some quantities for substrate current calculations
5197  mm1 = mm1_loc;
5198  dmm1vgs = dmm1vgs_loc;
5199  dmm1vds = dmm1vds_loc;
5200  dmm1vbs = dmm1vbs_loc;
5201 
5202  return true;
5203 }
5204 
5205 
5206 //-----------------------------------------------------------------------------
5207 // VDMOS Master functions:
5208 //-----------------------------------------------------------------------------
5209 
5210 //-----------------------------------------------------------------------------
5211 // Function : Master::updateState
5212 // Purpose :
5213 // Special Notes :
5214 // Scope : public
5215 // Creator : Eric Keiter, SNL
5216 // Creation Date : 11/26/08
5217 //-----------------------------------------------------------------------------
5218 bool Master::updateState (double * solVec, double * staVec, double * stoVec)
5219 {
5220  bool bsuccess = true;
5221 
5222  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
5223  {
5224  Instance & mi = *(*it);
5225 
5226  double * currStaVector = mi.extData.currStaVectorRawPtr;
5227 
5228  bool btmp = mi.updateIntermediateVars ();
5229  bsuccess = bsuccess && btmp;
5230 
5231  double vgs1(0.0), vgdd1(0.0), vbs1(0.0),vgb1(0.0), vdds1(0.0);
5232 
5233  // voltage drops:
5234  staVec[mi.li_state_vbdd] = mi.vbdd;
5235  staVec[mi.li_state_vbs] = mi.vbs;
5236  staVec[mi.li_state_vgps] = mi.vgps;
5237  staVec[mi.li_state_vdds] = mi.vdds;
5238  staVec[mi.li_state_D1vd] = mi.D1vd;
5239 
5240  // now the meyer capacitances
5241  // we didn't calculate these charges in update IntermediateVars
5242  // but we did calculate the voltage drops and capacitances.
5243  // first store the capacitances themselves:
5244  staVec[mi.li_state_capgs] = mi.capgs;
5245  staVec[mi.li_state_capgdd] = mi.capgdd;
5246  staVec[mi.li_state_capgb] = mi.capgb;
5247 
5248  // now the charges
5249  // BE CAREFUL! We can only do Q=CV for DCOP! Otherwise it's
5250  // supposed to be *INTEGRATED*:
5251  // Q = int(t0,t1)C(V)*dV --- and we approximate that by
5252  // Q(t1)-Q(t0) = CBar*(V(t1)-V(t0)) where CBar is the average.
5253  // Now with Meyer back averaging, Capxx is the average between the last
5254  // time step and this one. So we gotta do the right thing for non-DCOP
5255  // when backaverage is on.
5256 
5257  if((getSolverState().dcopFlag))
5258  {
5259  mi.qgs = mi.Capgs *mi.vgps;
5260  mi.qgdd = mi.Capgdd*mi.vgpdd;
5261  mi.qgb = mi.Capgb *mi.Vgpb;
5262  }
5263  else
5264  {
5265  // get the ones from last time step
5266  mi.qgs = currStaVector[mi.li_state_qgs];
5267  mi.qgdd = currStaVector[mi.li_state_qgdd];
5268  mi.qgb = currStaVector[mi.li_state_qgb];
5269  // get the voltage drops, too
5270  vgs1 = currStaVector[mi.li_state_vgps];
5271  vbs1 = currStaVector[mi.li_state_vbs];
5272  vdds1 = currStaVector[mi.li_state_vdds];
5273 
5274  vgb1 = vgs1-vbs1;
5275  vgdd1 = vgs1-vdds1;
5276 
5277  // NOW we can calculate the charge update
5278  mi.qgs += mi.Capgs*(mi.vgps-vgs1);
5279  mi.qgdd += mi.Capgdd*(mi.vgpdd-vgdd1);
5280  mi.qgb += mi.Capgb*((mi.vgps-mi.vbs)-vgb1);
5281  }
5282 
5283  staVec[mi.li_state_qgs] = mi.qgs;
5284  staVec[mi.li_state_qgdd] = mi.qgdd;
5285  staVec[mi.li_state_qgb] = mi.qgb;
5286 
5287  // and the diode parasitic capacitors
5288  // these charges were set in updateIntermediateVars
5289  staVec[mi.li_state_qbd] = mi.qbd;
5290  staVec[mi.li_state_qbs] = mi.qbs;
5291 
5292  // diode #1 quatntities
5293 
5294  staVec[mi.li_state_D1DIOcapCharge] = mi.D1DIOcapCharge;
5295 
5296  // In the case of a charge, need to handle this:
5297  // Ensure dQ/dt = 0 for initial step after dcOP by setting the history equal
5298  // to the current value
5300  {
5301  currStaVector[mi.li_state_D1DIOcapCharge] = mi.D1DIOcapCharge;
5302  }
5303 
5304  staVec[ mi.li_state_von ] = mi.von;
5305 
5306  }
5307 // }
5308 
5309  return bsuccess;
5310 }
5311 
5312 //-----------------------------------------------------------------------------
5313 // Function : Master::updateSecondaryState
5314 // Purpose :
5315 // Special Notes :
5316 // Scope : public
5317 // Creator : Eric Keiter, SNL
5318 // Creation Date : 11/26/08
5319 //-----------------------------------------------------------------------------
5320 bool Master::updateSecondaryState ( double * staDerivVec, double * stoVec )
5321 {
5322  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
5323  {
5324  Instance & mi = *(*it);
5325 
5326  mi.cqgs = staDerivVec[mi.li_state_qgs];
5327  mi.cqgdd = staDerivVec[mi.li_state_qgdd];
5328  mi.cqgb = staDerivVec[mi.li_state_qgb];
5329  mi.cqbd = staDerivVec[mi.li_state_qbd];
5330  mi.cqbs = staDerivVec[mi.li_state_qbs];
5331 
5332  mi.D1DIOcapCurrent = staDerivVec[mi.li_state_D1DIOcapCharge];
5333 
5334  double Dtype =mi.getModel().dtype;
5335  double D1current = Dtype*(mi.D1cdeq + mi.D1DIOcapCurrent);
5336 
5337  // Even though these do not look exactly like the RHS loads, they do account
5338  // for the currents into the leads even with the lead resistors in place
5339  // --- since the lead resistances are the only thing between the primed node
5340  // and the lead, the currents are still correct AT CONVERGED SOLUTIONS
5341  /*
5342  Lead current load has been moved to F & Q loads
5343  mi.LeadCurrentD = -mi.getModel().dtype*(-mi.mode*mi.cdrain+mi.cbd+mi.cqgdd+mi.cqbd)-D1current-mi.Irdsshunt;
5344  mi.LeadCurrentG = mi.getModel().dtype*(mi.cqgs+mi.cqgb+mi.cqgdd);
5345  mi.LeadCurrentS = mi.getModel().dtype*(-mi.mode*mi.cdrain+mi.cbs-mi.cqgs-mi.cqbs) + D1current+mi.Irdsshunt;
5346  mi.LeadCurrentB = mi.getModel().dtype*(mi.cbd-mi.cbs-mi.cqgb+mi.cqbd+mi.cqbs);
5347  */
5348  }
5349 // }
5350 
5351  return true;
5352 }
5353 
5354 //-----------------------------------------------------------------------------
5355 // Function : Master::loadDAEVectors
5356 // Purpose :
5357 // Special Notes :
5358 // Scope : public
5359 // Creator : Eric Keiter, SNL
5360 // Creation Date : 11/26/08
5361 //-----------------------------------------------------------------------------
5362 bool Master::loadDAEVectors (double * solVec, double * fVec, double *qVec, double * bVec, double * storeLeadF, double * storeLeadQ, double * leadF, double * leadQ, double * junctionV)
5363 {
5364  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
5365  {
5366  Instance & mi = *(*it);
5367 
5368  // F-vector:
5369  double * dFdxdVp = mi.extData.dFdxdVpVectorRawPtr;
5370 
5371  double coef_Jdxp(0.0), gd_Jdxp(0.0);
5372  double gmin1 = getDeviceOptions().gmin;
5373  double Dtype =mi.getModel().dtype;
5374 
5375  double ceqbs = Dtype*(mi.cbs+mi.cqbs);
5376  double ceqbd = Dtype*(mi.cbd+mi.cqbd);
5377 
5378  double D1current = Dtype*mi.D1cdeq; // don't add in the capacitor stuff
5379 
5380 
5381  fVec[mi.li_Drain] += (mi.Idraindrift + mi.Irdsshunt - D1current);
5382 
5383  if (mi.Igate != 0.0)
5384  {
5385  fVec[mi.li_Gate] += mi.Igate;
5386  fVec[mi.li_GatePrime] += -mi.Igate;
5387  }
5388 
5389 
5390  fVec[mi.li_Source] += (mi.Isource - mi.Irdsshunt + mi.Ird1rs);
5391  fVec[mi.li_Bulk] += (ceqbs + ceqbd);
5392  fVec[mi.li_DrainPrime] += (-mi.Idrain-(ceqbd - mi.cdreq));
5393  fVec[mi.li_SourcePrime] += (-mi.Isource-(ceqbs + mi.cdreq));
5394  fVec[mi.li_DrainDrift] += (mi.Idrain - mi.Idraindrift);
5395  fVec[mi.li_D1Prime] += (D1current - mi.Ird1rs);
5396 
5397  // limiter terms:
5398  if (!mi.origFlag)
5399  {
5400  // bulk
5401  coef_Jdxp = Dtype*( + ((mi.gbd-gmin1))*(mi.vbdd-mi.vbdd_orig)
5402  + ((mi.gbs-gmin1))*(mi.vbs-mi.vbs_orig));
5403  dFdxdVp[mi.li_Bulk] += coef_Jdxp;
5404 
5405  // drain-prime
5406  coef_Jdxp = Dtype*(-((mi.gbd-gmin1))*
5407  (mi.vbdd-mi.vbdd_orig)+mi.gdds*(mi.vdds-mi.vdds_orig)
5408  +mi.Gm*((mi.mode>0)?(mi.vgps-mi.vgps_orig):(mi.vgpdd-mi.vgpdd_orig))
5409  +mi.Gmbs*((mi.mode>0)?(mi.vbs-mi.vbs_orig):(mi.vbdd-mi.vbdd_orig)));
5410  dFdxdVp[mi.li_DrainPrime] += coef_Jdxp;
5411 
5412  // source prime
5413  coef_Jdxp = Dtype*(-((mi.gbs-gmin1))*(mi.vbs-mi.vbs_orig)
5414  -mi.gdds*(mi.vdds-mi.vdds_orig)
5415  -mi.Gm*((mi.mode>0)?(mi.vgps-mi.vgps_orig):(mi.vgpdd-mi.vgpdd_orig))
5416  -mi.Gmbs*((mi.mode>0)?(mi.vbs-mi.vbs_orig):(mi.vbdd-mi.vbdd_orig)));
5417 
5418 
5419  dFdxdVp[mi.li_SourcePrime] += coef_Jdxp;
5420  gd_Jdxp = -mi.D1gd * (mi.D1vd-mi.D1vd_orig);
5421  dFdxdVp[mi.li_Drain] += gd_Jdxp;
5422  dFdxdVp[mi.li_D1Prime] -= gd_Jdxp;
5423  }
5424 
5425  // Q-vector:
5426  double * dQdxdVp = mi.extData.dQdxdVpVectorRawPtr;
5427 
5428  double qcoef_Jdxp(0.0);
5429 
5430  double Qeqbs = Dtype*mi.qbs;
5431  double Qeqbd = Dtype*mi.qbd;
5432  double Qeqgb = Dtype*mi.qgb;
5433  double Qeqgs = Dtype*mi.qgs;
5434  double Qeqgdd = Dtype*mi.qgdd;
5435 
5436 
5437  qVec[mi.li_Bulk] += (Qeqbs + Qeqbd - Qeqgb);
5438  qVec[mi.li_DrainPrime] += -(Qeqbd + Qeqgdd);
5439  qVec[mi.li_GatePrime] += (Qeqgs+Qeqgdd+Qeqgb);
5440  qVec[mi.li_SourcePrime] += (-(Qeqbs + Qeqgs));
5441  qVec[mi.li_D1Prime] += mi.D1DIOcapCharge;
5442  qVec[mi.li_Drain] += -mi.D1DIOcapCharge;
5443 
5444  // voltlim section:
5445  if (!mi.origFlag)
5446  {
5447  // bulk
5448  qcoef_Jdxp = Dtype*(-mi.Capgb*(mi.vgps-mi.vgps_orig-mi.vbs+mi.vbs_orig)
5449  + (+mi.Capgb)*(mi.vbdd-mi.vbdd_orig)
5450  + (+mi.capbs)*(mi.vbs-mi.vbs_orig));
5451  dQdxdVp[mi.li_Bulk] += qcoef_Jdxp;
5452 
5453  // drain-prime
5454  qcoef_Jdxp = Dtype*(-mi.Capgdd*(mi.vgpdd-mi.vgpdd_orig)-mi.capbd*(mi.vbdd-mi.vbdd_orig));
5455  dQdxdVp[mi.li_DrainPrime] += qcoef_Jdxp;
5456 
5457  // gate-prime
5458  qcoef_Jdxp = Dtype*(mi.Capgdd*(mi.vgpdd-mi.vgpdd_orig)+mi.Capgs*(mi.vgps-mi.vgps_orig)+
5459  mi.Capgb*(mi.vgps-mi.vgps_orig-mi.vbs+mi.vbs_orig));
5460  dQdxdVp[mi.li_GatePrime] += qcoef_Jdxp;
5461 
5462  // source-prime
5463  qcoef_Jdxp = Dtype*(-mi.Capgs*(mi.vgps-mi.vgps_orig)-mi.capbs*(mi.vbs-mi.vbs_orig));
5464  dQdxdVp[mi.li_SourcePrime] += qcoef_Jdxp;
5465  qcoef_Jdxp = -mi.D1capd*(mi.D1vd-mi.D1vd_orig);
5466  dQdxdVp[mi.li_D1Prime] -= qcoef_Jdxp;
5467  dQdxdVp[mi.li_Drain] += qcoef_Jdxp;
5468  }
5469 
5470  if( mi.loadLeadCurrent )
5471  {
5472  storeLeadF[mi.li_store_dev_id] = (mi.Idraindrift + mi.Irdsshunt - D1current);
5473  storeLeadF[mi.li_store_dev_is] = (mi.Isource - mi.Irdsshunt + mi.Ird1rs);
5474  storeLeadF[mi.li_store_dev_ig] = 0;
5475  storeLeadF[mi.li_store_dev_ib] = (ceqbs + ceqbd);
5476  // case where optional nodes become external nodes:
5477  if (mi.Igate != 0.0)
5478  {
5479  storeLeadF[mi.li_store_dev_ig] += mi.Igate;
5480  }
5481  if( !mi.gateCond )
5482  {
5483  // G' is G
5484  storeLeadF[mi.li_store_dev_ig] += -mi.Igate;
5485  }
5486 
5487  if( !mi.sourceCond )
5488  {
5489  // S' is S
5490  storeLeadF[mi.li_store_dev_is] += (-mi.Isource-(ceqbs + mi.cdreq));
5491  }
5492 
5493  if( !mi.drainCond )
5494  {
5495  // Ddrift is D'
5496  }
5497 
5498  if( !mi.model_.D1DIOconductance )
5499  {
5500  // D1' is S
5501  storeLeadF[mi.li_store_dev_is] += (D1current - mi.Ird1rs);
5502  }
5503 
5504  storeLeadQ[mi.li_store_dev_id] = -mi.D1DIOcapCharge;
5505  storeLeadQ[mi.li_store_dev_is] = 0;
5506  storeLeadQ[mi.li_store_dev_ig] = 0;
5507  storeLeadQ[mi.li_store_dev_ib] = (Qeqbs + Qeqbd - Qeqgb);
5508  // case where optional nodes become external nodes:
5509  if( !mi.gateCond )
5510  {
5511  // G' is G
5512  storeLeadQ[mi.li_store_dev_ig] += (Qeqgs+Qeqgdd+Qeqgb);
5513  }
5514 
5515  if( !mi.sourceCond )
5516  {
5517  // S' is S
5518  storeLeadQ[mi.li_store_dev_is] += (-(Qeqbs + Qeqgs));
5519  }
5520 
5521  if( !mi.drainCond )
5522  {
5523  // Ddrift is D'
5524  }
5525 
5526  if( !mi.model_.D1DIOconductance )
5527  {
5528  // D1' is S
5529  storeLeadQ[mi.li_store_dev_is] += mi.D1DIOcapCharge;
5530  }
5531  }
5532 
5533  }
5534 
5535  return true;
5536 }
5537 
5538 
5539 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
5540 //-----------------------------------------------------------------------------
5541 // Function : Master::loadDAEMatrices
5542 // Purpose :
5543 // Special Notes :
5544 // Scope : public
5545 // Creator : Eric Keiter, SNL
5546 // Creation Date : 11/26/08
5547 //-----------------------------------------------------------------------------
5548 bool Master::loadDAEMatrices (Linear::Matrix & dFdx, Linear::Matrix & dQdx)
5549 {
5550 #ifdef _OMP
5551 #pragma omp parallel for
5552 #endif
5553  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
5554  {
5555  Instance & mi = *(*it);
5556 
5557  // F-matrix:
5558 
5559  *mi.f_DrainEquDrainNodePtr += mi.dIdd_dVd + mi.gdsshunt + mi.D1gd;
5560 
5562 
5564 
5565  *mi.f_DrainEquD1PrimeNodePtr -= mi.D1gd;
5566 
5567  if (mi.gateCond != 0)
5568  {
5569 
5570  *mi.f_GateEquGateNodePtr += mi.gateCond;
5571 
5573  }
5574 
5575 
5577 
5578  *mi.f_SourceEquSourceNodePtr += mi.sourceCond + mi.gdsshunt + mi.D1gspr;
5579 
5580  *mi.f_SourceEquD1PrimeNodePtr += -mi.D1gspr;
5581 
5582  if (mi.sourceCond != 0.0)
5583  {
5584 
5586  }
5587 
5588 
5589  *mi.f_BulkEquBulkNodePtr += mi.gbs+mi.gbd;
5590 
5591  *mi.f_BulkEquDrainPrimeNodePtr -= mi.gbd;
5592 
5593 
5594  *mi.f_BulkEquSourcePrimeNodePtr -= mi.gbs;
5595 
5596 
5597  *mi.f_DrainPrimeEquBulkNodePtr += -mi.gbd+mi.Gmbs;
5598 
5600 
5602 
5604  if (mi.drainCond != 0.0)
5605  {
5606 
5608  }
5609 
5610  if (mi.gateCond != 0)
5611  {
5612 
5614 
5616  }
5617 
5618  if (mi.sourceCond != 0.0)
5619  {
5620 
5622  }
5623 
5624  *mi.f_SourcePrimeEquBulkNodePtr -= mi.gbs+mi.Gmbs;
5625 
5627 
5629 
5631 
5632 
5634 
5635  if (mi.drainCond != 0.0)
5636  {
5637 
5639  }
5640 
5642 
5643 
5644  *mi.f_D1PrimeEquDrainNodePtr -= mi.D1gd;
5645 
5646  *mi.f_D1PrimeEquSourceNodePtr += -mi.D1gspr;
5647 
5648  *mi.f_D1PrimeEquD1PrimeNodePtr += mi.D1gd + mi.D1gspr;
5649 
5650  // Q-matrix:
5651  if (!getSolverState().dcopFlag)
5652  {
5653 
5654  *mi.q_BulkEquBulkNodePtr += +mi.capbs+mi.capbd+mi.Capgb;
5655 
5656  *mi.q_BulkEquDrainPrimeNodePtr -= +mi.capbd;
5657 
5658  *mi.q_BulkEquGatePrimeNodePtr -= mi.Capgb;
5659 
5661 
5662 
5663  *mi.q_DrainPrimeEquBulkNodePtr += -mi.capbd;
5664 
5666 
5668 
5669 
5670  *mi.q_GatePrimeEquBulkNodePtr -= mi.Capgb;
5671 
5673 
5675 
5677 
5678 
5680 
5682 
5684 
5685 
5686  *mi.q_DrainEquDrainNodePtr += mi.D1capd;
5687 
5688  *mi.q_DrainEquD1PrimeNodePtr -= mi.D1capd;
5689 
5690 
5691  *mi.q_D1PrimeEquDrainNodePtr -= mi.D1capd;
5692 
5694  }
5695  }
5696  return true;
5697 }
5698 #else
5699 //-----------------------------------------------------------------------------
5700 // Function : Master::loadDAEMatrices
5701 // Purpose :
5702 // Special Notes :
5703 // Scope : public
5704 // Creator : Eric Keiter, SNL
5705 // Creation Date : 11/26/08
5706 //-----------------------------------------------------------------------------
5707 bool Master::loadDAEMatrices (Linear::Matrix & dFdx, Linear::Matrix & dQdx)
5708 {
5709  int sizeInstances = instanceContainer_.size();
5710 
5711 #ifdef _OMP
5712 #pragma omp parallel for
5713 #endif
5714  for (int i=0; i<sizeInstances; ++i)
5715  {
5716  Instance & mi = *(instanceContainer_.at(i));
5717 
5718  // F-matrix:
5719 
5720  dFdx[mi.li_Drain][mi.ADrainEquDrainNodeOffset] += mi.dIdd_dVd + mi.gdsshunt + mi.D1gd;
5721 
5722  dFdx[mi.li_Drain][mi.ADrainEquSourceNodeOffset] -= mi.gdsshunt;
5723 
5724  dFdx[mi.li_Drain][mi.ADrainEquDrainDriftNodeOffset] -= mi.dIdd_dVd;
5725 
5726  dFdx[mi.li_Drain][mi.ADrainEquD1PrimeNodeOffset] -= mi.D1gd;
5727 
5728  if (mi.gateCond != 0)
5729  {
5730 
5731  dFdx[mi.li_Gate][mi.AGateEquGateNodeOffset] += mi.gateCond;
5732 
5733  dFdx[mi.li_Gate][mi.AGateEquGatePrimeNodeOffset] -= mi.gateCond;
5734  }
5735 
5736 
5737  dFdx[mi.li_Source][mi.ASourceEquDrainNodeOffset] -= mi.gdsshunt;
5738 
5739  dFdx[mi.li_Source][mi.ASourceEquSourceNodeOffset] += mi.sourceCond + mi.gdsshunt + mi.D1gspr;
5740 
5741  dFdx[mi.li_Source][mi.ASourceEquD1PrimeNodeOffset] += -mi.D1gspr;
5742 
5743  if (mi.sourceCond != 0.0)
5744  {
5745 
5747  }
5748 
5749 
5750  dFdx[mi.li_Bulk][mi.ABulkEquBulkNodeOffset] += mi.gbs+mi.gbd;
5751 
5752  dFdx[mi.li_Bulk][mi.ABulkEquDrainPrimeNodeOffset] -= mi.gbd;
5753 
5754  dFdx[mi.li_Bulk][mi.ABulkEquSourcePrimeNodeOffset] -= mi.gbs;
5755 
5756 
5757  dFdx[mi.li_DrainPrime][mi.ADrainPrimeEquBulkNodeOffset] += -mi.gbd+mi.Gmbs;
5758 
5760 
5762 
5764  if (mi.drainCond != 0.0)
5765  {
5766 
5768  }
5769 
5770  if (mi.gateCond != 0)
5771  {
5772 
5774 
5776  }
5777 
5778  if (mi.sourceCond != 0.0)
5779  {
5780 
5782  }
5783 
5784  dFdx[mi.li_SourcePrime][mi.ASourcePrimeEquBulkNodeOffset] -= mi.gbs+mi.Gmbs;
5785 
5787 
5789 
5791 
5792 
5794 
5795  if (mi.drainCond != 0.0)
5796  {
5797 
5799  }
5800 
5802 
5803 
5804  dFdx[mi.li_D1Prime][mi.AD1PrimeEquDrainNodeOffset] -= mi.D1gd;
5805 
5806  dFdx[mi.li_D1Prime][mi.AD1PrimeEquSourceNodeOffset] += -mi.D1gspr;
5807 
5808  dFdx[mi.li_D1Prime][mi.AD1PrimeEquD1PrimeNodeOffset] += mi.D1gd + mi.D1gspr;
5809 
5810  // Q-matrix:
5811  if (!getSolverState().dcopFlag)
5812  {
5813 
5814  dQdx[mi.li_Bulk][mi.ABulkEquBulkNodeOffset] += +mi.capbs+mi.capbd+mi.Capgb;
5815 
5816  dQdx[mi.li_Bulk][mi.ABulkEquDrainPrimeNodeOffset] -= +mi.capbd;
5817 
5818  dQdx[mi.li_Bulk][mi.ABulkEquGatePrimeNodeOffset] -= mi.Capgb;
5819 
5820  dQdx[mi.li_Bulk][mi.ABulkEquSourcePrimeNodeOffset] -= +mi.capbs;
5821 
5822 
5823  dQdx[mi.li_DrainPrime][mi.ADrainPrimeEquBulkNodeOffset] += -mi.capbd;
5824 
5826 
5828 
5829 
5830  dQdx[mi.li_GatePrime][mi.AGatePrimeEquBulkNodeOffset] -= mi.Capgb;
5831 
5833 
5835 
5837 
5838 
5840 
5842 
5844 
5845 
5846  dQdx[mi.li_Drain][mi.ADrainEquDrainNodeOffset] += mi.D1capd;
5847 
5848  dQdx[mi.li_Drain][mi.ADrainEquD1PrimeNodeOffset] -= mi.D1capd;
5849 
5850 
5851  dQdx[mi.li_D1Prime][mi.AD1PrimeEquDrainNodeOffset] -= mi.D1capd;
5852 
5853  dQdx[mi.li_D1Prime][mi.AD1PrimeEquD1PrimeNodeOffset] += mi.D1capd;
5854  }
5855 
5856  }
5857  return true;
5858 }
5859 #endif
5860 
5861 Device *Traits::factory(const Configuration &configuration, const FactoryBlock &factory_block)
5862 {
5863 
5864  return new Master(configuration, factory_block, factory_block.solverState_, factory_block.deviceOptions_);
5865 }
5866 
5868 {
5870  .registerDevice("m", 18)
5871  .registerModelType("pmos", 18)
5872  .registerModelType("nmos", 18);
5873 }
5874 
5875 } // namespace VDMOS
5876 } // namespace Device
5877 } // namespace Xyce
const InstanceName & getName() const
#define CONSTPMOS
Definition: N_DEV_Const.h:77
bool UCCMqmeyer(double vgps, double vgpdd, double vgpb, double von_local, double vddsat_local, double &capgs_local, double &capgdd_local, double &capgb_local, double phi, double cox)
Definition: N_DEV_VDMOS.C:4170
static std::vector< std::vector< int > > jacStamp_D1C_DC_GC
Definition: N_DEV_VDMOS.h:651
static std::vector< std::vector< int > > jacStamp_D1C
Definition: N_DEV_VDMOS.h:657
static void initThermalModel(ParametricData< T > &parametric_data)
Add the parameter "TEMPMODEL" to the parametric_data.
static std::vector< std::vector< int > > jacMap2_DC_SC_GC
Definition: N_DEV_VDMOS.h:640
const DeviceOptions & deviceOptions_
Descriptor & addPar(const char *parName, T default_value, T U::*varPtr)
Adds the parameter description to the parameter map.
Definition: N_DEV_Pars.h:1429
double * q_SourcePrimeEquSourcePrimeNodePtr
Definition: N_DEV_VDMOS.h:522
static std::vector< int > jacMap_D1C_DC_SC
Definition: N_DEV_VDMOS.h:662
#define INVMC
static std::vector< std::vector< int > > jacMap2_GC
Definition: N_DEV_VDMOS.h:644
#define CONSTQ
Definition: N_DEV_Const.h:51
#define CONSTREFTEMP
Definition: N_DEV_Const.h:56
static std::vector< std::vector< int > > jacMap2_DC_GC
Definition: N_DEV_VDMOS.h:641
static std::vector< std::vector< int > > jacMap2_D1C
Definition: N_DEV_VDMOS.h:675
#define RD
double pnjlim(double vnew, double vold, double vt, double vcrit, int *icheck)
static std::vector< std::vector< int > > jacMap2_D1C_DC_SC
Definition: N_DEV_VDMOS.h:671
bool UCCMMeyercap(double vgps, double vgpdd, double vgpb, double &cgs, double &cgd, double &cgb)
Definition: N_DEV_VDMOS.C:4267
bool given(const std::string &parameter_name) const
Pure virtual class to augment a linear system.
Parameter may be specified as time dependent expression from netlist.
Definition: N_DEV_Pars.h:67
static std::vector< std::vector< int > > jacStamp_SC_GC
Definition: N_DEV_VDMOS.h:624
void addInternalNode(Util::SymbolTable &symbol_table, int index, const InstanceName &instance_name, const std::string &lead_name)
#define MC
#define VSIGMA
double * f_SourcePrimeEquSourcePrimeNodePtr
Definition: N_DEV_VDMOS.h:477
static std::vector< std::vector< int > > jacMap2_D1C_SC
Definition: N_DEV_VDMOS.h:673
#define CONSTMAX_EXP_ARG
Definition: N_DEV_Const.h:107
static std::vector< int > jacMap_D1C_SC
Definition: N_DEV_VDMOS.h:664
#define ISUBMOD
bool UCCMmosa2(double vgps, double vdds, double dvonvbs, double &cdraindrift_loc, double &vsate)
Definition: N_DEV_VDMOS.C:4826
static std::vector< int > jacMap
Definition: N_DEV_VDMOS.h:638
virtual bool updateState(double *solVec, double *staVec, double *stoVec)
Updates the devices state information.
Definition: N_DEV_VDMOS.C:5218
static std::vector< std::vector< int > > jacStamp_DC_SC_GC
Definition: N_DEV_VDMOS.h:622
InstanceVector::const_iterator getInstanceEnd() const
Returns an iterator to the ending of the vector of all instances created for this device...
static std::vector< int > jacMap_SC_GC
Definition: N_DEV_VDMOS.h:633
static void loadModelParameters(ParametricData< Model > &model_parameters)
Definition: N_DEV_VDMOS.C:128
static std::vector< std::vector< int > > jacStamp_D1C_GC
Definition: N_DEV_VDMOS.h:654
#define CONSTEPSSI
Definition: N_DEV_Const.h:80
#define ETA
#define AssertLIDs(cmp)
#define DELTA
static std::vector< std::vector< int > > jacStamp_DC
Definition: N_DEV_VDMOS.h:628
bool processParams()
processParams
Definition: N_DEV_VDMOS.C:1056
#define CONSTNMOS
Definition: N_DEV_Const.h:76
#define AI
bool updateTemperature(const double &temp_tmp)
Definition: N_DEV_VDMOS.C:3631
Parameter is subject to being set to minimum junction capacitance.
Definition: N_DEV_Pars.h:71
#define CONSTe
Definition: N_DEV_Const.h:60
static std::vector< std::vector< int > > jacStamp_D1C_SC_GC
Definition: N_DEV_VDMOS.h:652
virtual void registerJacLIDs(const JacobianStamp &jacLIDVec)
Parameter is subject to being set to minimum lead resistance.
Definition: N_DEV_Pars.h:70
static std::vector< std::vector< int > > jacStamp_D1C_DC
Definition: N_DEV_VDMOS.h:656
#define LS
double fetlim(double vnew, double vold, double vto)
#define BI
static std::vector< std::vector< int > > jacStamp_D1C_DC_SC
Definition: N_DEV_VDMOS.h:653
#define DELMAX
#define EXP_MAX
bool UCCMcvon(double vbs_local, double &von_local, double &dvonvbs_local)
Definition: N_DEV_VDMOS.C:4574
#define INVM
static std::vector< int > jacMap_D1C
Definition: N_DEV_VDMOS.h:666
static std::vector< int > jacMap_DC
Definition: N_DEV_VDMOS.h:637
static std::vector< std::vector< int > > jacMap2_SC
Definition: N_DEV_VDMOS.h:645
InstanceVector::const_iterator getInstanceBegin() const
Returns an iterator to the beginning of the vector of all instances created for this device...
static std::vector< std::vector< int > > jacMap2_D1C_DC_SC_GC
Definition: N_DEV_VDMOS.h:668
std::vector< Param > params
Parameters from the line.
#define THETA
static std::vector< std::vector< int > > jacMap2_SC_GC
Definition: N_DEV_VDMOS.h:642
void registerStoreLIDs(const std::vector< int > &stoLIDVecRef)
Definition: N_DEV_VDMOS.C:2006
void setParams(const std::vector< Param > &params)
const std::string & getName() const
static std::vector< std::vector< int > > jacMap2
Definition: N_DEV_VDMOS.h:647
double * q_SourcePrimeEquDrainPrimeNodePtr
Definition: N_DEV_VDMOS.h:520
The FactoryBlock contains parameters needed by the device, instance and model creation functions...
static std::vector< std::vector< int > > jacMap2_D1C_GC
Definition: N_DEV_VDMOS.h:672
static std::vector< std::vector< int > > jacMap2_D1C_DC
Definition: N_DEV_VDMOS.h:674
bool UCCMCharges(double vgps, double vgpdd, double vgpb, double &qD, double &qS, double &qB)
Definition: N_DEV_VDMOS.C:4435
const DeviceOptions & getDeviceOptions() const
static void loadInstanceParameters(ParametricData< Instance > &instance_parameters)
Definition: N_DEV_VDMOS.C:73
#define CONSTroot2
Definition: N_DEV_Const.h:50
static std::vector< std::vector< int > > jacStamp_D1C_DC_SC_GC
Definition: N_DEV_VDMOS.h:650
virtual bool loadDAEMatrices(Linear::Matrix &dFdx, Linear::Matrix &dQdx)
Populates the device's Jacobian object with these pointers.
Definition: N_DEV_VDMOS.C:5548
const std::vector< std::vector< int > > & jacobianStamp() const
Definition: N_DEV_VDMOS.C:2029
#define CONSTEPSOX
Definition: N_DEV_Const.h:79
virtual bool loadDAEVectors(double *solVec, double *fVec, double *qVec, double *bVec, double *storeLeadF, double *storeLeadQ, double *leadF, double *leadQ, double *junctionV)
Populates the device's ExternData object with these pointers.
Definition: N_DEV_VDMOS.C:5362
#define VMAX
static std::vector< std::vector< int > > jacStamp_GC
Definition: N_DEV_VDMOS.h:626
double limvds(double vnew, double vold)
Linear::Vector * nextStaVectorPtr
static std::vector< int > jacMap_DC_GC
Definition: N_DEV_VDMOS.h:632
virtual void forEachInstance(DeviceInstanceOp &op) const
Apply a device instance "op" to all instances associated with this model.
Definition: N_DEV_VDMOS.C:1041
static Config< T > & addConfiguration()
Adds the device to the Xyce device configuration.
static std::vector< std::vector< int > > jacMap2_DC_SC
Definition: N_DEV_VDMOS.h:643
Linear::Matrix * dFdxMatrixPtr
static std::vector< int > jacMap_GC
Definition: N_DEV_VDMOS.h:635
const DeviceOptions & getDeviceOptions() const
Returns the device options given during device construction.
The Device class is an interface for device implementations.
Definition: N_DEV_Device.h:101
static std::vector< std::vector< int > > jacMap2_D1C_SC_GC
Definition: N_DEV_VDMOS.h:670
#define SIGMA0
#define TOX
static std::vector< int > jacMap_D1C_DC_GC
Definition: N_DEV_VDMOS.h:660
void addStoreNode(Util::SymbolTable &symbol_table, int index, const InstanceName &instance_name, const std::string &lead_name)
static std::vector< int > jacMap_SC
Definition: N_DEV_VDMOS.h:636
const SolverState & solverState_
void registerJacLIDs(const std::vector< std::vector< int > > &jacLIDVec)
Definition: N_DEV_VDMOS.C:2059
Class Configuration contains device configuration data.
#define L
bool UCCMmosa1(double vgps, double vdds, double dvonvbs, double &cdraindrift_loc, double &vsate)
Definition: N_DEV_VDMOS.C:4615
double * f_DrainPrimeEquSourcePrimeNodePtr
Definition: N_DEV_VDMOS.h:464
#define M
void jacStampMap(const JacobianStamp &stamp_parent, IdVector &map_parent, JacobianStamp &map2_parent, JacobianStamp &stamp, IdVector &map, JacobianStamp &map2, int from, int to, int original_size)
const SolverState & getSolverState() const
double * f_SourcePrimeEquDrainPrimeNodePtr
Definition: N_DEV_VDMOS.h:475
static std::vector< std::vector< int > > jacStamp_D1C_SC
Definition: N_DEV_VDMOS.h:655
Linear::Vector * currStaVectorPtr
virtual std::ostream & printOutInstances(std::ostream &os) const
Definition: N_DEV_VDMOS.C:1006
void registerLIDs(const std::vector< int > &intLIDVecRef, const std::vector< int > &extLIDVecRef)
Definition: N_DEV_VDMOS.C:1811
static std::vector< int > jacMap_DC_SC_GC
Definition: N_DEV_VDMOS.h:631
#define Xyce_NONPOINTER_MATRIX_LOAD
Definition: N_DEV_Bsrc.C:92
static std::vector< int > jacMap_D1C_DC
Definition: N_DEV_VDMOS.h:665
const std::string & getType() const
static std::vector< int > jacMap_D1C_DC_SC_GC
Definition: N_DEV_VDMOS.h:659
static std::vector< int > jacMap_D1C_GC
Definition: N_DEV_VDMOS.h:663
static std::vector< std::vector< int > > jacStamp
Definition: N_DEV_VDMOS.h:629
#define RS
double * q_DrainPrimeEquSourcePrimeNodePtr
Definition: N_DEV_VDMOS.h:509
static std::vector< int > jacMap_D1C_SC_GC
Definition: N_DEV_VDMOS.h:661
#define CONSTKoverQ
Definition: N_DEV_Const.h:58
#define VSIGMAT
std::vector< Instance * > instanceContainer
Definition: N_DEV_VDMOS.h:722
#define W
ModelBlock represents a .MODEL line from the netlist.
Manages parameter binding for class C.
Definition: N_DEV_Pars.h:214
static std::vector< std::vector< int > > jacStamp_DC_GC
Definition: N_DEV_VDMOS.h:623
InstanceBlock represent a device instance line from the netlist.
static std::vector< std::vector< int > > jacStamp_DC_SC
Definition: N_DEV_VDMOS.h:625
std::vector< Param > params
Linear::Matrix * dQdxMatrixPtr
static Device * factory(const Configuration &configuration, const FactoryBlock &factory_block)
Definition: N_DEV_VDMOS.C:5861
static std::vector< std::vector< int > > jacMap2_DC
Definition: N_DEV_VDMOS.h:646
static std::vector< std::vector< int > > jacStamp_SC
Definition: N_DEV_VDMOS.h:627
#define MD
Linear::Vector * flagSolVectorPtr
#define INVMD
bool processInstanceParams()
processInstanceParams
Definition: N_DEV_VDMOS.C:1160
virtual bool updateSecondaryState(double *staDeriv, double *stoVec)
Updates the devices secondary state information.
Definition: N_DEV_VDMOS.C:5320
#define LAMBDA
static std::vector< int > jacMap_DC_SC
Definition: N_DEV_VDMOS.h:634
const SolverState & getSolverState() const
Returns the solver state given during device construction.
void setModParams(const std::vector< Param > &params)
static std::vector< std::vector< int > > jacMap2_D1C_DC_GC
Definition: N_DEV_VDMOS.h:669
void registerStateLIDs(const std::vector< int > &staLIDVecRef)
Definition: N_DEV_VDMOS.C:1921
Instance(const Configuration &configuration, const InstanceBlock &IB, Model &Miter, const FactoryBlock &factory_block)
Definition: N_DEV_VDMOS.C:1186
#define RSUB
#define CONSTboltz
Definition: N_DEV_Const.h:53
#define ALPHA
void loadNodeSymbols(Util::SymbolTable &symbol_table) const
Populates and returns the store name map.
Definition: N_DEV_VDMOS.C:1888