Xyce  6.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
N_DEV_MESFET.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-2014 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_MESFET.C,v $
27 //
28 // Purpose :
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.116.2.3 $
40 //
41 // Revision Date : $Date: 2014/03/06 23:33:43 $
42 //
43 // Current Owner : $Author: tvrusso $
44 //-------------------------------------------------------------------------
45 
46 #include <Xyce_config.h>
47 
48 
49 // ---------- Standard Includes ----------
50 #ifdef HAVE_CMATH
51 #include <cmath>
52 #else
53 #include <math.h>
54 #endif
55 
56 #ifdef HAVE_CSTDIO
57 #include <cstdio>
58 #else
59 #include <stdio.h>
60 #endif
61 
62 // ---------- Xyce Includes ----------
63 #include <N_DEV_Const.h>
64 #include <N_DEV_DeviceOptions.h>
65 #include <N_DEV_ExternData.h>
66 #include <N_DEV_MESFET.h>
67 #include <N_DEV_MatrixLoadData.h>
68 #include <N_DEV_SolverState.h>
69 #include <N_DEV_Message.h>
70 #include <N_ERH_ErrorMgr.h>
71 
72 #include <N_LAS_Matrix.h>
73 #include <N_LAS_Vector.h>
74 
75 namespace Xyce {
76 namespace Device {
77 namespace MESFET {
78 
80 {
81  p.addPar("TEMP", 0.0, &MESFET::Instance::temp)
82  .setExpressionAccess(ParameterType::TIME_DEP)
83  .setDescription("Device temperature");
84 
85  p.addPar("AREA", 1.0, &MESFET::Instance::area)
86  .setUnit(U_METER2)
87  .setCategory(CAT_GEOMETRY)
88  .setDescription("device area");
89 }
90 
92 {
93  p.addPar("AF", 1.0, &MESFET::Model::AF)
94  .setUnit(U_NONE)
95  .setCategory(CAT_FLICKER)
96  .setDescription("Flicker noise exponent");
97 
98  p.addPar("B", 0.3, &MESFET::Model::B)
99  .setUnit(U_VOLTM1)
100  .setCategory(CAT_PROCESS)
101  .setDescription("Doping tail parameter");
102 
103  p.addPar("BETA", 2.5e-3, &MESFET::Model::BETA)
104  .setUnit(U_AMPVM2)
105  .setCategory(CAT_PROCESS)
106  .setDescription("Transconductance parameter");
107 
108  p.addPar("ALPHA", 2.0, &MESFET::Model::ALPHA)
109  .setUnit(U_VOLTM1)
110  .setCategory(CAT_PROCESS)
111  .setDescription("Saturation voltage parameter");
112 
113  p.addPar("CGS", 0.0, &MESFET::Model::CGS)
114  .setExpressionAccess(ParameterType::MIN_CAP)
115  .setUnit(U_FARAD)
116  .setCategory(CAT_CAP)
117  .setDescription("Zero-bias gate-source junction capacitance");
118 
119  p.addPar("CGD", 0.0, &MESFET::Model::CGD)
120  .setExpressionAccess(ParameterType::MIN_CAP)
121  .setUnit(U_FARAD)
122  .setCategory(CAT_CAP)
123  .setDescription("Zero-bias gate-drain junction capacitance");
124 
125  p.addPar("FC", 0.5, &MESFET::Model::FC)
126  .setUnit(U_FARAD)
127  .setCategory(CAT_CAP)
128  .setDescription("Coefficient for forward-bias depletion capacitance");
129 
130  p.addPar("IS", 1e-14, &MESFET::Model::IS)
131  .setUnit(U_AMP)
132  .setCategory(CAT_CURRENT)
133  .setDescription("Gate junction saturation current");
134 
135  p.addPar("KF", 0.05, &MESFET::Model::KF)
136  .setUnit(U_NONE)
137  .setCategory(CAT_FLICKER)
138  .setDescription("Flicker noise coefficient");
139 
140  p.addPar("LAMBDA",0.0, &MESFET::Model::LAMBDA)
141  .setUnit(U_VOLTM1)
142  .setCategory(CAT_VOLT)
143  .setDescription("Channel length modulation");
144 
145  p.addPar("PB", 1.0, &MESFET::Model::PB)
146  .setUnit(U_VOLT)
147  .setCategory(CAT_VOLT)
148  .setDescription("Gate junction potential");
149 
150  p.addPar("RD", 0.0, &MESFET::Model::RD)
151  .setExpressionAccess(ParameterType::MIN_RES)
152  .setUnit(U_OHM)
153  .setCategory(CAT_RES)
154  .setDescription("Drain ohmic resistance");
155 
156  p.addPar("RS", 0.0, &MESFET::Model::RS)
157  .setExpressionAccess(ParameterType::MIN_RES)
158  .setUnit(U_OHM)
159  .setCategory(CAT_RES)
160  .setDescription("Source ohmic resistance");
161 
162  p.addPar("TNOM", 0.0, &MESFET::Model::TNOM);
163 
164  p.addPar("VTO", 0.0, &MESFET::Model::VTO)
165  .setUnit(U_VOLT)
166  .setCategory(CAT_VOLT)
167  .setDescription("Threshold voltage");
168 
170 }
171 
172 std::vector< std::vector<int> > Instance::jacStamp_DC_SC;
173 std::vector< std::vector<int> > Instance::jacStamp_DC;
174 std::vector< std::vector<int> > Instance::jacStamp_SC;
175 std::vector< std::vector<int> > Instance::jacStamp;
176 
177 std::vector<int> Instance::jacMap_DC_SC;
178 std::vector<int> Instance::jacMap_DC;
179 std::vector<int> Instance::jacMap_SC;
180 std::vector<int> Instance::jacMap;
181 
182 std::vector< std::vector<int> > Instance::jacMap2_DC_SC;
183 std::vector< std::vector<int> > Instance::jacMap2_DC;
184 std::vector< std::vector<int> > Instance::jacMap2_SC;
185 std::vector< std::vector<int> > Instance::jacMap2;
186 
187 //------------------- Class Model ---------------------------------
188 //-----------------------------------------------------------------------------
189 // Function : Model::Model
190 // Purpose : model block constructor
191 // Special Notes :
192 // Scope : public
193 // Creator : pmc
194 // Creation Date : 11/16/2003
195 //-----------------------------------------------------------------------------
197  const Configuration & configuration,
198  const ModelBlock & MB,
199  const FactoryBlock & factory_block)
200  : DeviceModel(MB, configuration.getModelParameters(), factory_block),
201  AF(1.0),
202  B(0.3),
203  ALPHA(2.0),
204  BETA(2.5e-3),
205  CGS(0.0),
206  CGD(0.0),
207  FC(0.5),
208  IS(1.0e-14),
209  KF(0.0),
210  LAMBDA(0.0),
211  PB(1.0),
212  RD(0.0),
213  RS(0.0),
214  TNOM(CONSTREFTEMP),
215  VTO(-2.0),
216  fNcoef(0.0),
217  fNexp(1.0),
218  dtype(CONSTNMOS)
219 {
220  if (getType() != "")
221  {
222  if (getType() == "NMF") {
223  dtype = CONSTNMOS;
224  }
225  else if (getType() == "PMF") {
226  dtype = CONSTPMOS;
227  }
228  else
229  {
230  UserError0(*this) << "Could not recognize the type for model " << getName();
231  }
232  }
233 
234 
235  // Set params to constant default values:
236  setDefaultParams ();
237 
238  // Set params according to .model line and constant defaults from metadata:
239  setModParams (MB.params);
240 
241  // Set any non-constant parameter defaults:
242  if (!given("TNOM"))
244 
245  // Calculate any parameters specified as expressions:
247 
248  // calculate dependent (ie computed) params and check for errors:
249  processParams ();
250 }
251 
252 //-----------------------------------------------------------------------------
253 // Function : Model::~Model
254 // Purpose : destructor
255 // Special Notes :
256 // Scope : public
257 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
258 // Creation Date : 3/16/00
259 //-----------------------------------------------------------------------------
261 {
262  std::vector<Instance*>::iterator iter;
263  std::vector<Instance*>::iterator first = instanceContainer.begin();
264  std::vector<Instance*>::iterator last = instanceContainer.end();
265 
266  for (iter=first; iter!=last; ++iter)
267  {
268  delete (*iter);
269  }
270 
271 }
272 
273 
274 //-----------------------------------------------------------------------------
275 // Function : Model::printOutInstances
276 // Purpose : debugging tool.
277 // Special Notes :
278 // Scope : public
279 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
280 // Creation Date : 3/16/00
281 //-----------------------------------------------------------------------------
282 std::ostream &Model::printOutInstances(std::ostream &os) const
283 {
284  std::vector<Instance*>::const_iterator iter;
285  std::vector<Instance*>::const_iterator first = instanceContainer.begin();
286  std::vector<Instance*>::const_iterator last = instanceContainer.end();
287 
288  int i,isize;
289  isize = instanceContainer.size();
290  os << std::endl;
291  os << "Number of MESFET Instances: " << isize << std::endl;
292  os << " name model name Parameters" << std::endl;
293 
294  for (i=0, iter=first; iter!=last; ++iter, ++i)
295  {
296  os << " " << i << ": " << (*iter)->getName() << "\t";
297  os << getName();
298  os << std::endl;
299  }
300 
301  os << std::endl;
302 
303  return os;
304 }
305 
306 //-----------------------------------------------------------------------------
307 // Function : Model::forEachInstance
308 // Purpose :
309 // Special Notes :
310 // Scope : public
311 // Creator : David Baur
312 // Creation Date : 2/4/2014
313 //-----------------------------------------------------------------------------
314 /// Apply a device instance "op" to all instances associated with this
315 /// model
316 ///
317 /// @param[in] op Operator to apply to all instances.
318 ///
319 ///
320 void Model::forEachInstance(DeviceInstanceOp &op) const /* override */
321 {
322  for (std::vector<Instance *>::const_iterator it = instanceContainer.begin(); it != instanceContainer.end(); ++it)
323  op(*it);
324 }
325 
326 
327 //-----------------------------------------------------------------------------
328 // Function : Model::processParams
329 // Purpose :
330 // Special Notes :
331 // Scope : public
332 // Creator : pmc
333 // Creation Date : 11/16/2003
334 //-----------------------------------------------------------------------------
336 {
337 
338  return true;
339 }
340 
341 //----------------------------------------------------------------------------
342 // Function : Model::processInstanceParams
343 // Purpose :
344 // Special Notes :
345 // Scope : public
346 // Creator : Dave Shirely, PSSI
347 // Creation Date : 03/23/06
348 //----------------------------------------------------------------------------
350 {
351 
352  std::vector<Instance*>::iterator iter;
353  std::vector<Instance*>::iterator first = instanceContainer.begin();
354  std::vector<Instance*>::iterator last = instanceContainer.end();
355 
356  for (iter=first; iter!=last; ++iter)
357  {
358  (*iter)->processParams();
359  }
360 
361  return true;
362 }
363 
364 //------------------------ Class Instance -------------------------
365 //-----------------------------------------------------------------------------
366 // Function : Instance::Instance
367 // Purpose : instance block constructor
368 // Special Notes :
369 // Scope : public
370 // Creator : pmc
371 // Creation Date : 11/16/2003
372 //-----------------------------------------------------------------------------
374  const Configuration & configuration,
375  const InstanceBlock & instance_block,
376  Model & model,
377  const FactoryBlock & factory_block)
378  : DeviceInstance(instance_block, configuration.getInstanceParameters(), factory_block),
379  model_(model),
380  limitedFlag(false),
381  off(0),
382  ic(0),
383  area(1.0),
384  ic_vds(0.0),
385  ic_vgs(0.0),
386  temp(getDeviceOptions().temp.getImmutableValue<double>()),
387  sourceCond(0.0),
388  drainCond(0.0),
389  tCGS(0.0),
390  tCGD(0.0),
391  tIS(0.0),
392  tPB(0.0),
393  tBeta(0.0),
394  tvt0(0.0),
395  tLambda(0.0),
396  tAlpha(0.0),
397  tRD(0.0),
398  tRS(0.0),
399  tMESb(0.0),
400  Bfac(0.0),
401  dNode(0),
402  gNode(0),
403  sNode(0),
404  dpNode(0),
405  spNode(0),
406  Vgs(0.0),
407  Vgd(0.0),
408  gm(0.0),
409  gds(0.0),
410  ggs(0.0),
411  ggd(0.0),
412  p(0.0),
413  // Solution variables and intermediate quantities
414  // drain,source,gate, drainprime and sourceprime voltages
415  Vd(0.0),
416  Vs(0.0),
417  Vg(0.0),
418  Vdp(0.0),
419  Vsp(0.0),
420  // vector local indices
421  li_Drain(-1),
422  li_DrainPrime(-1),
423  li_Source(-1),
424  li_SourcePrime(-1),
425  li_Gate(-1),
426  // Jacobian Matrix
427  // Jacobian Matrix Offset:
428  // V_d Row:
429  ADrainEquDrainNodeOffset(-1),
430  ADrainEquDrainPrimeNodeOffset(-1),
431  // V_g Row:
432  AGateEquGateNodeOffset(-1),
433  AGateEquDrainPrimeNodeOffset(-1),
434  AGateEquSourcePrimeNodeOffset(-1),
435  // V_s Row:
436  ASourceEquSourceNodeOffset(-1),
437  ASourceEquSourcePrimeNodeOffset(-1),
438  // V_d' Row:
439  ADrainPrimeEquDrainNodeOffset(-1),
440  ADrainPrimeEquGateNodeOffset(-1),
441  ADrainPrimeEquDrainPrimeNodeOffset(-1),
442  ADrainPrimeEquSourcePrimeNodeOffset(-1),
443  // V_s' Row:
444  ASourcePrimeEquGateNodeOffset(-1),
445  ASourcePrimeEquSourceNodeOffset(-1),
446  ASourcePrimeEquDrainPrimeNodeOffset(-1),
447  ASourcePrimeEquSourcePrimeNodeOffset(-1),
448 
450  // dFdx Matrix Ptr:
451  // V_d Row:
452  f_DrainEquDrainNodePtr(0),
453  f_DrainEquDrainPrimeNodePtr(0),
454  // V_g Row:
455  f_GateEquGateNodePtr(0),
456  f_GateEquDrainPrimeNodePtr(0),
457  f_GateEquSourcePrimeNodePtr(0),
458  // V_s Row:
459  f_SourceEquSourceNodePtr(0),
460  f_SourceEquSourcePrimeNodePtr(0),
461  // V_d' Row:
462  f_DrainPrimeEquDrainNodePtr(0),
463  f_DrainPrimeEquGateNodePtr(0),
464  f_DrainPrimeEquDrainPrimeNodePtr(0),
465  f_DrainPrimeEquSourcePrimeNodePtr(0),
466  // V_s' Row:
467  f_SourcePrimeEquGateNodePtr(0),
468  f_SourcePrimeEquSourceNodePtr(0),
469  f_SourcePrimeEquDrainPrimeNodePtr(0),
470  f_SourcePrimeEquSourcePrimeNodePtr(0),
471 
472  // dQdx Matrix Ptr:
473  // V_d Row:
474  q_DrainEquDrainNodePtr(0),
475  q_DrainEquDrainPrimeNodePtr(0),
476  // V_g Row:
477  q_GateEquGateNodePtr(0),
478  q_GateEquDrainPrimeNodePtr(0),
479  q_GateEquSourcePrimeNodePtr(0),
480  // V_s Row:
481  q_SourceEquSourceNodePtr(0),
482  q_SourceEquSourcePrimeNodePtr(0),
483  // V_d' Row:
484  q_DrainPrimeEquDrainNodePtr(0),
485  q_DrainPrimeEquGateNodePtr(0),
486  q_DrainPrimeEquDrainPrimeNodePtr(0),
487  q_DrainPrimeEquSourcePrimeNodePtr(0),
488  // V_s' Row:
489  q_SourcePrimeEquGateNodePtr(0),
490  q_SourcePrimeEquSourceNodePtr(0),
491  q_SourcePrimeEquDrainPrimeNodePtr(0),
492  q_SourcePrimeEquSourcePrimeNodePtr(0),
493 #endif
494  vgs(0.0),
495  vgd(0.0),
496  vgs_old(0.0),
497  vgd_old(0.0),
498  vds_old(0.0),
499  vgs_orig(0.0),
500  vgd_orig(0.0),
501 
502  capgs(0.0),
503  qgs(0.0),
504  cqgs(0.0),
505  capgd(0.0),
506  qgd(0.0),
507  cqgd(0.0),
508  mode(1),
509  // local indices
510  li_store_vgs(-1),
511  li_store_vgd(-1),
512  li_store_dev_id(-1),
513  li_store_dev_ig(-1),
514  li_store_dev_is(-1),
515  li_state_qgs(-1),
516  li_state_gcgs(-1),
517  li_state_qgd(-1),
518  li_state_gcgd(-1)
519 {
520  numIntVars = 2;
521  numExtVars = 3;
522  numStateVars = 4;
523  setNumStoreVars(2);
524  numLeadCurrentStoreVars = 3; // lead currents drain, gate and source
525 
526  devConMap.resize(3);
527  devConMap[0] = 1;
528  devConMap[1] = 2;
529  devConMap[2] = 1;
530 
531  if( jacStamp.empty() )
532  {
533  // stamp for RS!=0, RD!=0
534  jacStamp_DC_SC.resize(5);
535  jacStamp_DC_SC[0].resize(2); // Drain row
536  jacStamp_DC_SC[0][0]=0; // d-d
537  jacStamp_DC_SC[0][1]=3; // d-d'
538  jacStamp_DC_SC[1].resize(3); // Gate row
539  jacStamp_DC_SC[1][0]=1; // g-g
540  jacStamp_DC_SC[1][1]=3; // g-d'
541  jacStamp_DC_SC[1][2]=4; // g-s'
542  jacStamp_DC_SC[2].resize(2); // Source row
543  jacStamp_DC_SC[2][0]=2; // s-s
544  jacStamp_DC_SC[2][1]=4; // s-s'
545  jacStamp_DC_SC[3].resize(4); // Drain' row
546  jacStamp_DC_SC[3][0]=0; // d'-d
547  jacStamp_DC_SC[3][1]=1; // d'-g
548  jacStamp_DC_SC[3][2]=3; // d'-d'
549  jacStamp_DC_SC[3][3]=4; // d'-s'
550  jacStamp_DC_SC[4].resize(4); // Source' row
551  jacStamp_DC_SC[4][0]=1; // s'-g
552  jacStamp_DC_SC[4][1]=2; // s'-s
553  jacStamp_DC_SC[4][2]=3; // s'-d'
554  jacStamp_DC_SC[4][3]=4; // s'-s'
555 
556  jacMap_DC_SC.clear();
558  jacStamp_DC, jacMap_DC, jacMap2_DC, 4, 2, 5);
559 
561  jacStamp_SC, jacMap_SC, jacMap2_SC, 3, 0, 5);
562 
564  jacStamp, jacMap, jacMap2, 3, 0, 5);
565 
566  }
567 
568  // Set params to constant default values:
569  setDefaultParams ();
570 
571  // Set params according to instance line and constant defaults from metadata:
572  setParams(instance_block.params);
573 
574  // Set any non-constant parameter defaults:
575  if (!given("TEMP"))
576  temp = getDeviceOptions().temp.getImmutableValue<double>();
577 
579 
580  // Calculate any parameters specified as expressions:
581  processParams ();
582 
583  // process source/drain series resistance
584  drainCond = 0;
585  if (model_.RD != 0)
587  sourceCond = 0;
588  if (model_.RS != 0)
590 
591  numIntVars = (((sourceCond == 0.0)?0:1)+((drainCond == 0.0)?0:1));
592 
593 }
594 
595 //-----------------------------------------------------------------------------
596 // Function : Instance::~Instance
597 // Purpose : destructor
598 // Special Notes :
599 // Scope : public
600 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
601 // Creation Date : 3/16/00
602 //-----------------------------------------------------------------------------
604 {
605 }
606 
607 //----------------------------------------------------------------------------
608 // Function : Instance::registerLIDs
609 // Purpose :
610 // Special Notes :
611 // Scope : public
612 // Creator : pmc
613 // Creation Date : 11/16/2003
614 //----------------------------------------------------------------------------
615 void Instance::registerLIDs( const std::vector<int> & intLIDVecRef,
616  const std::vector<int> & extLIDVecRef )
617 {
618  numIntVars = (((sourceCond == 0.0)?0:1)+((drainCond == 0.0)?0:1));
619 
620  AssertLIDs(intLIDVecRef.size() == numIntVars);
621  AssertLIDs(extLIDVecRef.size() == numExtVars);
622 
623 #ifdef Xyce_DEBUG_DEVICE
624  if (getDeviceOptions().debugLevel > 0)
625  {
626  Xyce::dout() << std::endl << section_divider << std::endl;
627  Xyce::dout() << " Instance::registerLIDs" << std::endl;
628  Xyce::dout() << " name = " << getName() << std::endl;
629  Xyce::dout() << " number of internal variables: " << numIntVars << std::endl;
630  Xyce::dout() << " number of external variables: " << numExtVars << std::endl;
631  }
632 #endif
633 
634  // copy over the global ID lists.
635  intLIDVec = intLIDVecRef;
636  extLIDVec = extLIDVecRef;
637 
638  // now use these lists to obtain the indices into the
639  // linear algebra entities. This assumes an order.
640  // For the matrix indices, first do the rows.
641 
642  li_Drain = extLIDVec[0];
643  li_Gate = extLIDVec[1];
644  li_Source = extLIDVec[2];
645 
646  int intLoc = 0;
647 
648  if( drainCond )
649  li_DrainPrime = intLIDVec[intLoc++];
650  else
652 
653  if( sourceCond )
654  li_SourcePrime = intLIDVec[intLoc];
655  else
657 
658 #ifdef Xyce_DEBUG_DEVICE
659  if (getDeviceOptions().debugLevel > 0)
660  {
661  Xyce::dout() << "\n variable local indices:\n";
662  Xyce::dout() << " li_Drain = " << li_Drain << std::endl;
663  Xyce::dout() << " li_DrainPrime = " << li_DrainPrime << std::endl;
664  Xyce::dout() << " li_Source = " << li_Source << std::endl;
665  Xyce::dout() << " li_SourcePrime = " << li_SourcePrime << std::endl;
666  Xyce::dout() << " li_Gate = " << li_Gate << std::endl;
667 
668  Xyce::dout() << section_divider << std::endl;
669  }
670 #endif
671 }
672 
673 //-----------------------------------------------------------------------------
674 // Function : Instance::getIntNameMap
675 // Purpose :
676 // Special Notes :
677 // Scope : public
678 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
679 // Creation Date : 05/13/05
680 //-----------------------------------------------------------------------------
681 std::map<int,std::string> & Instance::getIntNameMap ()
682 {
683  // set up the internal name map, if it hasn't been already.
684  if (intNameMap.empty ())
685  {
686  // set up the internal name map:
687  std::string tmpstr;
688 
689  if (drainCond != 0.0)
690  {
691  tmpstr = getName()+"_drainprime";
692  spiceInternalName (tmpstr);
693  intNameMap[ li_DrainPrime ] = tmpstr;
694  }
695 
696  if (sourceCond != 0.0)
697  {
698  tmpstr = getName()+"_sourceprime";
699  spiceInternalName (tmpstr);
700  intNameMap[ li_SourcePrime ] = tmpstr;
701  }
702  }
703 
704  return intNameMap;
705 }
706 
707 
708 //-----------------------------------------------------------------------------
709 // Function : N_DEV_MESFETInstance::getStoreNameMap
710 // Purpose :
711 // Special Notes :
712 // Scope : public
713 // Creator : Richard Schiek, Electrical Systems Modeling
714 // Creation Date : 4/4/2013
715 //-----------------------------------------------------------------------------
716 std::map<int,std::string> & N_DEV_MESFETInstance::getStoreNameMap ()
717 {
718  // set up the internal name map, if it hasn't been already.
719  if( loadLeadCurrent && storeNameMap.empty ())
720  {
721  // change subcircuitname:devicetype_deviceName to
722  // devicetype:subcircuitName:deviceName
723  std::string modName(getName());
724  spiceInternalName(modName);
725  std::string tmpstr;
726  tmpstr = modName+":DEV_ID";
727  storeNameMap[ li_store_dev_id ] = tmpstr;
728  tmpstr = modName+":DEV_IG";
729  storeNameMap[ li_store_dev_ig ] = tmpstr;
730  tmpstr = modName+":DEV_IS";
731  storeNameMap[ li_store_dev_is ] = tmpstr;
732  }
733 
734  return storeNameMap;
735 }
736 
737 
738 //----------------------------------------------------------------------------
739 // Function : Instance::registerStateLIDs
740 // Purpose :
741 // Special Notes :
742 // Scope : public
743 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
744 // Creation Date : 3/16/00
745 //----------------------------------------------------------------------------
746 void Instance::registerStateLIDs(const std::vector<int> & staLIDVecRef)
747 {
748  AssertLIDs(staLIDVecRef.size() == numStateVars);
749 
750 #ifdef Xyce_DEBUG_DEVICE
751  if (getDeviceOptions().debugLevel > 0)
752  {
753  Xyce::dout() << std::endl;
754  Xyce::dout() << section_divider << std::endl;
755  Xyce::dout() << " In Instance::registerStateLIDs\n\n";
756  Xyce::dout() << " name = " << getName() << std::endl;
757  Xyce::dout() << " Number of State LIDs: " << numStateVars << std::endl;
758  }
759 #endif
760 
761  // Copy over the global ID lists:
762  staLIDVec = staLIDVecRef;
763 
764  int lid=0;
765  li_state_qgs = staLIDVec[lid++];
766  li_state_gcgs = staLIDVec[lid++];
767 
768  li_state_qgd = staLIDVec[lid++];
769  li_state_gcgd = staLIDVec[lid++];
770 
771 #ifdef Xyce_DEBUG_DEVICE
772  if (getDeviceOptions().debugLevel > 0)
773  {
774  Xyce::dout() << " State local indices:" << std::endl;
775  Xyce::dout() << std::endl;
776 
777  Xyce::dout() << " li_state_qgs = " << li_state_qgs << std::endl;
778  Xyce::dout() << " li_state_gcgs = " << li_state_gcgs;
779  Xyce::dout() << " li_state_qgd = " << li_state_qgd;
780  Xyce::dout() << " li_state_gcgd = " << li_state_gcgd << std::endl;;
781 
782  Xyce::dout() << section_divider << std::endl;
783  }
784 #endif
785 
786 }
787 
788 //----------------------------------------------------------------------------
789 // Function : Instance::registerStoreLIDs
790 // Purpose :
791 // Special Notes :
792 // Scope : public
793 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
794 // Creation Date : 12/9/11
795 //----------------------------------------------------------------------------
796 void Instance::registerStoreLIDs(const std::vector<int> & stoLIDVecRef)
797 {
798  AssertLIDs(stoLIDVecRef.size() == getNumStoreVars());
799 
800  // Copy over the global ID lists:
801  stoLIDVec = stoLIDVecRef;
802 
803  int lid=0;
804  li_store_vgs = stoLIDVec[lid++];
805  li_store_vgd = stoLIDVec[lid++];
806 
807  if( loadLeadCurrent )
808  {
809  li_store_dev_id = stoLIDVec[lid++];
810  li_store_dev_ig = stoLIDVec[lid++];
811  li_store_dev_is = stoLIDVec[lid++];
812  }
813 }
814 
815 //----------------------------------------------------------------------------
816 // Function : Instance::jacobianStamp
817 // Purpose :
818 // Special Notes :
819 // Scope : public
820 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
821 // Creation Date : 3/16/00
822 //----------------------------------------------------------------------------
823 const std::vector< std::vector<int> > & Instance::jacobianStamp() const
824 {
825  if( drainCond != 0.0 && sourceCond != 0.0 )
826  return jacStamp_DC_SC;
827  else if( drainCond != 0.0 && sourceCond == 0.0 )
828  return jacStamp_DC;
829  else if( drainCond == 0.0 && sourceCond != 0.0 )
830  return jacStamp_SC;
831  else if( drainCond == 0.0 && sourceCond == 0.0 )
832  return jacStamp;
833  else
834  return jacStamp;
835 }
836 
837 //----------------------------------------------------------------------------
838 // Function : Instance::registerJacLIDs
839 // Special Notes :
840 // Scope : public
841 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
842 // Creation Date : 3/16/00
843 //----------------------------------------------------------------------------
844 void Instance::registerJacLIDs( const std::vector< std::vector<int> > & jacLIDVec )
845 {
846  DeviceInstance::registerJacLIDs( jacLIDVec );
847  std::vector<int> map;
848  std::vector< std::vector<int> > map2;
849 
850  if (drainCond != 0.0)
851  {
852  if (sourceCond != 0.0)
853  {
854  map = jacMap_DC_SC;
855  map2 = jacMap2_DC_SC;
856  }
857  else
858  {
859  map = jacMap_DC;
860  map2 = jacMap2_DC;
861  }
862  }
863  else
864  {
865  if (sourceCond != 0.0)
866  {
867  map = jacMap_SC;
868  map2 = jacMap2_SC;
869  }
870  else
871  {
872  map = jacMap;
873  map2 = jacMap2;
874  }
875  }
876 
877  ADrainEquDrainNodeOffset = jacLIDVec[map[0]][map2[0][0]];
878  ADrainEquDrainPrimeNodeOffset = jacLIDVec[map[0]][map2[0][1]];
879 
880  AGateEquGateNodeOffset = jacLIDVec[map[1]][map2[1][0]];
881  AGateEquDrainPrimeNodeOffset = jacLIDVec[map[1]][map2[1][1]];
882  AGateEquSourcePrimeNodeOffset = jacLIDVec[map[1]][map2[1][2]];
883 
884  ASourceEquSourceNodeOffset = jacLIDVec[map[2]][map2[2][0]];
885  ASourceEquSourcePrimeNodeOffset = jacLIDVec[map[2]][map2[2][1]];
886 
887  ADrainPrimeEquDrainNodeOffset = jacLIDVec[map[3]][map2[3][0]];
888  ADrainPrimeEquGateNodeOffset = jacLIDVec[map[3]][map2[3][1]];
889  ADrainPrimeEquDrainPrimeNodeOffset = jacLIDVec[map[3]][map2[3][2]];
890  ADrainPrimeEquSourcePrimeNodeOffset = jacLIDVec[map[3]][map2[3][3]];
891 
892  ASourcePrimeEquGateNodeOffset = jacLIDVec[map[4]][map2[4][0]];
893  ASourcePrimeEquSourceNodeOffset = jacLIDVec[map[4]][map2[4][1]];
894  ASourcePrimeEquDrainPrimeNodeOffset = jacLIDVec[map[4]][map2[4][2]];
895  ASourcePrimeEquSourcePrimeNodeOffset = jacLIDVec[map[4]][map2[4][3]];
896 }
897 
898 //-----------------------------------------------------------------------------
899 // Function : Instance::setupPointers
900 // Purpose :
901 // Special Notes :
902 // Scope : public
903 // Creator : Eric Keiter, SNL
904 // Creation Date : 11/30/08
905 //-----------------------------------------------------------------------------
907 {
908 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
909  N_LAS_Matrix & dFdx = *(extData.dFdxMatrixPtr);
910  N_LAS_Matrix & dQdx = *(extData.dQdxMatrixPtr);
911 
912  // F-matrix:
915 
919 
922 
927 
932 
933  // Q-matrix:
936 
940 
943 
948 
953 
954 #endif
955 }
956 
957 //----------------------------------------------------------------------------
958 // Function : Instance::updatePrimaryState
959 // Purpose :
960 // Special Notes :
961 // Scope : public
962 // Creator : pmc
963 // Creation Date : 11/16/2003
964 //----------------------------------------------------------------------------
966 {
967  double * staVec = extData.nextStaVectorRawPtr;
968  bool bsuccess = updateIntermediateVars ();
969 
970  double * stoVec = extData.nextStoVectorRawPtr;
971  stoVec[li_store_vgs] = vgs;
972  stoVec[li_store_vgd] = vgd;
973  staVec[li_state_qgs] = qgs;
974  staVec[li_state_qgd] = qgd;
975 
976  return bsuccess;
977 }
978 
979 //-----------------------------------------------------------------------------
980 // Function : Instance::updateIntermediateVars
981 // Purpose :
982 // Special Notes :
983 // Scope : public
984 // Creator : pmc
985 // Creation Date : 11/16/2003
986 //-----------------------------------------------------------------------------
988 {
989  double * solVec = extData.nextSolVectorRawPtr;
990  double * currStaVec = extData.currStaVectorRawPtr;
991 
992  int dtype;
993  double csat, betap;
994  double vgst, vgdt;
995  double evgs, evgd;
996  double sarg, vtf;;
997 
998 // from the spice jfet
999  double czgd, czgs;
1000  double czgdf2, czgsf2;
1001  double fcpb2;
1002  double twop;
1003  int icheck, ichk1;
1004 
1005 // for the Shockley version
1006 // double A, B, C, B12, C12, D, Vdsat;
1007 // double delta;
1008  double prod, denom, invdenom, afact, lfact;
1009 
1010 #ifdef Xyce_DEBUG_DEVICE
1011  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1012  {
1013  Xyce::dout() << subsection_divider << std::endl;
1014  Xyce::dout() <<" Instance::updateIntermediateVars.\n"<<std::endl;
1015  Xyce::dout() <<" name = " << getName() << std::endl;
1016  Xyce::dout() <<" Model name = " << model_.getName() << std::endl;
1017  Xyce::dout() <<" dtype is " << model_.dtype << std::endl;
1018  Xyce::dout() << std::endl;
1019  Xyce::dout().width(25); Xyce::dout().precision(17); Xyce::dout().setf(std::ios::scientific);
1020  }
1021 #endif
1022 
1023  icheck = 1;
1024  dtype = model_.dtype;
1025 
1026  // we need our solution variables for any of this stuff
1027  Vd = 0.0;
1028  Vs = 0.0;
1029  Vg = 0.0;
1030  Vdp = 0.0;
1031  Vsp = 0.0;
1032 
1033  Vd = solVec[li_Drain];
1034  Vg = solVec[li_Gate];
1035  Vs = solVec[li_Source];
1036  Vsp = solVec[li_SourcePrime];
1037  Vdp = solVec[li_DrainPrime];
1038 
1039 #ifdef Xyce_DEBUG_DEVICE
1040  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1041  {
1042  Xyce::dout() << " " << std::endl;
1043  Xyce::dout() << " Vg = " << Vg << std::endl;
1044  Xyce::dout() << " Vd = " << Vd << std::endl;
1045  Xyce::dout() << " Vs = " << Vs << std::endl;
1046  Xyce::dout() << " Vdp = " << Vdp << std::endl;
1047  Xyce::dout() << " Vsp = " << Vsp << std::endl;
1048  }
1049 #endif
1050 
1051  // now we need voltage drops
1052  Vddp = Vd - Vdp;
1053  Vssp = Vs - Vsp;
1054  Vgsp = Vg - Vsp;
1055  Vgdp = Vg - Vdp;
1056  Vdpsp = Vdp - Vsp;
1057 
1058  // Now the things that the 3f5 code really uses
1059  vgs = dtype * Vgsp;
1060  vgd = dtype * Vgdp;
1061  vds = vgs-vgd;
1062 
1063  origFlag = 1;
1064  limitedFlag = false;
1065  vgs_orig = vgs;
1066  vgd_orig = vgd;
1067  vds_orig = vds;
1068 
1069  if (getSolverState().newtonIter == 0)
1070  {
1071  newtonIterOld=0;
1072  if (getSolverState().initJctFlag && getDeviceOptions().voltageLimiterFlag)
1073  {
1074  if (getSolverState().inputOPFlag)
1075  {
1076  N_LAS_Vector * flagSolVectorPtr = extData.flagSolVectorPtr;
1077  if ((*flagSolVectorPtr)[li_Drain] == 0 || (*flagSolVectorPtr)[li_Gate] == 0 ||
1078  (*flagSolVectorPtr)[li_Source] == 0 || (*flagSolVectorPtr)[li_SourcePrime] ||
1079  (*flagSolVectorPtr)[li_DrainPrime] )
1080  {
1081  vgs = 0;
1082  vgd = 0;
1083  vds = vgs-vgd;
1084  }
1085  }
1086  else
1087  {
1088  vgs = 0;
1089  vgd = 0;
1090  vds = vgs-vgd;
1091  }
1092  }
1093  if (!(getSolverState().dcopFlag)||(getSolverState().locaEnabledFlag && getSolverState().dcopFlag))
1094  {
1095  double * currStoVec = extData.currStoVectorRawPtr;
1096  vgs_old = currStoVec[li_store_vgs];
1097  vgd_old = currStoVec[li_store_vgd];
1098  }
1099  else
1100  { // there is no history
1101  vgs_old = vgs;
1102  vgd_old = vgd;
1103  }
1104  }
1105  else
1106  {
1107  double *stoVec = extData.nextStoVectorRawPtr;
1108  vgs_old = stoVec[li_store_vgs];
1109  vgd_old = stoVec[li_store_vgd];
1110  }
1111 
1112  // SPICE-type Voltage Limiting
1113  ///////////////////////////////
1114 
1115  if (getDeviceOptions().voltageLimiterFlag)
1116  {
1117 #ifdef Xyce_DEBUG_DEVICE
1118  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1119  {
1120  Xyce::dout() << " before limiting: " << std::endl;
1121  Xyce::dout() << " vgs = " << vgs << " vgs_old = " << vgs_old << std::endl;
1122  Xyce::dout() << " vgd = " << vgd << " vgd_old = " << vgd_old << std::endl;
1123  }
1124 #endif
1125 
1126  ichk1=1;
1127  vgs = devSupport.pnjlim(vgs, vgs_old, vt, vcrit, &icheck);
1128  vgd = devSupport.pnjlim(vgd, vgd_old, vt, vcrit, &ichk1);
1129 
1130  if (ichk1 == 1) {icheck=1;}
1131  if (icheck == 1) limitedFlag=true;
1132 
1134  vgd = devSupport.fetlim(vgd, vgd_old, tvt0);
1135  vds = vgs-vgd;
1136 
1137 #ifdef Xyce_DEBUG_DEVICE
1138  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1139  {
1140  Xyce::dout() << " After limiting: " << std::endl;
1141  Xyce::dout() << " vgs = " << vgs << std::endl;
1142  Xyce::dout() << " vgd = " << vgd << std::endl;
1143  Xyce::dout() << " " << std::endl;
1144  }
1145 #endif
1146  }
1147 
1148 #ifdef Xyce_DEBUG_DEVICE
1149  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1150  {
1151  Xyce::dout() << "vgs = " << vgs << std::endl;
1152  Xyce::dout() << "vgd = " << vgd << std::endl;
1153  Xyce::dout() << "vds = " << vds << std::endl;
1154  Xyce::dout() << "Vddp = " << Vddp << std::endl;
1155  Xyce::dout() << "Vssp = " << Vssp << std::endl;
1156  Xyce::dout() << "Vgsp = " << Vgsp << std::endl;
1157  Xyce::dout() << "Vgdp = " << Vgdp << std::endl;
1158  Xyce::dout() << "Vdpsp = " << Vdpsp << std::endl;
1159  Xyce::dout() << " " << std::endl;
1160  }
1161 #endif
1162  // Now set the origFlag
1163  if (vgs_orig != vgs || vds_orig != vds || vgd_orig != vgd) origFlag = 0;
1164 
1165 #ifdef Xyce_DEBUG_DEVICE
1166  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1167  {
1168  if (origFlag == 0)
1169  {
1170  Xyce::dout() << " Something modified the voltages. " << std::endl;
1171  Xyce::dout() << " Voltage before after diff " << std::endl;
1172  Xyce::dout() << " vgs " << vgs_orig << " " << vgs << " " << vgs-vgs_orig << std::endl;
1173  Xyce::dout() << " vgd " << vgd_orig << " " << vgd << " " << vgd-vgd_orig << std::endl;
1174  Xyce::dout() << " vds " << vds_orig << " " << vds << " " << vds-vds_orig << std::endl;
1175  Xyce::dout() << " " << std::endl;
1176  }
1177  }
1178 #endif
1179 
1180  // update the "old" variables:
1181  if (getSolverState().newtonIter != 0 && getSolverState().newtonIter != newtonIterOld)
1182  {
1184  }
1185 
1186  //
1187  // the following block of code evaluates the dc current and its
1188  // derivatives and the charges associated with the gate and
1189  // channel
1190  //
1191 
1192  // vt set in updateTemperature
1193  vtf = 5.0*vt;
1194  csat = tIS;
1195  if (vgs <= -vtf)
1196  {
1197  ggs = -csat/vgs + getDeviceOptions().gmin;
1198  cg = ggs*vgs;
1199  }
1200  else
1201  {
1202  evgs = exp(vgs/vt);
1203  ggs = csat*evgs/vt + getDeviceOptions().gmin;
1204  cg = csat*(evgs-1) + getDeviceOptions().gmin*vgs;
1205  }
1206  if (vgd <= -vtf)
1207  {
1208  ggd = -csat/vgd + getDeviceOptions().gmin;
1209  cgd = ggd*vgd;
1210  }
1211  else
1212  {
1213  evgd = exp(vgd/vt);
1214  ggd = csat*evgd/vt + getDeviceOptions().gmin;
1215  cgd = csat*(evgd-1) + getDeviceOptions().gmin*vgd;
1216  }
1217  cg = cg + cgd;
1218 
1219  // 3f5 does this simple stuff
1220  if (vds >= 0)
1221  mode = 1;
1222  else
1223  mode = -1;
1224 
1225  if (vds >= 0) // normal mode
1226  {
1227  vgst = vgs-tvt0;
1228  if (vgst <= 0)
1229  {
1230  //
1231  // normal mode, cutoff region
1232  //
1233  cdrain = 0;
1234  gm = 0;
1235  gds = 0;
1236  }
1237  else
1238  {
1239  prod = 1 + tLambda*vds;
1240  betap = tBeta*prod;
1241  denom = 1 + tMESb*vgst;
1242  invdenom = 1/denom;
1243  if (vds >= ( 3/tAlpha ) )
1244  {
1245  //
1246  // normal mode, saturation region
1247  //
1248  cdrain = betap*vgst*vgst*invdenom;
1249  gm = betap*vgst*(1 + denom)*invdenom*invdenom;
1250  gds = tLambda*tBeta*vgst*vgst*invdenom;
1251  }
1252  else
1253  {
1254  //
1255  // normal mode, linear region
1256  //
1257  afact = 1 - tAlpha*vds/3;
1258  lfact = 1 - afact*afact*afact;
1259  cdrain = betap*vgst*vgst*invdenom*lfact;
1260  gm = betap*vgst*(1 + denom)*invdenom*invdenom*lfact;
1261  gds = tBeta*vgst*vgst*invdenom*(tAlpha*afact*afact*prod + lfact*tLambda);
1262  }
1263  }
1264  }
1265  else // inverse mode
1266  {
1267  vgdt = vgd - tvt0;
1268  if (vgdt <= 0)
1269  {
1270  //
1271  // inverse mode, cutoff region
1272  //
1273  cdrain = 0;
1274  gm = 0;
1275  gds = 0;
1276  }
1277  else
1278  {
1279  //
1280  // inverse mode, saturation region
1281  //
1282  prod = 1 - tLambda*vds;
1283  betap = tBeta*prod;
1284  denom = 1 + tMESb*vgdt;
1285  invdenom = 1/denom;
1286  if ( -vds >= 3/tAlpha )
1287  {
1288  cdrain = -betap*vgdt*vgdt*invdenom;
1289  gm = -betap*vgdt*(1 + denom)*invdenom*invdenom;
1290  gds = tLambda*tBeta*vgdt*vgdt*invdenom - gm;
1291  }
1292  else
1293  {
1294  //
1295  // inverse mode, linear region
1296  //
1297  afact = 1 + tAlpha*vds/3;
1298  lfact = 1 - afact*afact*afact;
1299  cdrain = -betap*vgdt*vgdt*invdenom*lfact;
1300  gm = -betap*vgdt*(1 + denom)*invdenom*invdenom*lfact;
1301  gds = tBeta*vgdt*vgdt*invdenom*(tAlpha*afact*afact*prod
1302  + lfact*tLambda) - gm;
1303  }
1304  }
1305  }
1306  cd = cdrain-cgd;
1307 
1308  //
1309  // charge storage elements
1310  //
1311  twop = 2.0*tPB;
1312  fcpb2 = corDepCap*corDepCap;
1313  czgs = tCGS;
1314  czgd = tCGD;
1315  if(czgs != 0)
1316  {
1317  czgsf2=czgs/f2;
1318  if (vgs < corDepCap)
1319  {
1320  sarg=sqrt(1-vgs/tPB);
1321  qgs = twop*czgs*(1-sarg);
1322  capgs=czgs/sarg;
1323  }
1324  else
1325  {
1326  qgs = czgs*f1 + czgsf2*(f3 *(vgs - corDepCap)
1327  +(vgs*vgs - fcpb2)/(2*twop));
1328  capgs=czgsf2*(f3 + vgs/twop);
1329  }
1330  }
1331  else
1332  {
1333  qgs=0.0;
1334  capgs=0.0;
1335  }
1336 
1337  if(czgd != 0)
1338  {
1339  czgdf2=czgd/f2;
1340  if (vgd < corDepCap)
1341  {
1342  sarg=sqrt(1-vgd/tPB);
1343  qgd = twop*czgd*(1-sarg);
1344  capgd=czgd/sarg;
1345  }
1346  else
1347  {
1348  qgd = czgd*f1 + czgdf2*( f3*(vgd - corDepCap)
1349  +(vgd*vgd - fcpb2)/(2*twop) );
1350  capgd=czgdf2*(f3 + vgd/twop);
1351  }
1352  }
1353  else
1354  {
1355  qgd=0.0;
1356  capgd=0.0;
1357  }
1358 
1359  Idrain = drainCond * Vddp;
1360  Isource = sourceCond * Vssp;
1361 
1362 #ifdef Xyce_DEBUG_DEVICE
1364  {
1365  Xyce::dout() << " Done with Instance::updateIntermediateVars." << std::endl;
1366  Xyce::dout() << " mode = " << mode << std::endl;
1367  Xyce::dout() << " tBeta = " << tBeta << std::endl;
1368  Xyce::dout() << " Idrain = " << Idrain << std::endl;
1369  Xyce::dout() << " Isource = " << Isource << std::endl;
1370  Xyce::dout() << " gds = " << gds << std::endl;
1371  Xyce::dout() << " gm = " << gm << std::endl;
1372  }
1373 #endif
1374 
1375  /// CURRENTS to load into RHS:
1376 
1377  // so at this point:
1378 
1379  // current out of drain is
1380  // Idrain
1381 
1382  // current out of gate:
1383  // dtype*( d/dt(qgs) + d/dt(qgd) )
1384 
1385  // the current *out of* the source should be simply
1386  // Isource
1387 
1388  // current out of drain' is
1389  // -Idrain - dtype*( d/dt(qgd) - cdrain )
1390 
1391  // the current out of the source' is
1392  // -Isource - dtype*( d/dt(qgs) + cdrain )
1393 
1394  return true;
1395 }
1396 
1397 //-----------------------------------------------------------------------------
1398 // Function : Instance::loadDAEQVector
1399 //
1400 // Purpose : Loads the Q-vector contributions for a single
1401 // voltage source instance.
1402 //
1403 // Special Notes : The "Q" vector is part of a standard DAE formalism in
1404 // which the system of equations is represented as:
1405 //
1406 // f(x) = dQ(x)/dt + F(x) + B(t) = 0
1407 //
1408 // The "Q" vector contains charges and fluxes, mostly.
1409 // The voltage source will not make any contributions to Q,
1410 // so this function does nothing.
1411 //
1412 // from updateSecondaryState:
1413 //
1414 // ggd = ggd + capgd*(getSolverState().pdt);
1415 // ggs = ggs + capgs*(getSolverState().pdt);
1416 //
1417 // // Sum the capacitor currents into the DC currents.
1418 // cg = cg + cqgs + cqgd;
1419 // cd = cd - cqgd;
1420 // cgd = cgd + cqgd;
1421 //
1422 // So:
1423 //
1424 // replace ggd with capgd.
1425 // replace ggs with capgs
1426 //
1427 // replace cg with qgs+qgd
1428 // replace cd with -qgd
1429 // replace cgd with qgd.
1430 //
1431 //
1432 // Scope : public
1433 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
1434 // Creation Date : 06/01/05
1435 //-----------------------------------------------------------------------------
1437 {
1438  double * qVec = extData.daeQVectorRawPtr;
1439  double * dQdxdVp = extData.dQdxdVpVectorRawPtr;
1440 
1441  // set up the final load variables:
1442  int Dtype = model_.dtype;
1443  double ceqgd = Dtype*(qgd);
1444  double ceqgs = Dtype*(((qgs+qgd)-qgd));
1445  double cdreq = Dtype*(((-qgd)+qgd));
1446 
1447  double ceqgd_Jdxp = -Dtype*(capgd*(vgd-vgd_orig));
1448  double ceqgs_Jdxp = -Dtype*(capgs*(vgs-vgs_orig));
1449  double cdreq_Jdxp = 0.0;
1450 
1451  qVec[li_Gate ] += ( ceqgs+ceqgd);
1452  qVec[li_DrainPrime ] -= (-cdreq+ceqgd);
1453  qVec[li_SourcePrime] -= ( cdreq+ceqgs);
1454 
1455  if (!origFlag)
1456  {
1457  dQdxdVp[li_Gate ] -= ( ceqgs_Jdxp+ceqgd_Jdxp);
1458  dQdxdVp[li_DrainPrime ] += (-cdreq_Jdxp+ceqgd_Jdxp);
1459  dQdxdVp[li_SourcePrime] += ( cdreq_Jdxp+ceqgs_Jdxp);
1460  }
1461 
1462  return true;
1463 }
1464 
1465 //-----------------------------------------------------------------------------
1466 // Function : Instance::loadDAEFVector
1467 //
1468 // Purpose : Loads the F-vector contributions for a single
1469 // MESFET instance.
1470 //
1471 // Special Notes :
1472 //
1473 // Scope : public
1474 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
1475 // Creation Date : 06/01/05
1476 //-----------------------------------------------------------------------------
1478 {
1479  double * fVec = extData.daeFVectorRawPtr;
1480  double * dFdxdVp = extData.dFdxdVpVectorRawPtr;
1481 
1482  // set up the final load variables:
1483  int Dtype = model_.dtype;
1484  double ceqgd = Dtype*(cgd);
1485  double ceqgs = Dtype*((cg-cgd));
1486  double cdreq = Dtype*((cd+cgd));
1487 
1488  double ceqgd_Jdxp = -Dtype*(ggd*(vgd-vgd_orig));
1489  double ceqgs_Jdxp = -Dtype*(ggs*(vgs-vgs_orig));
1490  double cdreq_Jdxp = -Dtype*(gds*(vds-vds_orig)+gm*(vgs-vgs_orig));
1491 
1492  // optional load resistors:
1493  if (drainCond != 0.0)
1494  {
1495  fVec[li_Drain ] += Idrain;
1496  }
1497  if (sourceCond != 0.0)
1498  {
1499  fVec[li_Source] += Isource;
1500  }
1501 
1502  fVec[li_Gate ] += (ceqgs+ceqgd);
1503  fVec[li_DrainPrime ] -= (Idrain +(-cdreq+ceqgd));
1504  fVec[li_SourcePrime] -= (Isource+(cdreq+ceqgs));
1505 
1506  if (!origFlag)
1507  {
1508  dFdxdVp[li_Gate ] -= ( ceqgs_Jdxp+ceqgd_Jdxp);
1509  dFdxdVp[li_DrainPrime ] += (-cdreq_Jdxp+ceqgd_Jdxp);
1510  dFdxdVp[li_SourcePrime] += ( cdreq_Jdxp+ceqgs_Jdxp);
1511  }
1512 
1513  return true;
1514 }
1515 
1516 //-----------------------------------------------------------------------------
1517 // Function : Instance::loadDAEdQdx
1518 //
1519 // Purpose : Loads the Q-vector contributions for a single
1520 // MESFET instance.
1521 //
1522 // Special Notes : The "Q" vector is part of a standard DAE formalism in
1523 // which the system of equations is represented as:
1524 //
1525 // f(x) = dQ(x)/dt + F(x) - B(t) = 0
1526 //
1527 // Scope : public
1528 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
1529 // Creation Date : 06/01/05
1530 //-----------------------------------------------------------------------------
1532 {
1533  N_LAS_Matrix & dQdx = *(extData.dQdxMatrixPtr);
1534 
1542 
1543  return true;
1544 }
1545 
1546 //-----------------------------------------------------------------------------
1547 // Function : Instance::loadDAEdFdx ()
1548 //
1549 // Purpose : Loads the F-vector contributions for a single
1550 // resistor instance.
1551 //
1552 // Special Notes : The F-vector is an algebraic constaint.
1553 //
1554 // Scope : public
1555 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
1556 // Creation Date : 06/01/05
1557 //-----------------------------------------------------------------------------
1559 {
1560  N_LAS_Matrix & dFdx = *(extData.dFdxMatrixPtr);
1561 
1564 
1568 
1571 
1575  drainCond+gds+ggd;
1577 
1582  += sourceCond+gds+gm+ggs;
1583 
1584  return true;
1585 }
1586 
1587 //-----------------------------------------------------------------------------
1588 // Function : Instance::updateTemperature
1589 // Purpose :
1590 // Special Notes :
1591 // Scope : public
1592 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1593 // Creation Date : 3/16/00
1594 //-----------------------------------------------------------------------------
1595 bool Instance::updateTemperature ( const double & temp_tmp)
1596 {
1597  bool bsuccess = true;
1598  double tnom, ratio;
1599  double arg, arg1;
1600  double ratio1;
1601  double fact1, fact2;
1602  double kt, kt1;
1603  double vtnom;
1604  double egfet, egfet1;
1605  double pbfact;
1606  double cjfact, cjfact1;
1607  double gmanew, gmaold;
1608  double pbo;
1609  double xfc;
1610  double Pb;
1611 
1612 #ifdef Xyce_DEBUG_DEVICE
1613  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1614  {
1615 // Xyce::dout() << subsection_divider << std::endl;
1616  Xyce::dout() << " Instance::Begin of updateTemperature. \n";
1617  Xyce::dout() << " name = " << getName() << std::endl;
1618  Xyce::dout() << std::endl;
1619  }
1620 #endif
1621 
1622  // first set the instance temperature to the new temperature:
1623  if (temp_tmp != -999.0) temp = temp_tmp;
1625  {
1626  // make sure interpolation doesn't take any resistance negative
1627  if(model_.RD < 0) model_.RD = 0;
1628  if(model_.RS < 0) model_.RS = 0;
1629 
1630  // some params may have changed during interpolation
1631  // model_.processParams();
1632  }
1633 
1634  Pb = model_.PB;
1635  tnom = model_.TNOM;
1636  ratio = temp/tnom;
1637 
1638  // first do the model stuff
1639 
1640  vtnom = tnom*CONSTKoverQ;
1641  fact1 = tnom/CONSTREFTEMP;
1642  kt1 = CONSTboltz*tnom;
1643  egfet1 = 1.16 - (7.02e-4*tnom*tnom)/(tnom + 1108);
1644  arg1 = -egfet1/(2.0*kt1) + 1.1150877/(CONSTboltz*2.0*CONSTREFTEMP);
1645  pbfact = -2.0*vtnom*(1.5*log(fact1) + CONSTQ*arg1);
1646  pbo = (Pb - pbfact)/fact1;
1647  gmaold = (Pb - pbo)/pbo;
1648  cjfact = 1.0/(1.0 + 0.5*(4e-4*(tnom - CONSTREFTEMP) - gmaold));
1649 
1650  if(model_.FC >.95) {
1651  Xyce::dout() << "Depletion cap. coeff. FC too large, limited to .95" <<
1652  Xyce::dout() << std::endl;
1653  model_.FC = .95;
1654  }
1655  xfc = log(1.0 - model_.FC);
1656  f2 = exp(1.5*xfc);
1657  f3 = 1.0 - 1.5*model_.FC;
1658  // skip bFac
1659 
1660  // now do the instance stuff
1661 
1662  vt = temp*CONSTKoverQ;
1663  kt = temp*CONSTboltz;
1664  fact2 = temp/CONSTREFTEMP;
1665  ratio1 = ratio - 1.0;
1666  tIS = model_.IS*exp(ratio1*1.11/vt)*area;
1667 
1668  tCGS = model_.CGS*cjfact*area;
1669  tCGD = model_.CGD*cjfact*area;
1670  egfet = 1.16 - (7.02e-4*temp*temp)/(temp + 1108);
1671  arg = -egfet/(2.0*kt) + 1.1150877/(CONSTboltz*2.0*CONSTREFTEMP);
1672  pbfact = -2.0*vt*(1.5*log(fact2) + CONSTQ*arg);
1673  tPB = fact2*pbo + pbfact;
1674  gmanew = (tPB - pbo)/pbo;
1675  cjfact1 = 1.0 + 0.5*(4e-4*(temp - CONSTREFTEMP) - gmanew);
1676  tCGS *= cjfact1;
1677  tCGD *= cjfact1;
1678 
1679  corDepCap = model_.FC*tPB;
1680  f1 = tPB*(1.0 - exp((0.5)*xfc))/(0.5);
1681  vcrit = vt * log(vt/(CONSTroot2 * tIS));
1682 
1683  // the following parameters have no temperature dependence in Spice 3f5
1684  //
1685  tBeta = model_.BETA*area; // transconductance parameter
1686  tvt0 = model_.VTO; // threshold voltage
1687  tLambda = model_.LAMBDA; // channel-length modulation
1688  tAlpha = model_.ALPHA; // saturation voltage parameter
1689  tRD = model_.RD/area; // drain ohmic resistance
1690  tRS = model_.RS/area; // source ohmic resistance
1691  tMESb = model_.B; // dopinng tail parameter
1692 
1693 #ifdef Xyce_DEBUG_DEVICE
1694  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1695  {
1696  Xyce::dout() << "temp = "<< temp << std::endl;
1697  Xyce::dout() << "tnom = " << tnom << std::endl;
1698  Xyce::dout() << "ratio = " << ratio << std::endl;
1699  Xyce::dout() << "vt = " << vt << std::endl;
1700  Xyce::dout() << "kt = " << kt << std::endl;
1701  Xyce::dout() << "fact2 = " << fact2 << std::endl;
1702  Xyce::dout() << "egfet = " << egfet << std::endl;
1703  Xyce::dout() << "arg = " << arg << std::endl;
1704  Xyce::dout() << "pbfact = " << pbfact << std::endl;
1705  Xyce::dout() << "PB = " << Pb << std::endl;
1706  Xyce::dout() << "pbo = " << pbo << std::endl;
1707  Xyce::dout() << "f2 = " << f2 << std::endl;
1708  Xyce::dout() << "f3 = " << f3 << std::endl;
1709  Xyce::dout() << "corDepCap= " << corDepCap << std::endl;
1710  Xyce::dout() << "tBeta = " << tBeta << std::endl;
1711  Xyce::dout() << "tvt0 = " << tvt0 << std::endl;
1712  Xyce::dout() << "tPB = " << tPB << std::endl;
1713  Xyce::dout() << "tMESb = " << tMESb << std::endl;
1714  Xyce::dout() << "tLambda = " << tLambda << std::endl;
1715  //Xyce::dout() << "tTheta = " << tTheta << std::endl;
1716  Xyce::dout() << "tRD = " << tRD << std::endl;
1717  Xyce::dout() << "tRS = " << tRS << std::endl;
1718  Xyce::dout() << " " << std::endl;
1719  }
1720 #endif
1721 
1722  return bsuccess;
1723 }
1724 
1725 //-----------------------------------------------------------------------------
1726 // Function : Instance::processParams
1727 // Purpose :
1728 // Special Notes :
1729 // Scope : public
1730 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1731 // Creation Date : 3/16/00
1732 //-----------------------------------------------------------------------------
1734 {
1736 
1737  return true;
1738 }
1739 
1740 // MESFET Master functions:
1741 
1742 //-----------------------------------------------------------------------------
1743 // Function : Master::updateState
1744 // Purpose :
1745 // Special Notes :
1746 // Scope : public
1747 // Creator : Eric Keiter, SNL
1748 // Creation Date : 11/26/08
1749 //-----------------------------------------------------------------------------
1750 bool Master::updateState (double * solVec, double * staVec, double * stoVec)
1751 {
1752  bool bsuccess = true;
1753 
1754  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
1755  {
1756  Instance & ji = *(*it);
1757 
1758  bool btmp = ji.updateIntermediateVars ();
1759  bsuccess = bsuccess && btmp;
1760 
1761  double * stoVec = ji.extData.nextStoVectorRawPtr;
1762  stoVec[ji.li_store_vgs] = ji.vgs;
1763  stoVec[ji.li_store_vgd] = ji.vgd;
1764  staVec[ji.li_state_qgs] = ji.qgs;
1765  staVec[ji.li_state_qgd] = ji.qgd;
1766  }
1767 
1768  return bsuccess;
1769 }
1770 
1771 //-----------------------------------------------------------------------------
1772 // Function : Master::loadDAEVectors
1773 // Purpose :
1774 // Special Notes :
1775 // Scope : public
1776 // Creator : Eric Keiter, SNL
1777 // Creation Date : 11/26/08
1778 //-----------------------------------------------------------------------------
1779 bool Master::loadDAEVectors (double * solVec, double * fVec, double *qVec, double * storeLeadF, double * storeLeadQ)
1780 {
1781  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
1782  {
1783  Instance & ji = *(*it);
1784 
1785  // F-vector:
1786  double * dFdxdVp = ji.extData.dFdxdVpVectorRawPtr;
1787 
1788  // set up the final load variables:
1789  int Dtype = ji.getModel().dtype;
1790  double f_ceqgd = Dtype*(ji.cgd);
1791  double f_ceqgs = Dtype*((ji.cg-ji.cgd));
1792  double f_cdreq = Dtype*((ji.cd+ji.cgd));
1793 
1794  double f_ceqgd_Jdxp = -Dtype*(ji.ggd*(ji.vgd-ji.vgd_orig));
1795  double f_ceqgs_Jdxp = -Dtype*(ji.ggs*(ji.vgs-ji.vgs_orig));
1796  double f_cdreq_Jdxp = -Dtype*(ji.gds*(ji.vds-ji.vds_orig)+ji.gm*(ji.vgs-ji.vgs_orig));
1797 
1798  // optional load resistors:
1799  if (ji.drainCond != 0.0)
1800  {
1801  fVec[ji.li_Drain ] += ji.Idrain;
1802  }
1803  if (ji.sourceCond != 0.0)
1804  {
1805  fVec[ji.li_Source] += ji.Isource;
1806  }
1807  fVec[ji.li_Gate ] += (f_ceqgs+f_ceqgd);
1808  fVec[ji.li_DrainPrime ] -= (ji.Idrain +(-f_cdreq+f_ceqgd));
1809  fVec[ji.li_SourcePrime] -= (ji.Isource+(f_cdreq+f_ceqgs));
1810 
1811  if (!ji.origFlag)
1812  {
1813  dFdxdVp[ji.li_Gate ] -= ( f_ceqgs_Jdxp+f_ceqgd_Jdxp);
1814  dFdxdVp[ji.li_DrainPrime ] += (-f_cdreq_Jdxp+f_ceqgd_Jdxp);
1815  dFdxdVp[ji.li_SourcePrime] += ( f_cdreq_Jdxp+f_ceqgs_Jdxp);
1816  }
1817 
1818  // Q-vector:
1819  double * dQdxdVp = ji.extData.dQdxdVpVectorRawPtr;
1820 
1821  // set up the final load variables:
1822  double q_ceqgd = Dtype*(ji.qgd);
1823  double q_ceqgs = Dtype*(((ji.qgs+ji.qgd)-ji.qgd));
1824  double q_cdreq = Dtype*(((-ji.qgd)+ji.qgd));
1825 
1826  double q_ceqgd_Jdxp = -Dtype*(ji.capgd*(ji.vgd-ji.vgd_orig));
1827  double q_ceqgs_Jdxp = -Dtype*(ji.capgs*(ji.vgs-ji.vgs_orig));
1828  double q_cdreq_Jdxp = 0.0;
1829 
1830  qVec[ji.li_Gate ] += ( q_ceqgs+q_ceqgd);
1831  qVec[ji.li_DrainPrime ] -= (-q_cdreq+q_ceqgd);
1832  qVec[ji.li_SourcePrime] -= ( q_cdreq+q_ceqgs);
1833 
1834  if (!ji.origFlag)
1835  {
1836  dQdxdVp[ji.li_Gate ] -= ( q_ceqgs_Jdxp+q_ceqgd_Jdxp);
1837  dQdxdVp[ji.li_DrainPrime ] += (-q_cdreq_Jdxp+q_ceqgd_Jdxp);
1838  dQdxdVp[ji.li_SourcePrime] += ( q_cdreq_Jdxp+q_ceqgs_Jdxp);
1839  }
1840 
1841  if( ji.loadLeadCurrent )
1842  {
1843  if (ji.drainCond != 0.0)
1844  {
1845  storeLeadF[ji.li_store_dev_id] = ji.Idrain;
1846  }
1847  else
1848  {
1849  storeLeadF[ji.li_store_dev_id] = -(ji.Idrain +(-f_cdreq+f_ceqgd));
1850  storeLeadQ[ji.li_store_dev_id] = -(-q_cdreq+q_ceqgd);
1851  }
1852  if (ji.sourceCond != 0.0)
1853  {
1854  storeLeadF[ji.li_store_dev_is] = ji.Isource;
1855  }
1856  else
1857  {
1858  storeLeadF[ji.li_store_dev_is] = -(ji.Isource+(f_cdreq+f_ceqgs));
1859  storeLeadQ[ji.li_store_dev_is] = -( q_cdreq+q_ceqgs);
1860  }
1861  storeLeadF[ji.li_store_dev_ig] = (f_ceqgs+f_ceqgd);
1862  storeLeadQ[ji.li_store_dev_ig] = (q_ceqgs+q_ceqgd);
1863  }
1864 
1865  }
1866 
1867  return true;
1868 }
1869 
1870 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
1871 //-----------------------------------------------------------------------------
1872 // Function : Master::loadDAEMatrices
1873 // Purpose :
1874 // Special Notes :
1875 // Scope : public
1876 // Creator : Eric Keiter, SNL
1877 // Creation Date : 12/12/08
1878 //-----------------------------------------------------------------------------
1879 bool Master::loadDAEMatrices (N_LAS_Matrix & dFdx, N_LAS_Matrix & dQdx)
1880 {
1881  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
1882  {
1883  Instance & ji = *(*it);
1884 
1885  // F-matrix:
1886 
1888 
1890 
1891 
1892  *ji.f_GateEquGateNodePtr += ji.ggd+ji.ggs;
1893 
1894  *ji.f_GateEquDrainPrimeNodePtr -= ji.ggd;
1895 
1896  *ji.f_GateEquSourcePrimeNodePtr -= ji.ggs;
1897 
1898 
1900 
1902 
1903 
1905 
1906  *ji.f_DrainPrimeEquGateNodePtr += ji.gm-ji.ggd;
1907 
1909 
1911 
1912 
1913  *ji.f_SourcePrimeEquGateNodePtr -= ji.gm+ji.ggs;
1914 
1916 
1918 
1920 
1921  // Q-matrix:
1922 
1923  *ji.q_GateEquGateNodePtr += ji.capgd+ji.capgs;
1924 
1926 
1928 
1930 
1932 
1934 
1936  }
1937 
1938  return true;
1939 }
1940 
1941 #else
1942 //-----------------------------------------------------------------------------
1943 // Function : Master::loadDAEMatrices
1944 // Purpose :
1945 // Special Notes :
1946 // Scope : public
1947 // Creator : Eric Keiter, SNL
1948 // Creation Date : 12/12/08
1949 //-----------------------------------------------------------------------------
1950 bool Master::loadDAEMatrices (N_LAS_Matrix & dFdx, N_LAS_Matrix & dQdx)
1951 {
1952  int sizeInstances = instanceContainer_.size();
1953  for (int i=0; i<sizeInstances; ++i)
1954  {
1955  Instance & ji = *(instanceContainer_.at(i));
1956 
1957  // dFdx matrix:
1958 
1959  dFdx[ji.li_Drain][ji.ADrainEquDrainNodeOffset] += ji.drainCond;
1960 
1962 
1963 
1964  dFdx[ji.li_Gate][ji.AGateEquGateNodeOffset] += ji.ggd+ji.ggs;
1965 
1966  dFdx[ji.li_Gate][ji.AGateEquDrainPrimeNodeOffset] -= ji.ggd;
1967 
1968  dFdx[ji.li_Gate][ji.AGateEquSourcePrimeNodeOffset] -= ji.ggs;
1969 
1970 
1971  dFdx[ji.li_Source][ji.ASourceEquSourceNodeOffset] += ji.sourceCond;
1972 
1974 
1975 
1977 
1978  dFdx[ji.li_DrainPrime][ji.ADrainPrimeEquGateNodeOffset] += ji.gm-ji.ggd;
1979 
1981  ji.drainCond+ji.gds+ji.ggd;
1982 
1984 
1985 
1986  dFdx[ji.li_SourcePrime][ji.ASourcePrimeEquGateNodeOffset] -= ji.gm+ji.ggs;
1987 
1989 
1991 
1993  += ji.sourceCond+ji.gds+ji.gm+ji.ggs;
1994 
1995  // dQdx matrix:
1996 
1997  dQdx[ji.li_Gate ][ji.AGateEquGateNodeOffset ] += ji.capgd+ji.capgs;
1998 
1999  dQdx[ji.li_Gate ][ji.AGateEquDrainPrimeNodeOffset ] -= ji.capgd;
2000 
2001  dQdx[ji.li_Gate ][ji.AGateEquSourcePrimeNodeOffset ] -= ji.capgs;
2002 
2003  dQdx[ji.li_DrainPrime ][ji.ADrainPrimeEquGateNodeOffset ] -= ji.capgd;
2004 
2006 
2008 
2010  }
2011 
2012  return true;
2013 }
2014 #endif
2015 
2016 Device *Traits::factory(const Configuration &configuration, const FactoryBlock &factory_block)
2017 {
2018 
2019  return new Master(configuration, factory_block, factory_block.solverState_, factory_block.deviceOptions_);
2020 }
2021 
2023 {
2025  .registerDevice("z", 1)
2026  .registerModelType("nmf", 1)
2027  .registerModelType("pmf", 1);
2028 }
2029 
2030 } // namespace MESFET
2031 } // namespace Device
2032 } // namespace Xyce