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.119 $
40 //
41 // Revision Date : $Date: 2014/05/22 17:40:29 $
42 //
43 // Current Owner : $Author: erkeite $
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  if (drainCond != 0.0)
687  {
688  intNameMap[ li_DrainPrime ] = spiceInternalName(getName(), "drainprime");
689  }
690 
691  if (sourceCond != 0.0)
692  {
693  intNameMap[ li_SourcePrime ] = spiceInternalName(getName(), "sourceprime");
694  }
695  }
696 
697  return intNameMap;
698 }
699 
700 
701 //-----------------------------------------------------------------------------
702 // Function : N_DEV_MESFETInstance::getStoreNameMap
703 // Purpose :
704 // Special Notes :
705 // Scope : public
706 // Creator : Richard Schiek, Electrical Systems Modeling
707 // Creation Date : 4/4/2013
708 //-----------------------------------------------------------------------------
709 std::map<int,std::string> & N_DEV_MESFETInstance::getStoreNameMap ()
710 {
711  // set up the internal name map, if it hasn't been already.
712  if( loadLeadCurrent && storeNameMap.empty ())
713  {
714  storeNameMap[ li_store_dev_id ] = spiceStoreName(getName(), "DEV_ID");
715  storeNameMap[ li_store_dev_is ] = spiceStoreName(getName(), "DEV_IS");
716  storeNameMap[ li_store_dev_ig ] = spiceStoreName(getName(), "DEV_IG");
717  }
718 
719  return storeNameMap;
720 }
721 
722 
723 //----------------------------------------------------------------------------
724 // Function : Instance::registerStateLIDs
725 // Purpose :
726 // Special Notes :
727 // Scope : public
728 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
729 // Creation Date : 3/16/00
730 //----------------------------------------------------------------------------
731 void Instance::registerStateLIDs(const std::vector<int> & staLIDVecRef)
732 {
733  AssertLIDs(staLIDVecRef.size() == numStateVars);
734 
735 #ifdef Xyce_DEBUG_DEVICE
736  if (getDeviceOptions().debugLevel > 0)
737  {
738  Xyce::dout() << std::endl;
739  Xyce::dout() << section_divider << std::endl;
740  Xyce::dout() << " In Instance::registerStateLIDs\n\n";
741  Xyce::dout() << " name = " << getName() << std::endl;
742  Xyce::dout() << " Number of State LIDs: " << numStateVars << std::endl;
743  }
744 #endif
745 
746  // Copy over the global ID lists:
747  staLIDVec = staLIDVecRef;
748 
749  int lid=0;
750  li_state_qgs = staLIDVec[lid++];
751  li_state_gcgs = staLIDVec[lid++];
752 
753  li_state_qgd = staLIDVec[lid++];
754  li_state_gcgd = staLIDVec[lid++];
755 
756 #ifdef Xyce_DEBUG_DEVICE
757  if (getDeviceOptions().debugLevel > 0)
758  {
759  Xyce::dout() << " State local indices:" << std::endl;
760  Xyce::dout() << std::endl;
761 
762  Xyce::dout() << " li_state_qgs = " << li_state_qgs << std::endl;
763  Xyce::dout() << " li_state_gcgs = " << li_state_gcgs;
764  Xyce::dout() << " li_state_qgd = " << li_state_qgd;
765  Xyce::dout() << " li_state_gcgd = " << li_state_gcgd << std::endl;;
766 
767  Xyce::dout() << section_divider << std::endl;
768  }
769 #endif
770 
771 }
772 
773 //----------------------------------------------------------------------------
774 // Function : Instance::registerStoreLIDs
775 // Purpose :
776 // Special Notes :
777 // Scope : public
778 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
779 // Creation Date : 12/9/11
780 //----------------------------------------------------------------------------
781 void Instance::registerStoreLIDs(const std::vector<int> & stoLIDVecRef)
782 {
783  AssertLIDs(stoLIDVecRef.size() == getNumStoreVars());
784 
785  // Copy over the global ID lists:
786  stoLIDVec = stoLIDVecRef;
787 
788  int lid=0;
789  li_store_vgs = stoLIDVec[lid++];
790  li_store_vgd = stoLIDVec[lid++];
791 
792  if( loadLeadCurrent )
793  {
794  li_store_dev_id = stoLIDVec[lid++];
795  li_store_dev_ig = stoLIDVec[lid++];
796  li_store_dev_is = stoLIDVec[lid++];
797  }
798 }
799 
800 //----------------------------------------------------------------------------
801 // Function : Instance::jacobianStamp
802 // Purpose :
803 // Special Notes :
804 // Scope : public
805 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
806 // Creation Date : 3/16/00
807 //----------------------------------------------------------------------------
808 const std::vector< std::vector<int> > & Instance::jacobianStamp() const
809 {
810  if( drainCond != 0.0 && sourceCond != 0.0 )
811  return jacStamp_DC_SC;
812  else if( drainCond != 0.0 && sourceCond == 0.0 )
813  return jacStamp_DC;
814  else if( drainCond == 0.0 && sourceCond != 0.0 )
815  return jacStamp_SC;
816  else if( drainCond == 0.0 && sourceCond == 0.0 )
817  return jacStamp;
818  else
819  return jacStamp;
820 }
821 
822 //----------------------------------------------------------------------------
823 // Function : Instance::registerJacLIDs
824 // Special Notes :
825 // Scope : public
826 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
827 // Creation Date : 3/16/00
828 //----------------------------------------------------------------------------
829 void Instance::registerJacLIDs( const std::vector< std::vector<int> > & jacLIDVec )
830 {
831  DeviceInstance::registerJacLIDs( jacLIDVec );
832  std::vector<int> map;
833  std::vector< std::vector<int> > map2;
834 
835  if (drainCond != 0.0)
836  {
837  if (sourceCond != 0.0)
838  {
839  map = jacMap_DC_SC;
840  map2 = jacMap2_DC_SC;
841  }
842  else
843  {
844  map = jacMap_DC;
845  map2 = jacMap2_DC;
846  }
847  }
848  else
849  {
850  if (sourceCond != 0.0)
851  {
852  map = jacMap_SC;
853  map2 = jacMap2_SC;
854  }
855  else
856  {
857  map = jacMap;
858  map2 = jacMap2;
859  }
860  }
861 
862  ADrainEquDrainNodeOffset = jacLIDVec[map[0]][map2[0][0]];
863  ADrainEquDrainPrimeNodeOffset = jacLIDVec[map[0]][map2[0][1]];
864 
865  AGateEquGateNodeOffset = jacLIDVec[map[1]][map2[1][0]];
866  AGateEquDrainPrimeNodeOffset = jacLIDVec[map[1]][map2[1][1]];
867  AGateEquSourcePrimeNodeOffset = jacLIDVec[map[1]][map2[1][2]];
868 
869  ASourceEquSourceNodeOffset = jacLIDVec[map[2]][map2[2][0]];
870  ASourceEquSourcePrimeNodeOffset = jacLIDVec[map[2]][map2[2][1]];
871 
872  ADrainPrimeEquDrainNodeOffset = jacLIDVec[map[3]][map2[3][0]];
873  ADrainPrimeEquGateNodeOffset = jacLIDVec[map[3]][map2[3][1]];
874  ADrainPrimeEquDrainPrimeNodeOffset = jacLIDVec[map[3]][map2[3][2]];
875  ADrainPrimeEquSourcePrimeNodeOffset = jacLIDVec[map[3]][map2[3][3]];
876 
877  ASourcePrimeEquGateNodeOffset = jacLIDVec[map[4]][map2[4][0]];
878  ASourcePrimeEquSourceNodeOffset = jacLIDVec[map[4]][map2[4][1]];
879  ASourcePrimeEquDrainPrimeNodeOffset = jacLIDVec[map[4]][map2[4][2]];
880  ASourcePrimeEquSourcePrimeNodeOffset = jacLIDVec[map[4]][map2[4][3]];
881 }
882 
883 //-----------------------------------------------------------------------------
884 // Function : Instance::setupPointers
885 // Purpose :
886 // Special Notes :
887 // Scope : public
888 // Creator : Eric Keiter, SNL
889 // Creation Date : 11/30/08
890 //-----------------------------------------------------------------------------
892 {
893 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
894  N_LAS_Matrix & dFdx = *(extData.dFdxMatrixPtr);
895  N_LAS_Matrix & dQdx = *(extData.dQdxMatrixPtr);
896 
897  // F-matrix:
900 
904 
907 
912 
917 
918  // Q-matrix:
921 
925 
928 
933 
938 
939 #endif
940 }
941 
942 //----------------------------------------------------------------------------
943 // Function : Instance::updatePrimaryState
944 // Purpose :
945 // Special Notes :
946 // Scope : public
947 // Creator : pmc
948 // Creation Date : 11/16/2003
949 //----------------------------------------------------------------------------
951 {
952  double * staVec = extData.nextStaVectorRawPtr;
953  bool bsuccess = updateIntermediateVars ();
954 
955  double * stoVec = extData.nextStoVectorRawPtr;
956  stoVec[li_store_vgs] = vgs;
957  stoVec[li_store_vgd] = vgd;
958  staVec[li_state_qgs] = qgs;
959  staVec[li_state_qgd] = qgd;
960 
961  return bsuccess;
962 }
963 
964 //-----------------------------------------------------------------------------
965 // Function : Instance::updateIntermediateVars
966 // Purpose :
967 // Special Notes :
968 // Scope : public
969 // Creator : pmc
970 // Creation Date : 11/16/2003
971 //-----------------------------------------------------------------------------
973 {
974  double * solVec = extData.nextSolVectorRawPtr;
975  double * currStaVec = extData.currStaVectorRawPtr;
976 
977  int dtype;
978  double csat, betap;
979  double vgst, vgdt;
980  double evgs, evgd;
981  double sarg, vtf;;
982 
983 // from the spice jfet
984  double czgd, czgs;
985  double czgdf2, czgsf2;
986  double fcpb2;
987  double twop;
988  int icheck, ichk1;
989 
990 // for the Shockley version
991 // double A, B, C, B12, C12, D, Vdsat;
992 // double delta;
993  double prod, denom, invdenom, afact, lfact;
994 
995 #ifdef Xyce_DEBUG_DEVICE
996  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
997  {
998  Xyce::dout() << subsection_divider << std::endl;
999  Xyce::dout() <<" Instance::updateIntermediateVars.\n"<<std::endl;
1000  Xyce::dout() <<" name = " << getName() << std::endl;
1001  Xyce::dout() <<" Model name = " << model_.getName() << std::endl;
1002  Xyce::dout() <<" dtype is " << model_.dtype << std::endl;
1003  Xyce::dout() << std::endl;
1004  Xyce::dout().width(25); Xyce::dout().precision(17); Xyce::dout().setf(std::ios::scientific);
1005  }
1006 #endif
1007 
1008  icheck = 1;
1009  dtype = model_.dtype;
1010 
1011  // we need our solution variables for any of this stuff
1012  Vd = 0.0;
1013  Vs = 0.0;
1014  Vg = 0.0;
1015  Vdp = 0.0;
1016  Vsp = 0.0;
1017 
1018  Vd = solVec[li_Drain];
1019  Vg = solVec[li_Gate];
1020  Vs = solVec[li_Source];
1021  Vsp = solVec[li_SourcePrime];
1022  Vdp = solVec[li_DrainPrime];
1023 
1024 #ifdef Xyce_DEBUG_DEVICE
1025  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1026  {
1027  Xyce::dout() << " " << std::endl;
1028  Xyce::dout() << " Vg = " << Vg << std::endl;
1029  Xyce::dout() << " Vd = " << Vd << std::endl;
1030  Xyce::dout() << " Vs = " << Vs << std::endl;
1031  Xyce::dout() << " Vdp = " << Vdp << std::endl;
1032  Xyce::dout() << " Vsp = " << Vsp << std::endl;
1033  }
1034 #endif
1035 
1036  // now we need voltage drops
1037  Vddp = Vd - Vdp;
1038  Vssp = Vs - Vsp;
1039  Vgsp = Vg - Vsp;
1040  Vgdp = Vg - Vdp;
1041  Vdpsp = Vdp - Vsp;
1042 
1043  // Now the things that the 3f5 code really uses
1044  vgs = dtype * Vgsp;
1045  vgd = dtype * Vgdp;
1046  vds = vgs-vgd;
1047 
1048  origFlag = 1;
1049  limitedFlag = false;
1050  vgs_orig = vgs;
1051  vgd_orig = vgd;
1052  vds_orig = vds;
1053 
1054  if (getSolverState().newtonIter == 0)
1055  {
1056  newtonIterOld=0;
1057  if (getSolverState().initJctFlag && getDeviceOptions().voltageLimiterFlag)
1058  {
1059  if (getSolverState().inputOPFlag)
1060  {
1061  N_LAS_Vector * flagSolVectorPtr = extData.flagSolVectorPtr;
1062  if ((*flagSolVectorPtr)[li_Drain] == 0 || (*flagSolVectorPtr)[li_Gate] == 0 ||
1063  (*flagSolVectorPtr)[li_Source] == 0 || (*flagSolVectorPtr)[li_SourcePrime] ||
1064  (*flagSolVectorPtr)[li_DrainPrime] )
1065  {
1066  vgs = 0;
1067  vgd = 0;
1068  vds = vgs-vgd;
1069  }
1070  }
1071  else
1072  {
1073  vgs = 0;
1074  vgd = 0;
1075  vds = vgs-vgd;
1076  }
1077  }
1078  if (!(getSolverState().dcopFlag)||(getSolverState().locaEnabledFlag && getSolverState().dcopFlag))
1079  {
1080  double * currStoVec = extData.currStoVectorRawPtr;
1081  vgs_old = currStoVec[li_store_vgs];
1082  vgd_old = currStoVec[li_store_vgd];
1083  }
1084  else
1085  { // there is no history
1086  vgs_old = vgs;
1087  vgd_old = vgd;
1088  }
1089  }
1090  else
1091  {
1092  double *stoVec = extData.nextStoVectorRawPtr;
1093  vgs_old = stoVec[li_store_vgs];
1094  vgd_old = stoVec[li_store_vgd];
1095  }
1096 
1097  // SPICE-type Voltage Limiting
1098  ///////////////////////////////
1099 
1100  if (getDeviceOptions().voltageLimiterFlag)
1101  {
1102 #ifdef Xyce_DEBUG_DEVICE
1103  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1104  {
1105  Xyce::dout() << " before limiting: " << std::endl;
1106  Xyce::dout() << " vgs = " << vgs << " vgs_old = " << vgs_old << std::endl;
1107  Xyce::dout() << " vgd = " << vgd << " vgd_old = " << vgd_old << std::endl;
1108  }
1109 #endif
1110 
1111  ichk1=1;
1112  vgs = devSupport.pnjlim(vgs, vgs_old, vt, vcrit, &icheck);
1113  vgd = devSupport.pnjlim(vgd, vgd_old, vt, vcrit, &ichk1);
1114 
1115  if (ichk1 == 1) {icheck=1;}
1116  if (icheck == 1) limitedFlag=true;
1117 
1119  vgd = devSupport.fetlim(vgd, vgd_old, tvt0);
1120  vds = vgs-vgd;
1121 
1122 #ifdef Xyce_DEBUG_DEVICE
1123  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1124  {
1125  Xyce::dout() << " After limiting: " << std::endl;
1126  Xyce::dout() << " vgs = " << vgs << std::endl;
1127  Xyce::dout() << " vgd = " << vgd << std::endl;
1128  Xyce::dout() << " " << std::endl;
1129  }
1130 #endif
1131  }
1132 
1133 #ifdef Xyce_DEBUG_DEVICE
1134  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1135  {
1136  Xyce::dout() << "vgs = " << vgs << std::endl;
1137  Xyce::dout() << "vgd = " << vgd << std::endl;
1138  Xyce::dout() << "vds = " << vds << std::endl;
1139  Xyce::dout() << "Vddp = " << Vddp << std::endl;
1140  Xyce::dout() << "Vssp = " << Vssp << std::endl;
1141  Xyce::dout() << "Vgsp = " << Vgsp << std::endl;
1142  Xyce::dout() << "Vgdp = " << Vgdp << std::endl;
1143  Xyce::dout() << "Vdpsp = " << Vdpsp << std::endl;
1144  Xyce::dout() << " " << std::endl;
1145  }
1146 #endif
1147  // Now set the origFlag
1148  if (vgs_orig != vgs || vds_orig != vds || vgd_orig != vgd) origFlag = 0;
1149 
1150 #ifdef Xyce_DEBUG_DEVICE
1151  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1152  {
1153  if (origFlag == 0)
1154  {
1155  Xyce::dout() << " Something modified the voltages. " << std::endl;
1156  Xyce::dout() << " Voltage before after diff " << std::endl;
1157  Xyce::dout() << " vgs " << vgs_orig << " " << vgs << " " << vgs-vgs_orig << std::endl;
1158  Xyce::dout() << " vgd " << vgd_orig << " " << vgd << " " << vgd-vgd_orig << std::endl;
1159  Xyce::dout() << " vds " << vds_orig << " " << vds << " " << vds-vds_orig << std::endl;
1160  Xyce::dout() << " " << std::endl;
1161  }
1162  }
1163 #endif
1164 
1165  // update the "old" variables:
1166  if (getSolverState().newtonIter != 0 && getSolverState().newtonIter != newtonIterOld)
1167  {
1169  }
1170 
1171  //
1172  // the following block of code evaluates the dc current and its
1173  // derivatives and the charges associated with the gate and
1174  // channel
1175  //
1176 
1177  // vt set in updateTemperature
1178  vtf = 5.0*vt;
1179  csat = tIS;
1180  if (vgs <= -vtf)
1181  {
1182  ggs = -csat/vgs + getDeviceOptions().gmin;
1183  cg = ggs*vgs;
1184  }
1185  else
1186  {
1187  evgs = exp(vgs/vt);
1188  ggs = csat*evgs/vt + getDeviceOptions().gmin;
1189  cg = csat*(evgs-1) + getDeviceOptions().gmin*vgs;
1190  }
1191  if (vgd <= -vtf)
1192  {
1193  ggd = -csat/vgd + getDeviceOptions().gmin;
1194  cgd = ggd*vgd;
1195  }
1196  else
1197  {
1198  evgd = exp(vgd/vt);
1199  ggd = csat*evgd/vt + getDeviceOptions().gmin;
1200  cgd = csat*(evgd-1) + getDeviceOptions().gmin*vgd;
1201  }
1202  cg = cg + cgd;
1203 
1204  // 3f5 does this simple stuff
1205  if (vds >= 0)
1206  mode = 1;
1207  else
1208  mode = -1;
1209 
1210  if (vds >= 0) // normal mode
1211  {
1212  vgst = vgs-tvt0;
1213  if (vgst <= 0)
1214  {
1215  //
1216  // normal mode, cutoff region
1217  //
1218  cdrain = 0;
1219  gm = 0;
1220  gds = 0;
1221  }
1222  else
1223  {
1224  prod = 1 + tLambda*vds;
1225  betap = tBeta*prod;
1226  denom = 1 + tMESb*vgst;
1227  invdenom = 1/denom;
1228  if (vds >= ( 3/tAlpha ) )
1229  {
1230  //
1231  // normal mode, saturation region
1232  //
1233  cdrain = betap*vgst*vgst*invdenom;
1234  gm = betap*vgst*(1 + denom)*invdenom*invdenom;
1235  gds = tLambda*tBeta*vgst*vgst*invdenom;
1236  }
1237  else
1238  {
1239  //
1240  // normal mode, linear region
1241  //
1242  afact = 1 - tAlpha*vds/3;
1243  lfact = 1 - afact*afact*afact;
1244  cdrain = betap*vgst*vgst*invdenom*lfact;
1245  gm = betap*vgst*(1 + denom)*invdenom*invdenom*lfact;
1246  gds = tBeta*vgst*vgst*invdenom*(tAlpha*afact*afact*prod + lfact*tLambda);
1247  }
1248  }
1249  }
1250  else // inverse mode
1251  {
1252  vgdt = vgd - tvt0;
1253  if (vgdt <= 0)
1254  {
1255  //
1256  // inverse mode, cutoff region
1257  //
1258  cdrain = 0;
1259  gm = 0;
1260  gds = 0;
1261  }
1262  else
1263  {
1264  //
1265  // inverse mode, saturation region
1266  //
1267  prod = 1 - tLambda*vds;
1268  betap = tBeta*prod;
1269  denom = 1 + tMESb*vgdt;
1270  invdenom = 1/denom;
1271  if ( -vds >= 3/tAlpha )
1272  {
1273  cdrain = -betap*vgdt*vgdt*invdenom;
1274  gm = -betap*vgdt*(1 + denom)*invdenom*invdenom;
1275  gds = tLambda*tBeta*vgdt*vgdt*invdenom - gm;
1276  }
1277  else
1278  {
1279  //
1280  // inverse mode, linear region
1281  //
1282  afact = 1 + tAlpha*vds/3;
1283  lfact = 1 - afact*afact*afact;
1284  cdrain = -betap*vgdt*vgdt*invdenom*lfact;
1285  gm = -betap*vgdt*(1 + denom)*invdenom*invdenom*lfact;
1286  gds = tBeta*vgdt*vgdt*invdenom*(tAlpha*afact*afact*prod
1287  + lfact*tLambda) - gm;
1288  }
1289  }
1290  }
1291  cd = cdrain-cgd;
1292 
1293  //
1294  // charge storage elements
1295  //
1296  twop = 2.0*tPB;
1297  fcpb2 = corDepCap*corDepCap;
1298  czgs = tCGS;
1299  czgd = tCGD;
1300  if(czgs != 0)
1301  {
1302  czgsf2=czgs/f2;
1303  if (vgs < corDepCap)
1304  {
1305  sarg=sqrt(1-vgs/tPB);
1306  qgs = twop*czgs*(1-sarg);
1307  capgs=czgs/sarg;
1308  }
1309  else
1310  {
1311  qgs = czgs*f1 + czgsf2*(f3 *(vgs - corDepCap)
1312  +(vgs*vgs - fcpb2)/(2*twop));
1313  capgs=czgsf2*(f3 + vgs/twop);
1314  }
1315  }
1316  else
1317  {
1318  qgs=0.0;
1319  capgs=0.0;
1320  }
1321 
1322  if(czgd != 0)
1323  {
1324  czgdf2=czgd/f2;
1325  if (vgd < corDepCap)
1326  {
1327  sarg=sqrt(1-vgd/tPB);
1328  qgd = twop*czgd*(1-sarg);
1329  capgd=czgd/sarg;
1330  }
1331  else
1332  {
1333  qgd = czgd*f1 + czgdf2*( f3*(vgd - corDepCap)
1334  +(vgd*vgd - fcpb2)/(2*twop) );
1335  capgd=czgdf2*(f3 + vgd/twop);
1336  }
1337  }
1338  else
1339  {
1340  qgd=0.0;
1341  capgd=0.0;
1342  }
1343 
1344  Idrain = drainCond * Vddp;
1345  Isource = sourceCond * Vssp;
1346 
1347 #ifdef Xyce_DEBUG_DEVICE
1349  {
1350  Xyce::dout() << " Done with Instance::updateIntermediateVars." << std::endl;
1351  Xyce::dout() << " mode = " << mode << std::endl;
1352  Xyce::dout() << " tBeta = " << tBeta << std::endl;
1353  Xyce::dout() << " Idrain = " << Idrain << std::endl;
1354  Xyce::dout() << " Isource = " << Isource << std::endl;
1355  Xyce::dout() << " gds = " << gds << std::endl;
1356  Xyce::dout() << " gm = " << gm << std::endl;
1357  }
1358 #endif
1359 
1360  /// CURRENTS to load into RHS:
1361 
1362  // so at this point:
1363 
1364  // current out of drain is
1365  // Idrain
1366 
1367  // current out of gate:
1368  // dtype*( d/dt(qgs) + d/dt(qgd) )
1369 
1370  // the current *out of* the source should be simply
1371  // Isource
1372 
1373  // current out of drain' is
1374  // -Idrain - dtype*( d/dt(qgd) - cdrain )
1375 
1376  // the current out of the source' is
1377  // -Isource - dtype*( d/dt(qgs) + cdrain )
1378 
1379  return true;
1380 }
1381 
1382 //-----------------------------------------------------------------------------
1383 // Function : Instance::loadDAEQVector
1384 //
1385 // Purpose : Loads the Q-vector contributions for a single
1386 // voltage source instance.
1387 //
1388 // Special Notes : The "Q" vector is part of a standard DAE formalism in
1389 // which the system of equations is represented as:
1390 //
1391 // f(x) = dQ(x)/dt + F(x) + B(t) = 0
1392 //
1393 // The "Q" vector contains charges and fluxes, mostly.
1394 // The voltage source will not make any contributions to Q,
1395 // so this function does nothing.
1396 //
1397 // from updateSecondaryState:
1398 //
1399 // ggd = ggd + capgd*(getSolverState().pdt);
1400 // ggs = ggs + capgs*(getSolverState().pdt);
1401 //
1402 // // Sum the capacitor currents into the DC currents.
1403 // cg = cg + cqgs + cqgd;
1404 // cd = cd - cqgd;
1405 // cgd = cgd + cqgd;
1406 //
1407 // So:
1408 //
1409 // replace ggd with capgd.
1410 // replace ggs with capgs
1411 //
1412 // replace cg with qgs+qgd
1413 // replace cd with -qgd
1414 // replace cgd with qgd.
1415 //
1416 //
1417 // Scope : public
1418 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
1419 // Creation Date : 06/01/05
1420 //-----------------------------------------------------------------------------
1422 {
1423  double * qVec = extData.daeQVectorRawPtr;
1424  double * dQdxdVp = extData.dQdxdVpVectorRawPtr;
1425 
1426  // set up the final load variables:
1427  int Dtype = model_.dtype;
1428  double ceqgd = Dtype*(qgd);
1429  double ceqgs = Dtype*(((qgs+qgd)-qgd));
1430  double cdreq = Dtype*(((-qgd)+qgd));
1431 
1432  double ceqgd_Jdxp = -Dtype*(capgd*(vgd-vgd_orig));
1433  double ceqgs_Jdxp = -Dtype*(capgs*(vgs-vgs_orig));
1434  double cdreq_Jdxp = 0.0;
1435 
1436  qVec[li_Gate ] += ( ceqgs+ceqgd);
1437  qVec[li_DrainPrime ] -= (-cdreq+ceqgd);
1438  qVec[li_SourcePrime] -= ( cdreq+ceqgs);
1439 
1440  if (!origFlag)
1441  {
1442  dQdxdVp[li_Gate ] -= ( ceqgs_Jdxp+ceqgd_Jdxp);
1443  dQdxdVp[li_DrainPrime ] += (-cdreq_Jdxp+ceqgd_Jdxp);
1444  dQdxdVp[li_SourcePrime] += ( cdreq_Jdxp+ceqgs_Jdxp);
1445  }
1446 
1447  return true;
1448 }
1449 
1450 //-----------------------------------------------------------------------------
1451 // Function : Instance::loadDAEFVector
1452 //
1453 // Purpose : Loads the F-vector contributions for a single
1454 // MESFET instance.
1455 //
1456 // Special Notes :
1457 //
1458 // Scope : public
1459 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
1460 // Creation Date : 06/01/05
1461 //-----------------------------------------------------------------------------
1463 {
1464  double * fVec = extData.daeFVectorRawPtr;
1465  double * dFdxdVp = extData.dFdxdVpVectorRawPtr;
1466 
1467  // set up the final load variables:
1468  int Dtype = model_.dtype;
1469  double ceqgd = Dtype*(cgd);
1470  double ceqgs = Dtype*((cg-cgd));
1471  double cdreq = Dtype*((cd+cgd));
1472 
1473  double ceqgd_Jdxp = -Dtype*(ggd*(vgd-vgd_orig));
1474  double ceqgs_Jdxp = -Dtype*(ggs*(vgs-vgs_orig));
1475  double cdreq_Jdxp = -Dtype*(gds*(vds-vds_orig)+gm*(vgs-vgs_orig));
1476 
1477  // optional load resistors:
1478  if (drainCond != 0.0)
1479  {
1480  fVec[li_Drain ] += Idrain;
1481  }
1482  if (sourceCond != 0.0)
1483  {
1484  fVec[li_Source] += Isource;
1485  }
1486 
1487  fVec[li_Gate ] += (ceqgs+ceqgd);
1488  fVec[li_DrainPrime ] -= (Idrain +(-cdreq+ceqgd));
1489  fVec[li_SourcePrime] -= (Isource+(cdreq+ceqgs));
1490 
1491  if (!origFlag)
1492  {
1493  dFdxdVp[li_Gate ] -= ( ceqgs_Jdxp+ceqgd_Jdxp);
1494  dFdxdVp[li_DrainPrime ] += (-cdreq_Jdxp+ceqgd_Jdxp);
1495  dFdxdVp[li_SourcePrime] += ( cdreq_Jdxp+ceqgs_Jdxp);
1496  }
1497 
1498  return true;
1499 }
1500 
1501 //-----------------------------------------------------------------------------
1502 // Function : Instance::loadDAEdQdx
1503 //
1504 // Purpose : Loads the Q-vector contributions for a single
1505 // MESFET instance.
1506 //
1507 // Special Notes : The "Q" vector is part of a standard DAE formalism in
1508 // which the system of equations is represented as:
1509 //
1510 // f(x) = dQ(x)/dt + F(x) - B(t) = 0
1511 //
1512 // Scope : public
1513 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
1514 // Creation Date : 06/01/05
1515 //-----------------------------------------------------------------------------
1517 {
1518  N_LAS_Matrix & dQdx = *(extData.dQdxMatrixPtr);
1519 
1527 
1528  return true;
1529 }
1530 
1531 //-----------------------------------------------------------------------------
1532 // Function : Instance::loadDAEdFdx ()
1533 //
1534 // Purpose : Loads the F-vector contributions for a single
1535 // resistor instance.
1536 //
1537 // Special Notes : The F-vector is an algebraic constaint.
1538 //
1539 // Scope : public
1540 // Creator : Eric R. Keiter, SNL, Parallel Computational Sciences
1541 // Creation Date : 06/01/05
1542 //-----------------------------------------------------------------------------
1544 {
1545  N_LAS_Matrix & dFdx = *(extData.dFdxMatrixPtr);
1546 
1549 
1553 
1556 
1560  drainCond+gds+ggd;
1562 
1567  += sourceCond+gds+gm+ggs;
1568 
1569  return true;
1570 }
1571 
1572 //-----------------------------------------------------------------------------
1573 // Function : Instance::updateTemperature
1574 // Purpose :
1575 // Special Notes :
1576 // Scope : public
1577 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1578 // Creation Date : 3/16/00
1579 //-----------------------------------------------------------------------------
1580 bool Instance::updateTemperature ( const double & temp_tmp)
1581 {
1582  bool bsuccess = true;
1583  double tnom, ratio;
1584  double arg, arg1;
1585  double ratio1;
1586  double fact1, fact2;
1587  double kt, kt1;
1588  double vtnom;
1589  double egfet, egfet1;
1590  double pbfact;
1591  double cjfact, cjfact1;
1592  double gmanew, gmaold;
1593  double pbo;
1594  double xfc;
1595  double Pb;
1596 
1597 #ifdef Xyce_DEBUG_DEVICE
1598  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1599  {
1600 // Xyce::dout() << subsection_divider << std::endl;
1601  Xyce::dout() << " Instance::Begin of updateTemperature. \n";
1602  Xyce::dout() << " name = " << getName() << std::endl;
1603  Xyce::dout() << std::endl;
1604  }
1605 #endif
1606 
1607  // first set the instance temperature to the new temperature:
1608  if (temp_tmp != -999.0) temp = temp_tmp;
1610  {
1611  // make sure interpolation doesn't take any resistance negative
1612  if(model_.RD < 0) model_.RD = 0;
1613  if(model_.RS < 0) model_.RS = 0;
1614 
1615  // some params may have changed during interpolation
1616  // model_.processParams();
1617  }
1618 
1619  Pb = model_.PB;
1620  tnom = model_.TNOM;
1621  ratio = temp/tnom;
1622 
1623  // first do the model stuff
1624 
1625  vtnom = tnom*CONSTKoverQ;
1626  fact1 = tnom/CONSTREFTEMP;
1627  kt1 = CONSTboltz*tnom;
1628  egfet1 = 1.16 - (7.02e-4*tnom*tnom)/(tnom + 1108);
1629  arg1 = -egfet1/(2.0*kt1) + 1.1150877/(CONSTboltz*2.0*CONSTREFTEMP);
1630  pbfact = -2.0*vtnom*(1.5*log(fact1) + CONSTQ*arg1);
1631  pbo = (Pb - pbfact)/fact1;
1632  gmaold = (Pb - pbo)/pbo;
1633  cjfact = 1.0/(1.0 + 0.5*(4e-4*(tnom - CONSTREFTEMP) - gmaold));
1634 
1635  if(model_.FC >.95) {
1636  Xyce::dout() << "Depletion cap. coeff. FC too large, limited to .95" <<
1637  Xyce::dout() << std::endl;
1638  model_.FC = .95;
1639  }
1640  xfc = log(1.0 - model_.FC);
1641  f2 = exp(1.5*xfc);
1642  f3 = 1.0 - 1.5*model_.FC;
1643  // skip bFac
1644 
1645  // now do the instance stuff
1646 
1647  vt = temp*CONSTKoverQ;
1648  kt = temp*CONSTboltz;
1649  fact2 = temp/CONSTREFTEMP;
1650  ratio1 = ratio - 1.0;
1651  tIS = model_.IS*exp(ratio1*1.11/vt)*area;
1652 
1653  tCGS = model_.CGS*cjfact*area;
1654  tCGD = model_.CGD*cjfact*area;
1655  egfet = 1.16 - (7.02e-4*temp*temp)/(temp + 1108);
1656  arg = -egfet/(2.0*kt) + 1.1150877/(CONSTboltz*2.0*CONSTREFTEMP);
1657  pbfact = -2.0*vt*(1.5*log(fact2) + CONSTQ*arg);
1658  tPB = fact2*pbo + pbfact;
1659  gmanew = (tPB - pbo)/pbo;
1660  cjfact1 = 1.0 + 0.5*(4e-4*(temp - CONSTREFTEMP) - gmanew);
1661  tCGS *= cjfact1;
1662  tCGD *= cjfact1;
1663 
1664  corDepCap = model_.FC*tPB;
1665  f1 = tPB*(1.0 - exp((0.5)*xfc))/(0.5);
1666  vcrit = vt * log(vt/(CONSTroot2 * tIS));
1667 
1668  // the following parameters have no temperature dependence in Spice 3f5
1669  //
1670  tBeta = model_.BETA*area; // transconductance parameter
1671  tvt0 = model_.VTO; // threshold voltage
1672  tLambda = model_.LAMBDA; // channel-length modulation
1673  tAlpha = model_.ALPHA; // saturation voltage parameter
1674  tRD = model_.RD/area; // drain ohmic resistance
1675  tRS = model_.RS/area; // source ohmic resistance
1676  tMESb = model_.B; // dopinng tail parameter
1677 
1678 #ifdef Xyce_DEBUG_DEVICE
1679  if (getDeviceOptions().debugLevel > 0 && getSolverState().debugTimeFlag)
1680  {
1681  Xyce::dout() << "temp = "<< temp << std::endl;
1682  Xyce::dout() << "tnom = " << tnom << std::endl;
1683  Xyce::dout() << "ratio = " << ratio << std::endl;
1684  Xyce::dout() << "vt = " << vt << std::endl;
1685  Xyce::dout() << "kt = " << kt << std::endl;
1686  Xyce::dout() << "fact2 = " << fact2 << std::endl;
1687  Xyce::dout() << "egfet = " << egfet << std::endl;
1688  Xyce::dout() << "arg = " << arg << std::endl;
1689  Xyce::dout() << "pbfact = " << pbfact << std::endl;
1690  Xyce::dout() << "PB = " << Pb << std::endl;
1691  Xyce::dout() << "pbo = " << pbo << std::endl;
1692  Xyce::dout() << "f2 = " << f2 << std::endl;
1693  Xyce::dout() << "f3 = " << f3 << std::endl;
1694  Xyce::dout() << "corDepCap= " << corDepCap << std::endl;
1695  Xyce::dout() << "tBeta = " << tBeta << std::endl;
1696  Xyce::dout() << "tvt0 = " << tvt0 << std::endl;
1697  Xyce::dout() << "tPB = " << tPB << std::endl;
1698  Xyce::dout() << "tMESb = " << tMESb << std::endl;
1699  Xyce::dout() << "tLambda = " << tLambda << std::endl;
1700  //Xyce::dout() << "tTheta = " << tTheta << std::endl;
1701  Xyce::dout() << "tRD = " << tRD << std::endl;
1702  Xyce::dout() << "tRS = " << tRS << std::endl;
1703  Xyce::dout() << " " << std::endl;
1704  }
1705 #endif
1706 
1707  return bsuccess;
1708 }
1709 
1710 //-----------------------------------------------------------------------------
1711 // Function : Instance::processParams
1712 // Purpose :
1713 // Special Notes :
1714 // Scope : public
1715 // Creator : Eric Keiter, SNL, Parallel Computational Sciences
1716 // Creation Date : 3/16/00
1717 //-----------------------------------------------------------------------------
1719 {
1721 
1722  return true;
1723 }
1724 
1725 // MESFET Master functions:
1726 
1727 //-----------------------------------------------------------------------------
1728 // Function : Master::updateState
1729 // Purpose :
1730 // Special Notes :
1731 // Scope : public
1732 // Creator : Eric Keiter, SNL
1733 // Creation Date : 11/26/08
1734 //-----------------------------------------------------------------------------
1735 bool Master::updateState (double * solVec, double * staVec, double * stoVec)
1736 {
1737  bool bsuccess = true;
1738 
1739  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
1740  {
1741  Instance & ji = *(*it);
1742 
1743  bool btmp = ji.updateIntermediateVars ();
1744  bsuccess = bsuccess && btmp;
1745 
1746  double * stoVec = ji.extData.nextStoVectorRawPtr;
1747  stoVec[ji.li_store_vgs] = ji.vgs;
1748  stoVec[ji.li_store_vgd] = ji.vgd;
1749  staVec[ji.li_state_qgs] = ji.qgs;
1750  staVec[ji.li_state_qgd] = ji.qgd;
1751  }
1752 
1753  return bsuccess;
1754 }
1755 
1756 //-----------------------------------------------------------------------------
1757 // Function : Master::loadDAEVectors
1758 // Purpose :
1759 // Special Notes :
1760 // Scope : public
1761 // Creator : Eric Keiter, SNL
1762 // Creation Date : 11/26/08
1763 //-----------------------------------------------------------------------------
1764 bool Master::loadDAEVectors (double * solVec, double * fVec, double *qVec, double * bVec, double * storeLeadF, double * storeLeadQ)
1765 {
1766  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
1767  {
1768  Instance & ji = *(*it);
1769 
1770  // F-vector:
1771  double * dFdxdVp = ji.extData.dFdxdVpVectorRawPtr;
1772 
1773  // set up the final load variables:
1774  int Dtype = ji.getModel().dtype;
1775  double f_ceqgd = Dtype*(ji.cgd);
1776  double f_ceqgs = Dtype*((ji.cg-ji.cgd));
1777  double f_cdreq = Dtype*((ji.cd+ji.cgd));
1778 
1779  double f_ceqgd_Jdxp = -Dtype*(ji.ggd*(ji.vgd-ji.vgd_orig));
1780  double f_ceqgs_Jdxp = -Dtype*(ji.ggs*(ji.vgs-ji.vgs_orig));
1781  double f_cdreq_Jdxp = -Dtype*(ji.gds*(ji.vds-ji.vds_orig)+ji.gm*(ji.vgs-ji.vgs_orig));
1782 
1783  // optional load resistors:
1784  if (ji.drainCond != 0.0)
1785  {
1786  fVec[ji.li_Drain ] += ji.Idrain;
1787  }
1788  if (ji.sourceCond != 0.0)
1789  {
1790  fVec[ji.li_Source] += ji.Isource;
1791  }
1792  fVec[ji.li_Gate ] += (f_ceqgs+f_ceqgd);
1793  fVec[ji.li_DrainPrime ] -= (ji.Idrain +(-f_cdreq+f_ceqgd));
1794  fVec[ji.li_SourcePrime] -= (ji.Isource+(f_cdreq+f_ceqgs));
1795 
1796  if (!ji.origFlag)
1797  {
1798  dFdxdVp[ji.li_Gate ] -= ( f_ceqgs_Jdxp+f_ceqgd_Jdxp);
1799  dFdxdVp[ji.li_DrainPrime ] += (-f_cdreq_Jdxp+f_ceqgd_Jdxp);
1800  dFdxdVp[ji.li_SourcePrime] += ( f_cdreq_Jdxp+f_ceqgs_Jdxp);
1801  }
1802 
1803  // Q-vector:
1804  double * dQdxdVp = ji.extData.dQdxdVpVectorRawPtr;
1805 
1806  // set up the final load variables:
1807  double q_ceqgd = Dtype*(ji.qgd);
1808  double q_ceqgs = Dtype*(((ji.qgs+ji.qgd)-ji.qgd));
1809  double q_cdreq = Dtype*(((-ji.qgd)+ji.qgd));
1810 
1811  double q_ceqgd_Jdxp = -Dtype*(ji.capgd*(ji.vgd-ji.vgd_orig));
1812  double q_ceqgs_Jdxp = -Dtype*(ji.capgs*(ji.vgs-ji.vgs_orig));
1813  double q_cdreq_Jdxp = 0.0;
1814 
1815  qVec[ji.li_Gate ] += ( q_ceqgs+q_ceqgd);
1816  qVec[ji.li_DrainPrime ] -= (-q_cdreq+q_ceqgd);
1817  qVec[ji.li_SourcePrime] -= ( q_cdreq+q_ceqgs);
1818 
1819  if (!ji.origFlag)
1820  {
1821  dQdxdVp[ji.li_Gate ] -= ( q_ceqgs_Jdxp+q_ceqgd_Jdxp);
1822  dQdxdVp[ji.li_DrainPrime ] += (-q_cdreq_Jdxp+q_ceqgd_Jdxp);
1823  dQdxdVp[ji.li_SourcePrime] += ( q_cdreq_Jdxp+q_ceqgs_Jdxp);
1824  }
1825 
1826  if( ji.loadLeadCurrent )
1827  {
1828  if (ji.drainCond != 0.0)
1829  {
1830  storeLeadF[ji.li_store_dev_id] = ji.Idrain;
1831  }
1832  else
1833  {
1834  storeLeadF[ji.li_store_dev_id] = -(ji.Idrain +(-f_cdreq+f_ceqgd));
1835  storeLeadQ[ji.li_store_dev_id] = -(-q_cdreq+q_ceqgd);
1836  }
1837  if (ji.sourceCond != 0.0)
1838  {
1839  storeLeadF[ji.li_store_dev_is] = ji.Isource;
1840  }
1841  else
1842  {
1843  storeLeadF[ji.li_store_dev_is] = -(ji.Isource+(f_cdreq+f_ceqgs));
1844  storeLeadQ[ji.li_store_dev_is] = -( q_cdreq+q_ceqgs);
1845  }
1846  storeLeadF[ji.li_store_dev_ig] = (f_ceqgs+f_ceqgd);
1847  storeLeadQ[ji.li_store_dev_ig] = (q_ceqgs+q_ceqgd);
1848  }
1849 
1850  }
1851 
1852  return true;
1853 }
1854 
1855 #ifndef Xyce_NONPOINTER_MATRIX_LOAD
1856 //-----------------------------------------------------------------------------
1857 // Function : Master::loadDAEMatrices
1858 // Purpose :
1859 // Special Notes :
1860 // Scope : public
1861 // Creator : Eric Keiter, SNL
1862 // Creation Date : 12/12/08
1863 //-----------------------------------------------------------------------------
1864 bool Master::loadDAEMatrices (N_LAS_Matrix & dFdx, N_LAS_Matrix & dQdx)
1865 {
1866  for (InstanceVector::const_iterator it = getInstanceBegin(); it != getInstanceEnd(); ++it)
1867  {
1868  Instance & ji = *(*it);
1869 
1870  // F-matrix:
1871 
1873 
1875 
1876 
1877  *ji.f_GateEquGateNodePtr += ji.ggd+ji.ggs;
1878 
1879  *ji.f_GateEquDrainPrimeNodePtr -= ji.ggd;
1880 
1881  *ji.f_GateEquSourcePrimeNodePtr -= ji.ggs;
1882 
1883 
1885 
1887 
1888 
1890 
1891  *ji.f_DrainPrimeEquGateNodePtr += ji.gm-ji.ggd;
1892 
1894 
1896 
1897 
1898  *ji.f_SourcePrimeEquGateNodePtr -= ji.gm+ji.ggs;
1899 
1901 
1903 
1905 
1906  // Q-matrix:
1907 
1908  *ji.q_GateEquGateNodePtr += ji.capgd+ji.capgs;
1909 
1911 
1913 
1915 
1917 
1919 
1921  }
1922 
1923  return true;
1924 }
1925 
1926 #else
1927 //-----------------------------------------------------------------------------
1928 // Function : Master::loadDAEMatrices
1929 // Purpose :
1930 // Special Notes :
1931 // Scope : public
1932 // Creator : Eric Keiter, SNL
1933 // Creation Date : 12/12/08
1934 //-----------------------------------------------------------------------------
1935 bool Master::loadDAEMatrices (N_LAS_Matrix & dFdx, N_LAS_Matrix & dQdx)
1936 {
1937  int sizeInstances = instanceContainer_.size();
1938  for (int i=0; i<sizeInstances; ++i)
1939  {
1940  Instance & ji = *(instanceContainer_.at(i));
1941 
1942  // dFdx matrix:
1943 
1944  dFdx[ji.li_Drain][ji.ADrainEquDrainNodeOffset] += ji.drainCond;
1945 
1947 
1948 
1949  dFdx[ji.li_Gate][ji.AGateEquGateNodeOffset] += ji.ggd+ji.ggs;
1950 
1951  dFdx[ji.li_Gate][ji.AGateEquDrainPrimeNodeOffset] -= ji.ggd;
1952 
1953  dFdx[ji.li_Gate][ji.AGateEquSourcePrimeNodeOffset] -= ji.ggs;
1954 
1955 
1956  dFdx[ji.li_Source][ji.ASourceEquSourceNodeOffset] += ji.sourceCond;
1957 
1959 
1960 
1962 
1963  dFdx[ji.li_DrainPrime][ji.ADrainPrimeEquGateNodeOffset] += ji.gm-ji.ggd;
1964 
1966  ji.drainCond+ji.gds+ji.ggd;
1967 
1969 
1970 
1971  dFdx[ji.li_SourcePrime][ji.ASourcePrimeEquGateNodeOffset] -= ji.gm+ji.ggs;
1972 
1974 
1976 
1978  += ji.sourceCond+ji.gds+ji.gm+ji.ggs;
1979 
1980  // dQdx matrix:
1981 
1982  dQdx[ji.li_Gate ][ji.AGateEquGateNodeOffset ] += ji.capgd+ji.capgs;
1983 
1984  dQdx[ji.li_Gate ][ji.AGateEquDrainPrimeNodeOffset ] -= ji.capgd;
1985 
1986  dQdx[ji.li_Gate ][ji.AGateEquSourcePrimeNodeOffset ] -= ji.capgs;
1987 
1988  dQdx[ji.li_DrainPrime ][ji.ADrainPrimeEquGateNodeOffset ] -= ji.capgd;
1989 
1991 
1993 
1995  }
1996 
1997  return true;
1998 }
1999 #endif
2000 
2001 Device *Traits::factory(const Configuration &configuration, const FactoryBlock &factory_block)
2002 {
2003 
2004  return new Master(configuration, factory_block, factory_block.solverState_, factory_block.deviceOptions_);
2005 }
2006 
2008 {
2010  .registerDevice("z", 1)
2011  .registerModelType("nmf", 1)
2012  .registerModelType("pmf", 1);
2013 }
2014 
2015 } // namespace MESFET
2016 } // namespace Device
2017 } // namespace Xyce